Manybody_Executable_Template.h

#ifndef QNANO_NEW_MANYBODY_EXECUTABLE_TEMPLATE_DEFINED_H
#define QNANO_NEW_MANYBODY_EXECUTABLE_TEMPLATE_DEFINED_H

#include "manybody/Manybody_Hamiltonian.h"
#include "tools/Abstract_Solver.h"
#include "tools/Eigenvector_Printer.h"
#include "tools/Parameter_Map.h"

#include "manybody/Manybody_Perturbative_Correction_Template.h"

template <class Tsolver> void Postprocess_Manybody_Template(Tsolver *solver, Manybody_Hamiltonian *hamil){
}


template <class Tsolver, class Tprinter> 
void Manybody_Executable_Template(Parameter_Map &param_map, bool verbose){

  if(verbose)std::cout<<"Setting up manybody Hamiltonian..."<<std::endl;

  Manybody_Hamiltonian *hamil=new Manybody_Hamiltonian;
  hamil->setup(param_map, verbose);

  if(verbose)std::cout<<"Dimension: "<<hamil->get_DIM()<<std::endl;
  //TODO: Test if dimension smaller than number of threads
  // if(hamil->get_DIM()<
  bool dump=param_map.get_as_bool("dump",false);
  if(verbose && dump){
    std::cout<<"Dumping Hamiltonian"<<std::endl;
    hamil->dump_sparse_Matrix("hamil_IA.dat", "hamil_JA.dat", "hamil_values.dat");
  }
  
  if(verbose)std::cout<<"Setting up Eigen Solver"<<std::endl;
  Tsolver *solver=new Tsolver(hamil);
  solver->setup();

  if(verbose)std::cout<<"Solving"<<std::endl;
  solver->solve();

  if(verbose)std::cout<<"Printing"<<std::endl;
  if(!param_map.is_specified("eigenvalue_filename")){
    param_map.add("eigenvalue_filename","energies.dat");
  }



  Eigenvector_Printer *EVprinter=new Tprinter(solver, hamil);
  EVprinter->setup(param_map);
  EVprinter->print_sparse=true;
  EVprinter->sparse_print_index=Manybody_Hilbert_Space::print_index;
  Manybody_Hilbert_Space *fmhs=hamil->hilbert;
  EVprinter->spi_ctx=(void*)( fmhs);
  {
   std::stringstream ss; 
   ss<<"rank_e "<<fmhs->get_rank_e()<<std::endl;
   ss<<"rank_h "<<fmhs->get_rank_h()<<std::endl;
   ss<<"nr_e "<<fmhs->get_nr_e()<<std::endl;
   ss<<"nr_h "<<fmhs->get_nr_h()<<std::endl;
   ss<<"eigenvalue $EIGENVALUE"<<std::endl;
   ss<<"eigenvector:"<<std::endl;
   EVprinter->EV_header=ss.str();
  }

  EVprinter->print();

  Calculate_Manybody_Perturbation<Tsolver>(param_map,solver,EVprinter,hamil);

  Postprocess_Manybody_Template<Tsolver>(solver,hamil);

  delete EVprinter;
  delete solver;
  delete hamil;
  if(verbose)std::cout<<"Done"<<std::endl;
}
#endif