HPhi++  3.1.0
CalcByTPQ.cpp File Reference

File for givinvg functions of TPQ method. More...

#include "FirstMultiply.hpp"
#include "Multiply.hpp"
#include "expec_energy_flct.hpp"
#include "expec_cisajs.hpp"
#include "expec_cisajscktaltdc.hpp"
#include "CalcByTPQ.hpp"
#include "FileIO.hpp"
#include "wrapperMPI.hpp"
#include "CalcTime.hpp"
#include "common/setmemory.hpp"
#include "mltplyCommon.hpp"

Go to the source code of this file.

Functions

int CalcByTPQ (const int NumAve, const int ExpecInterval, struct EDMainCalStruct *X)
 A main function to calculate physical quqntities by TPQ method. More...
 

Detailed Description

File for givinvg functions of TPQ method.

Version
0.1, 0.2
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)

Definition in file CalcByTPQ.cpp.

Function Documentation

◆ CalcByTPQ()

int CalcByTPQ ( const int  NumAve,
const int  ExpecInterval,
struct EDMainCalStruct X 
)

A main function to calculate physical quqntities by TPQ method.

Parameters
[in]NumAveNumber of samples
[in]ExpecIntervalinterval steps between the steps to calculate physical quantities
[in,out]XCalcStruct list for getting and giving calculation information
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)
Return values
0normally finished
-1unnormally finished

Initialize v1 and compute v0 = H*v1

Compute v1=0, and compute v0 = H*v1

Definition at line 48 of file CalcByTPQ.cpp.

References EDMainCalStruct::Bind, BindStruct::Check, childfopenALL(), childfopenMPI(), BindStruct::Def, PhysList::doublon, PhysList::doublon2, PhysList::energy, exitMPI(), expec_cisajs(), expec_cisajscktaltdc(), expec_energy_flct(), FirstMultiply(), global_1st_norm, global_norm, CheckList::idim_max, DefineList::iReStart, DefineList::istep, DefineList::Lanczos_max, LargeValue, Multiply(), myrank, DefineList::NsiteMPI, PhysList::num, PhysList::num2, NumAve, BindStruct::Phys, DefineList::St, StartTimer(), stdoutMPI, step_i, step_spin, StopTimer(), PhysList::Sz, PhysList::Sz2, TimeKeepStruct::tend, TimeKeeperWithRandAndStep(), TimeKeepStruct::tstart, v0, v1, and PhysList::var.

Referenced by main().

52  {
53  char sdt[D_FileNameMax];
54  char **sdt_phys, **sdt_norm, **sdt_flct;
55  int rand_i, iret;
56  long int i_max;
57  int step_iO = 0;
58  FILE *fp;
59  double *inv_temp, Ns;
60  struct TimeKeepStruct tstruct;
61  size_t byte_size;
62 
63  tstruct.tstart = time(NULL);
64  inv_temp = d_1d_allocate(NumAve);
65 
66  step_spin = ExpecInterval;
67  X->Bind.Def.St = 0;
68  fprintf(stdoutMPI, "%s", "###### Start: TPQCalculation. ######\n\n");
69  global_norm = d_1d_allocate(NumAve);
70  global_1st_norm = d_1d_allocate(NumAve);
71 
72  //for rand_i =0, rand_i<NumAve
73  sdt_phys = (char**)malloc(sizeof(char*)*NumAve);
74  sdt_norm = (char**)malloc(sizeof(char*)*NumAve);
75  sdt_flct = (char**)malloc(sizeof(char*)*NumAve);
76  for (rand_i = 0; rand_i < NumAve; rand_i++) {
77  sdt_phys[rand_i] = (char*)malloc(sizeof(char)*D_FileNameMax);
78  sdt_norm[rand_i] = (char*)malloc(sizeof(char)*D_FileNameMax);
79  sdt_flct[rand_i] = (char*)malloc(sizeof(char)*D_FileNameMax);
80  sprintf(sdt_phys[rand_i], "SS_rand%d.dat", rand_i);
81  sprintf(sdt_norm[rand_i], "Norm_rand%d.dat", rand_i);
82  sprintf(sdt_flct[rand_i], "Flct_rand%d.dat", rand_i);
83  }
84  Ns = 1.0 * X->Bind.Def.NsiteMPI;
85  fprintf(stdoutMPI, " rand_i / rand_max = %d / %d\n", 1, NumAve);
86  iret = 0;
87 
88  //Make or Read initial vector
89  if (X->Bind.Def.iReStart == RESTART_INOUT || X->Bind.Def.iReStart == RESTART_IN) {
90  StartTimer(3600);
91  TimeKeeperWithRandAndStep(&(X->Bind), "%s_Time_TPQ_Step.dat", " set %d step %d:output vector starts: %s\n", "a", 0, step_i);
92  fprintf(stdoutMPI, "%s", " Start: Input vector.\n");
93  sprintf(sdt, "tmpvec_set%d_rank_%d.dat", rand_i, myrank);
94  childfopenALL(sdt, "rb", &fp);
95  if (fp == NULL) {
96  fprintf(stdout, "A file of Inputvector does not exist.\n");
97  fprintf(stdout, "Start to calculate in normal procedure.\n");
98  iret = 1;
99  }
100  byte_size = fread(&step_i, sizeof(step_i), 1, fp);
101  byte_size = fread(&i_max, sizeof(long int), 1, fp);
102  if (i_max != X->Bind.Check.idim_max) {
103  fprintf(stderr, "Error: A file of Inputvector is incorrect.\n");
104  exitMPI(-1);
105  }
106  byte_size = fread(v0, sizeof(std::complex<double>), (X->Bind.Check.idim_max + 1)*NumAve, fp);
107  TimeKeeperWithRandAndStep(&(X->Bind), "%s_Time_TPQ_Step.dat", " set %d step %d:output vector finishes: %s\n", "a", 0, step_i);
108  fprintf(stdoutMPI, "%s", " End : Input vector.\n");
109  fclose(fp);
110  StopTimer(3600);
111  X->Bind.Def.istep = step_i;
112  StartTimer(3200);
113  iret = expec_energy_flct(&(X->Bind), NumAve, v0, v1);
114  StopTimer(3200);
115  if (iret != 0) return -1;
116 
117  step_iO = step_i - 1;
118  if (byte_size == 0) printf("byte_size: %d \n", (int)byte_size);
119  }/*if (X->Bind.Def.iReStart == RESTART_INOUT || X->Bind.Def.iReStart == RESTART_IN)*/
120 
121  if (X->Bind.Def.iReStart == RESTART_NOT || X->Bind.Def.iReStart == RESTART_OUT || iret == 1) {
122  StartTimer(3600);
123  for (rand_i = 0; rand_i < NumAve; rand_i++) {
124  if (childfopenMPI(sdt_phys[rand_i], "w", &fp) == 0) {
125  fprintf(fp, "%s", " # inv_tmp, energy, phys_var, phys_doublon, phys_num, step_i\n");
126  fclose(fp);
127  }
128  else return -1;
129  // for norm
130  if (childfopenMPI(sdt_norm[rand_i], "w", &fp) == 0) {
131  fprintf(fp, "%s", " # inv_temp, global_norm, global_1st_norm, step_i \n");
132  fclose(fp);
133  }
134  else return -1;
135  // for fluctuations
136  if (childfopenMPI(sdt_flct[rand_i], "w", &fp) == 0) {
137  fprintf(fp, "%s", " # inv_temp, N, N^2, D, D^2, Sz, Sz^2, step_i \n");
138  fclose(fp);
139  }
140  else return -1;
141  }
142  StopTimer(3600);
143 
144  step_i = 0;
145 
146  StartTimer(3100);
147  if (rand_i == 0) {
148  TimeKeeperWithRandAndStep(&(X->Bind), "%s_Time_TPQ_Step.dat", "set %d step %d:TPQ begins: %s", "w", 0, step_i);
149  }
150  else {
151  TimeKeeperWithRandAndStep(&(X->Bind), "%s_Time_TPQ_Step.dat", "set %d step %d:TPQ begins: %s", "a", 0, step_i);
152  }
156  FirstMultiply(&(X->Bind));
157  StopTimer(3100);
158  for (rand_i = 0; rand_i < NumAve; rand_i++) {
159  inv_temp[rand_i] = 0.0;
160  if (childfopenMPI(sdt_phys[rand_i], "a", &fp) == 0) {
161  fprintf(fp, "%.16lf %.16lf %.16lf %.16lf %.16lf %d\n",
162  inv_temp[rand_i], X->Bind.Phys.energy[rand_i], X->Bind.Phys.var[rand_i],
163  X->Bind.Phys.doublon[rand_i], X->Bind.Phys.num[rand_i], step_i);
164  fclose(fp);
165  }
166  else return -1;
167  // for norm
168  if (childfopenMPI(sdt_norm[rand_i], "a", &fp) == 0) {
169  fprintf(fp, "%.16lf %.16lf %.16lf %d\n",
170  inv_temp[rand_i], global_1st_norm[rand_i], global_1st_norm[rand_i], step_i);
171  fclose(fp);
172  }
173  else return -1;
174  }
178  StartTimer(3200);
179  iret = expec_energy_flct(&(X->Bind), NumAve, v0, v1); //v0 = H*v1
180  StopTimer(3200);
181  if (iret != 0) return -1;
182  step_i += 1;
183  StartTimer(3600);
184  for (rand_i = 0; rand_i < NumAve; rand_i++) {
185  inv_temp[rand_i] = (2.0 / Ns) / (LargeValue - X->Bind.Phys.energy[rand_i] / Ns);
186  if (childfopenMPI(sdt_phys[rand_i], "a", &fp) == 0) {
187  fprintf(fp, "%.16lf %.16lf %.16lf %.16lf %.16lf %d\n",
188  inv_temp[rand_i], X->Bind.Phys.energy[rand_i], X->Bind.Phys.var[rand_i],
189  X->Bind.Phys.doublon[rand_i], X->Bind.Phys.num[rand_i], step_i);
190  fclose(fp);
191  }
192  else return -1;
193  // for norm
194  if (childfopenMPI(sdt_norm[rand_i], "a", &fp) == 0) {
195  fprintf(fp, "%.16lf %.16lf %.16lf %d\n",
196  inv_temp[rand_i], global_norm[rand_i], global_1st_norm[rand_i], step_i);
197  fclose(fp);
198  }
199  else return -1;
200  // for fluctuations
201  if (childfopenMPI(sdt_flct[rand_i], "a", &fp) == 0) {
202  fprintf(fp, "%.16lf %.16lf %.16lf %.16lf %.16lf %.16lf %.16lf %d\n",
203  inv_temp[rand_i], X->Bind.Phys.num[rand_i], X->Bind.Phys.num2[rand_i],
204  X->Bind.Phys.doublon[rand_i], X->Bind.Phys.doublon2[rand_i],
205  X->Bind.Phys.Sz[rand_i], X->Bind.Phys.Sz2[rand_i], step_i);
206  fclose(fp);
207  }
208  else return -1;
209  }/*for (rand_i = 0; rand_i < NumAve; rand_i++)*/
210  StopTimer(3600);
211  step_i += 1;
212  X->Bind.Def.istep = step_i;
213  step_iO = 0;
214  }/*if (X->Bind.Def.iReStart == RESTART_NOT || X->Bind.Def.iReStart == RESTART_OUT || iret == 1)*/
215 
216  for (step_i = X->Bind.Def.istep; step_i < X->Bind.Def.Lanczos_max; step_i++) {
217  X->Bind.Def.istep = step_i;
218  if (step_i % ((X->Bind.Def.Lanczos_max - step_iO) / 10) == 0) {
219  fprintf(stdoutMPI, " step_i/total_step = %d/%d \n", step_i, X->Bind.Def.Lanczos_max);
220  }
221  X->Bind.Def.istep = step_i;
222  StartTimer(3600);
223  TimeKeeperWithRandAndStep(&(X->Bind), "%s_Time_TPQ_Step.dat", "set %d step %d:TPQ begins: %s", "a", 0, step_i);
224  StopTimer(3600);
225  StartTimer(3500);
226  Multiply(&(X->Bind));
227  StopTimer(3500);
228 
229  if (step_i%step_spin == 0) {
230  StartTimer(3300);
231  iret = expec_cisajs(&(X->Bind), NumAve, v1, v0);
232  StopTimer(3300);
233  if (iret != 0) return -1;
234 
235  StartTimer(3400);
236  iret = expec_cisajscktaltdc(&(X->Bind), NumAve, v1, v0);
237  StopTimer(3400);
238  if (iret != 0) return -1;
239  }
240 
241  StartTimer(3200);
242  iret = expec_energy_flct(&(X->Bind), NumAve, v0, v1);
243  StopTimer(3200);
244  if (iret != 0) return -1;
245 
246  StartTimer(3600);
247  for (rand_i = 0; rand_i < NumAve; rand_i++) {
248  inv_temp[rand_i] = (2.0*step_i / Ns) / (LargeValue - X->Bind.Phys.energy[rand_i] / Ns);
249  if (childfopenMPI(sdt_phys[rand_i], "a", &fp) == 0) {
250  fprintf(fp, "%.16lf %.16lf %.16lf %.16lf %.16lf %d\n",
251  inv_temp[rand_i], X->Bind.Phys.energy[rand_i], X->Bind.Phys.var[rand_i],
252  X->Bind.Phys.doublon[rand_i], X->Bind.Phys.num[rand_i], step_i);
253  // for
254  fclose(fp);
255  }
256  else return FALSE;
257 
258  if (childfopenMPI(sdt_norm[rand_i], "a", &fp) == 0) {
259  fprintf(fp, "%.16lf %.16lf %.16lf %d\n",
260  inv_temp[rand_i], global_norm[rand_i], global_1st_norm[rand_i], step_i);
261  fclose(fp);
262  }
263  else return FALSE;
264 
265  // for fluctuations
266  if (childfopenMPI(sdt_flct[rand_i], "a", &fp) == 0) {
267  fprintf(fp, "%.16lf %.16lf %.16lf %.16lf %.16lf %.16lf %.16lf %d\n",
268  inv_temp[rand_i], X->Bind.Phys.num[rand_i], X->Bind.Phys.num2[rand_i],
269  X->Bind.Phys.doublon[rand_i], X->Bind.Phys.doublon2[rand_i],
270  X->Bind.Phys.Sz[rand_i], X->Bind.Phys.Sz2[rand_i], step_i);
271  fclose(fp);
272  }
273  else return -1;
274  }/*for (rand_i = 0; rand_i < NumAve; rand_i++)*/
275  StopTimer(3600);
276  }/*for (step_i = X->Bind.Def.istep; step_i < X->Bind.Def.Lanczos_max; step_i++)*/
277 
278  if (X->Bind.Def.iReStart == RESTART_OUT || X->Bind.Def.iReStart == RESTART_INOUT) {
279  TimeKeeperWithRandAndStep(&(X->Bind), "%s_Time_TPQ_Step.dat", " set %d step %d:output vector starts: %s\n", "a", 0, step_i);
280  fprintf(stdoutMPI, "%s", " Start: Output vector.\n");
281  sprintf(sdt, "tmpvec_set%d_rank_%d.dat", 0, myrank);
282  if (childfopenALL(sdt, "wb", &fp) != 0) {
283  exitMPI(-1);
284  }
285  fwrite(&step_i, sizeof(step_i), 1, fp);
286  fwrite(&X->Bind.Check.idim_max, sizeof(X->Bind.Check.idim_max), 1, fp);
287  fwrite(v1, sizeof(std::complex<double>), (X->Bind.Check.idim_max + 1)*NumAve, fp);
288  fclose(fp);
289  TimeKeeperWithRandAndStep(&(X->Bind), "%s_Time_TPQ_Step.dat", " set %d step %d:output vector finishes: %s\n", "a", 0, step_i);
290  fprintf(stdoutMPI, "%s", " End : Output vector.\n");
291  }
292 
293  fprintf(stdoutMPI, "%s", "###### End : TPQCalculation. ######\n\n");
294 
295  tstruct.tend = time(NULL);
296  fprintf(stdoutMPI, "Finish: Elapsed time is %d [s].\n", (int)(tstruct.tend - tstruct.tstart));
297  free_d_1d_allocate(inv_temp);
298 
299  for (rand_i = 0; rand_i < NumAve; rand_i++) {
300  free(sdt_phys[rand_i]);
301  free(sdt_norm[rand_i]);
302  free(sdt_flct[rand_i]);
303  }
304  free(sdt_phys);
305  free(sdt_norm);
306  free(sdt_flct);
307 
308  return TRUE;
309 }
double * doublon
Expectation value of the Doublon.
Definition: struct.hpp:357
void exitMPI(int errorcode)
MPI Abortation wrapper.
Definition: wrapperMPI.cpp:86
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.hpp:395
double * var
Expectation value of the Energy variance.
Definition: struct.hpp:367
int St
0 or 1, but it affects nothing.
Definition: struct.hpp:80
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.cpp:75
int expec_cisajscktaltdc(struct BindStruct *X, int nstate, std::complex< double > **Xvec, std::complex< double > **vec)
Parent function to calculate two-body green&#39;s functions.
int childfopenALL(const char *_cPathChild, const char *_cmode, FILE **_fp)
All processes open file in output/ directory.
Definition: FileIO.cpp:50
Time keeping.
Definition: struct.hpp:410
double * global_1st_norm
Definition: global.cpp:46
int expec_energy_flct(struct BindStruct *X, int nstate, std::complex< double > **tmp_v0, std::complex< double > **tmp_v1)
Parent function to calculate expected values of energy and physical quantities.
int iReStart
Definition: struct.hpp:223
int NumAve
Definition: global.cpp:43
int Multiply(struct BindStruct *X)
Function of calculating the i-th step norm as and update the i+1-th wave vector as for TPQ calculat...
Definition: Multiply.cpp:44
std::complex< double > ** v0
Definition: global.cpp:20
int step_spin
Definition: global.cpp:47
double * num2
Expectation value of the quare of the number of electrons.
Definition: struct.hpp:360
struct PhysList Phys
Physical quantities.
Definition: struct.hpp:398
double * Sz2
Expectation value of the Square of total Sz.
Definition: struct.hpp:362
std::complex< double > ** v1
Definition: global.cpp:21
double * Sz
Expectation value of the Total Sz.
Definition: struct.hpp:361
int NsiteMPI
Total number of sites, differ from DefineList::Nsite.
Definition: struct.hpp:57
int TimeKeeperWithRandAndStep(struct BindStruct *X, const char *cFileName, const char *cTimeKeeper_Message, const char *cWriteType, const int irand, const int istep)
Functions for writing a time log.
Definition: log.cpp:117
double * num
Expectation value of the Number of electrons.
Definition: struct.hpp:359
int Lanczos_max
Maximum number of iterations.
Definition: struct.hpp:74
double * energy
Expectation value of the total energy.
Definition: struct.hpp:356
double * global_norm
Definition: global.cpp:45
int FirstMultiply(struct BindStruct *X)
Multiplication at the first step for TPQ mode ( is the random or inputted vector).
int myrank
Process ID, defined in InitializeMPI()
Definition: global.cpp:73
int expec_cisajs(struct BindStruct *X, int nstate, std::complex< double > **Xvec, std::complex< double > **vec)
function of calculation for one body green&#39;s function
int step_i
Definition: global.cpp:44
double LargeValue
Definition: global.cpp:42
void StopTimer(int n)
function for calculating elapse time [elapse time=StartTimer-StopTimer]
Definition: time.cpp:83
int istep
Index of TPQ step ???
Definition: struct.hpp:78
double * doublon2
Expectation value of the Square of doublon.
Definition: struct.hpp:358
time_t tstart
Start time.
Definition: struct.hpp:411
struct CheckList Check
Size of the Hilbert space.
Definition: struct.hpp:396
long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.hpp:305
struct BindStruct Bind
Binded struct.
Definition: struct.hpp:405
void StartTimer(int n)
function for initializing elapse time [start]
Definition: time.cpp:71
int childfopenMPI(const char *_cPathChild, const char *_cmode, FILE **_fp)
Only the root process open file in output/ directory.
Definition: FileIO.cpp:27