HPhi++  3.1.0
HPhiMain.cpp File Reference
#include <sz.hpp>
#include <HPhiTrans.hpp>
#include <output_list.hpp>
#include <diagonalcalc.hpp>
#include <CalcByLOBPCG.hpp>
#include <CalcByFullDiag.hpp>
#include <CalcByTPQ.hpp>
#include <CalcSpectrum.hpp>
#include <check.hpp>
#include "CalcByTEM.hpp"
#include "readdef.hpp"
#include "StdFace_main.hpp"
#include "wrapperMPI.hpp"
#include "splash.hpp"
#include "CalcTime.hpp"
#include "common/setmemory.hpp"

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 Main program for HPhi. More...
 

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

Main program for HPhi.

Parameters
argc[in] argument count
argv[in] argument vector
Version
2.1 Add Time evolution mode.
1.2 Add calculation spectrum mode.
1.0
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)
Return values
-1fail the calculation.
0succeed the calculation.

Definition at line 177 of file HPhiMain.cpp.

References EDMainCalStruct::Bind, BindStruct::Boost, CalcByFullDiag(), CalcByLOBPCG(), CalcByTEM(), CalcByTPQ(), CalcSpectrum(), check(), BindStruct::Def, diagonalcalc(), exitMPI(), ParamList::ExpecInterval, FinalizeMPI(), HPhiTrans(), DefineList::iCalcType, DefineList::iFlgCalcSpec, DefineList::iFlgScaLAPACK, InitializeMPI(), InitTimer(), JudgeDefType(), DefineList::k_exct, list_1, list_2_1, list_2_2, myrank, nproc, NumAve, DefineList::nvec, output_list(), OutputTimer(), DefineList::Param, ReadDefFileIdxPara(), ReadDefFileNInt(), SetConvergenceFactor(), setmem_def(), setmem_HEAD(), setmem_large(), splash(), StartTimer(), StdFace_main(), stdoutMPI, StopTimer(), sz(), TimeKeeper(), and DefineList::WRITE.

177  {
178 
179  int mode=0;
180  char cFileListName[D_FileNameMax];
181  struct EDMainCalStruct X;
182 
183  stdoutMPI = stdout;
184  if(JudgeDefType(argc, argv, &mode)!=0){
185  exitMPI(-1);
186  }
187 
188  if (mode == STANDARD_DRY_MODE) {
189  myrank = 0;
190  nproc = 1;
191  stdoutMPI = stdout;
192  splash();
193  }
194  else InitializeMPI(argc, argv);
195 
196  //Timer
197  InitTimer();
198  if (mode != STANDARD_DRY_MODE) StartTimer(0);
199 
200  //MakeDirectory for output
201  struct stat tmpst;
202  if (myrank == 0) {
203  if (stat("./output/", &tmpst) != 0) {
204  if (mkdir("./output/", 0777) != 0) {
205  fprintf(stdoutMPI, "%s", "Error: Fail to make output folder in current directory. \n");
206  exitMPI(-1);
207  }
208  }
209  }/*if (myrank == 0)*/
210 
211  strcpy(cFileListName, argv[2]);
212 
213  if(mode==STANDARD_MODE || mode == STANDARD_DRY_MODE){
214  if (myrank == 0) StdFace_main(argv[2]);
215  strcpy(cFileListName, "namelist.def");
216  if (mode == STANDARD_DRY_MODE){
217  fprintf(stdout, "Dry run is Finished. \n\n");
218  return 0;
219  }
220  }
221 
222  setmem_HEAD(&X.Bind);
223  if(ReadDefFileNInt(cFileListName, &(X.Bind.Def), &(X.Bind.Boost))!=0){
224  fprintf(stdoutMPI, "%s", "Error: Definition files(*.def) are incomplete.\n");
225  exitMPI(-1);
226  }
227 
228  if (X.Bind.Def.nvec < X.Bind.Def.k_exct){
229  fprintf(stdoutMPI, "%s", "Error: nvec should be smaller than exct are incorrect.\n");
230  fprintf(stdoutMPI, "Error: nvec = %d, exct=%d.\n",
231  X.Bind.Def.nvec, X.Bind.Def.k_exct);
232  exitMPI(-1);
233  }
234  fprintf(stdoutMPI, "%s", "\n###### Definition files are correct. ######\n\n");
235 
236  /*ALLOCATE-------------------------------------------*/
237  setmem_def(&X.Bind, &X.Bind.Boost);
238  /*-----------------------------------------------------*/
239 
240  /*Read Def files.*/
241  TimeKeeper(&(X.Bind), "%s_TimeKeeper.dat", "Read File starts: %s", "w");
242  if(ReadDefFileIdxPara(&(X.Bind.Def), &(X.Bind.Boost))!=0){
243  fprintf(stdoutMPI,
244  "Error: Indices and Parameters of Definition files(*.def) are incomplete.\n");
245  exitMPI(-1);
246  }
247  TimeKeeper(&(X.Bind), "%s_TimeKeeper.dat", "Read File finishes: %s", "a");
248  fprintf(stdoutMPI, "%s", "\n###### Indices and Parameters of Definition files(*.def) are complete. ######\n\n");
249 
250  /*Set convergence Factor*/
251  SetConvergenceFactor(&(X.Bind.Def));
252 
253  /*---------------------------*/
254  if(HPhiTrans(&(X.Bind))!=0) {
255  exitMPI(-1);
256  }
257 
258  //Start Calculation
259  if(X.Bind.Def.iFlgCalcSpec == CALCSPEC_NOT ||
260  X.Bind.Def.iFlgCalcSpec == CALCSPEC_SCRATCH) {
261 
262  if(check(&(X.Bind))==MPIFALSE){
263  exitMPI(-1);
264  }
265 
266  /*LARGE VECTORS ARE ALLOCATED*/
267  if (setmem_large(&X.Bind) != 0) {
268  fprintf(stdoutMPI, "Error: Fail for memory allocation.");
269  exitMPI(-1);
270  }
271 
272  StartTimer(1000);
273  if(sz(&(X.Bind), list_1, list_2_1, list_2_2)!=0){
274  exitMPI(-1);
275  }
276 
277  StopTimer(1000);
278  if(X.Bind.Def.WRITE==1){
279  output_list(&(X.Bind));
280  exitMPI(-2);
281  }
282  StartTimer(2000);
283  diagonalcalc(&(X.Bind));
284  StopTimer(2000);
285 
286  switch (X.Bind.Def.iCalcType) {
287  case CG:
288  if (CalcByLOBPCG(&X) != TRUE) {
289  exitMPI(-3);
290  }
291  break;
292 
293  case FullDiag:
294  StartTimer(5000);
295  if (X.Bind.Def.iFlgScaLAPACK == 0 && nproc != 1) {
296  fprintf(stdoutMPI, "Error: Full Diagonalization by LAPACK is only allowed for one process.\n");
297  FinalizeMPI();
298  }
299  if (CalcByFullDiag(&X) != TRUE) {
300  FinalizeMPI();
301  }
302  StopTimer(5000);
303  break;
304 
305  case TPQCalc:
306  StartTimer(3000);
307  if (CalcByTPQ(NumAve, X.Bind.Def.Param.ExpecInterval, &X) != TRUE) {
308  StopTimer(3000);
309  exitMPI(-3);
310  }
311  StopTimer(3000);
312  break;
313 
314  case TimeEvolution:
315  if (CalcByTEM(X.Bind.Def.Param.ExpecInterval, &X) != 0) {
316  exitMPI(-3);
317  }
318  break;
319 
320  default:
321  StopTimer(0);
322  exitMPI(-3);
323  }
324  }
325 
326  if(X.Bind.Def.iFlgCalcSpec != CALCSPEC_NOT){
327  StartTimer(6000);
328  if (CalcSpectrum(&X) != TRUE) {
329  StopTimer(6000);
330  exitMPI(-3);
331  }
332  StopTimer(6000);
333  }
334 
335  StopTimer(0);
336  OutputTimer(&(X.Bind));
337  FinalizeMPI();
338  return 0;
339 }
void exitMPI(int errorcode)
MPI Abortation wrapper.
Definition: wrapperMPI.cpp:86
void SetConvergenceFactor(struct DefineList *X)
function to set convergence factors
Definition: readdef.cpp:2584
int nproc
Number of processors, defined in InitializeMPI()
Definition: global.cpp:72
int JudgeDefType(const int argc, char *argv[], int *mode)
function of judging a type of define files.
Definition: readdef.cpp:2453
int HPhiTrans(struct BindStruct *X)
Function of checking transfers not to count the same type of operators. .
Definition: HPhiTrans.cpp:45
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.cpp:75
long int * list_2_1
Definition: global.cpp:27
long int * list_2_2
Definition: global.cpp:28
int NumAve
Definition: global.cpp:43
void FinalizeMPI()
MPI Finitialization wrapper.
Definition: wrapperMPI.cpp:74
void splash()
Print logo mark and version number.
Definition: splash.cpp:25
int CalcByTEM(const int ExpecInterval, struct EDMainCalStruct *X)
main function of time evolution calculation
Definition: CalcByTEM.cpp:87
int output_list(struct BindStruct *X)
Output list_1 for canonical ensembles.
Definition: output_list.cpp:40
int CalcByLOBPCG(struct EDMainCalStruct *X)
Driver routine for LOB(P)CG method.
void setmem_HEAD(struct BindStruct *X)
Set size of memories headers of output files.
Definition: xsetmem.cpp:39
int diagonalcalc(struct BindStruct *X)
Calculate diagonal components and obtain the list, list_diagonal.
void StdFace_main(char *fname)
Main routine for the standard mode.
int sz(struct BindStruct *X, long int *list_1_, long int *list_2_1_, long int *list_2_2_)
generating Hilbert space
Definition: sz.cpp:908
void InitTimer()
function for initializing Timer[]
Definition: time.cpp:53
void InitializeMPI(int argc, char *argv[])
MPI initialization wrapper Process ID (myrank), Number of processes (nproc), Number of threads (nthre...
Definition: wrapperMPI.cpp:44
void setmem_def(struct BindStruct *X, struct BoostList *xBoost)
Set size of memories for Def and Phys in BindStruct.
Definition: xsetmem.cpp:53
int myrank
Process ID, defined in InitializeMPI()
Definition: global.cpp:73
void OutputTimer(struct BindStruct *X)
function for outputting elapse time for each function
Definition: time.cpp:95
int CalcByTPQ(const int NumAve, const int ExpecInterval, struct EDMainCalStruct *X)
A main function to calculate physical quqntities by TPQ method.
Definition: CalcByTPQ.cpp:48
void StopTimer(int n)
function for calculating elapse time [elapse time=StartTimer-StopTimer]
Definition: time.cpp:83
int CalcSpectrum(struct EDMainCalStruct *X)
A main function to calculate spectrum.
int check(struct BindStruct *X)
A program to check size of dimension for Hilbert-space.
Definition: check.cpp:51
int ReadDefFileIdxPara(struct DefineList *X, struct BoostList *xBoost)
function of reading def files to get keyword index
Definition: readdef.cpp:1399
int setmem_large(struct BindStruct *X)
Set size of memories for vectors(vg, v0, v1, v2, vec, alpha, beta), lists (list_1, list_2_1, list_2_2, list_Diagonal) and Phys(BindStruct.PhysList) struct in the case of Full Diag mode.
Definition: xsetmem.cpp:154
int TimeKeeper(struct BindStruct *X, const char *cFileName, const char *cTimeKeeper_Message, const char *cWriteType)
Functions for writing a time log.
Definition: log.cpp:42
int ReadDefFileNInt(char *xNameListFile, struct DefineList *X, struct BoostList *xBoost)
Function of reading information about "ModPara" file and total number of parameters from other def fi...
Definition: readdef.cpp:456
long int * list_1
Definition: global.cpp:25
void StartTimer(int n)
function for initializing elapse time [start]
Definition: time.cpp:71
int CalcByFullDiag(struct EDMainCalStruct *X)
Parent function for FullDiag mode.