HPhi++  3.1.0
output.cpp
Go to the documentation of this file.
1 /* HPhi - Quantum Lattice Model Simulator */
2 /* Copyright (C) 2015 The University of Tokyo */
3 
4 /* This program is free software: you can redistribute it and/or modify */
5 /* it under the terms of the GNU General Public License as published by */
6 /* the Free Software Foundation, either version 3 of the License, or */
7 /* (at your option) any later version. */
8 
9 /* This program is distributed in the hope that it will be useful, */
10 /* but WITHOUT ANY WARRANTY; without even the implied warranty of */
11 /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
12 /* GNU General Public License for more details. */
13 
14 /* You should have received a copy of the GNU General Public License */
15 /* along with this program. If not, see <http://www.gnu.org/licenses/>. */
16 #include "output.hpp"
17 #include "FileIO.hpp"
18 #include "wrapperMPI.hpp"
19 
20 
26 
27 int output(struct BindStruct *X) {
28 
29  FILE *fp;
30  char sdt[D_FileNameMax];
31  long int i, i_max;
32  i_max = X->Check.idim_max;
33 
34  if (X->Def.iCalcType == FullDiag) {
35  switch (X->Def.iCalcModel) {
36  case Spin:
37  case Hubbard:
38  case Kondo:
39  sprintf(sdt, "%s_phys_Nup%d_Ndown%d.dat", X->Def.CDataFileHead, X->Def.Nup, X->Def.Ndown);
40  break;
41  case SpinGC:
42  case HubbardGC:
43  case KondoGC:
44  sprintf(sdt, "%s_phys.dat", X->Def.CDataFileHead);
45  break;
46  default:
47  break;
48  }
49  if (childfopenMPI(sdt, "w", &fp) != 0) {
50  return -1;
51  }
52  fprintf(fp, " <H> <N> <Sz> <S2> <D> \n");
53  for (i = 0; i < i_max; i++) {
54  fprintf(fp, " %10lf %10lf %10lf %10lf %10lf\n", X->Phys.energy[i], X->Phys.num_up[i]+X->Phys.num_down[i], X->Phys.Sz[i],
55  X->Phys.s2[i], X->Phys.doublon[i]);
56  }
57  fclose(fp);
58  }
59  else{
60  fprintf(stdoutMPI, "Error: output function is used only for FullDiag mode.");
61  return -1;
62  }
63 
64  return 0;
65 }
66 
72 
73 int outputHam(struct BindStruct *X){
74  long int i=0;
75  long int j=0;
76  long int imax = X->Check.idim_max;
77  long int ihermite=0;
78  char cHeader[256];
79  FILE *fp;
80  char sdt[D_FileNameMax];
81 
82 #pragma omp parallel for default(none) reduction(+:ihermite) firstprivate(imax) private(i, j) shared(v0)
83  for (i=1; i<=imax; i++){
84  for (j=1; j<=i; j++){
85  if(abs(v0[i][j])>1.0e-13){
86  ihermite += 1;
87  }
88  }
89  }
90 
91  strcpy(cHeader, "%%%%MatrixMarket matrix coordinate complex hermitian\n");
92  sprintf(sdt,"%s_Ham.dat", X->Def.CDataFileHead);
93  if(childfopenMPI(sdt,"w",&fp)!=0){
94  return -1;
95  }
96  fprintf(fp, "%s", cHeader);
97  fprintf(fp, "%ld %ld %ld \n", imax, imax, ihermite);
98  for (i=1; i<=imax; i++){
99  for (j=1; j<=i; j++){
100  if(abs(v0[i][j])>1.0e-13){
101  fprintf(fp, "%ld %ld %lf %lf\n",i,j,real(v0[i][j]),imag(v0[i][j]));
102  }
103  }
104  }
105  fclose(fp);
106  return 0;
107 }
double * doublon
Expectation value of the Doublon.
Definition: struct.hpp:357
int Nup
Number of spin-up electrons in this process.
Definition: struct.hpp:58
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.hpp:395
double * num_down
Expectation value of the number of down-spin electtrons.
Definition: struct.hpp:364
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.cpp:75
std::complex< double > ** v0
Definition: global.cpp:20
double * s2
Expectation value of the square of the total S.
Definition: struct.hpp:365
struct PhysList Phys
Physical quantities.
Definition: struct.hpp:398
double * Sz
Expectation value of the Total Sz.
Definition: struct.hpp:361
int output(struct BindStruct *X)
output function for FullDiag mode
Definition: output.cpp:27
Bind.
Definition: struct.hpp:394
double * energy
Expectation value of the total energy.
Definition: struct.hpp:356
double * num_up
Expectation value of the number of up-spin electtrons.
Definition: struct.hpp:363
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.hpp:200
int outputHam(struct BindStruct *X)
output Hamiltonian only used for FullDiag mode
Definition: output.cpp:73
int Ndown
Number of spin-down electrons in this process.
Definition: struct.hpp:59
struct CheckList Check
Size of the Hilbert space.
Definition: struct.hpp:396
char * CDataFileHead
Read from Calcmod in readdef.h. Header of output file such as Green&#39;s function.
Definition: struct.hpp:42
long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.hpp:305
int iCalcType
Switch for calculation type. 0:Lanczos, 1:TPQCalc, 2:FullDiag.
Definition: struct.hpp:194
int childfopenMPI(const char *_cPathChild, const char *_cmode, FILE **_fp)
Only the root process open file in output/ directory.
Definition: FileIO.cpp:27