HΦ  3.2.0
Functions
diagonalcalc.c File Reference

Calculate diagonal components, i.e. \( H_d |\phi_0> = E_d |\phi_0> \). More...

#include <bitcalc.h>
#include "FileIO.h"
#include "diagonalcalc.h"
#include "mltplySpinCore.h"
#include "wrapperMPI.h"
+ Include dependency graph for diagonalcalc.c:

Go to the source code of this file.

Functions

int SetDiagonalTETransfer (long unsigned int isite1, double dtmp_V, long unsigned int spin, struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1)
 Update the vector by the general one-body diagonal interaction, \( \mu_{i\sigma_1} n_ {i\sigma_1}\).
(Using in Time Evolution mode). More...
 
int SetDiagonalTEInterAll (long unsigned int isite1, long unsigned int isite2, long unsigned int isigma1, long unsigned int isigma2, double dtmp_V, struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1)
 Update the vector by the general two-body diagonal interaction, \( H_{i\sigma_1 j\sigma_2} n_ {i\sigma_1}n_{j\sigma_2}\).
(Using in Time Evolution mode). More...
 
int SetDiagonalTEChemi (long unsigned int isite1, long unsigned int spin, double dtmp_V, struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1)
 Update the vector by the chemical potential \( \mu_{i \sigma_1} n_ {i \sigma_1} \)
generated by the commutation relation in terms of the general two-body interaction,
\( c_ {i \sigma_1} a_{j\sigma_2}c_ {j \sigma_2}a_ {i \sigma_1} = c_ {i \sigma_1}a_ {i \sigma_1}-c_ {i \sigma_1} a_ {i \sigma_1} c_ {j \sigma_2}a_{j\sigma_2}\) . (Using in Time Evolution mode). More...
 
int diagonalcalc (struct BindStruct *X)
 Calculate diagonal components and obtain the list, list_diagonal. More...
 
int diagonalcalcForTE (const int _istep, struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1)
 
int SetDiagonalCoulombIntra (long unsigned int isite1, double dtmp_V, struct BindStruct *X)
 Calculate the components for Coulombintra interaction, \( U_i n_ {i \uparrow}n_{i \downarrow} \). More...
 
int SetDiagonalChemi (long unsigned int isite1, double dtmp_V, long unsigned int spin, struct BindStruct *X)
 Calculate the components for the chemical potential \( \mu_{i \sigma_1} n_ {i \sigma_1} \). More...
 
int SetDiagonalCoulombInter (long unsigned int isite1, long unsigned int isite2, double dtmp_V, struct BindStruct *X)
 Calculate the components for Coulombinter interaction, \( V_{ij} n_ {i}n_{j} \). More...
 
int SetDiagonalHund (long unsigned int isite1, long unsigned int isite2, double dtmp_V, struct BindStruct *X)
 Calculate the components for Hund interaction, \( H_{ij}(n_ {i\uparrow}n_{j\uparrow}+ n_ {i\downarrow}n_{j\downarrow})\). More...
 
int SetDiagonalInterAll (long unsigned int isite1, long unsigned int isite2, long unsigned int isigma1, long unsigned int isigma2, double dtmp_V, struct BindStruct *X)
 Calculate the components for general two-body diagonal interaction, \( H_{i\sigma_1 j\sigma_2} n_ {i\sigma_1}n_{j\sigma_2}\). More...
 

Detailed Description

Calculate diagonal components, i.e. \( H_d |\phi_0> = E_d |\phi_0> \).

Version
2.1

add functions to calculate diagonal components for Time evolution.

Author
Kazuyoshi Yoshimi (The University of Tokyo)
Version
0.2

modify functions to calculate diagonal components for general spin.

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

Definition in file diagonalcalc.c.

Function Documentation

◆ diagonalcalc()

int diagonalcalc ( struct BindStruct *  X)

Calculate diagonal components and obtain the list, list_diagonal.

Parameters
X[in] Struct to get the information of the diagonal operators.
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)
Return values
-1fail to calculate diagonal components.
0succeed to calculate diagonal components.

Definition at line 85 of file diagonalcalc.c.

References cDiagonalCalcFinish, cDiagonalCalcStart, cFileNameCheckChemi, cFileNameCheckCoulombIntra, cFileNameCheckHund, cFileNameCheckInterAll, cFileNameCheckInterU, cFileNameTimeKeep, BindStruct::Check, childfopenMPI(), DefineList::CoulombInter, DefineList::CoulombIntra, cProEndCalcDiag, cProStartCalcDiag, BindStruct::Def, diagonalcalcForTE(), DefineList::EDChemi, DefineList::EDNChemi, DefineList::EDParaChemi, DefineList::EDSpinChemi, DefineList::HundCoupling, CheckList::idim_max, DefineList::InterAll_Diagonal, list_Diagonal, DefineList::NCoulombInter, DefineList::NCoulombIntra, DefineList::NHundCoupling, DefineList::NInterAll_Diagonal, DefineList::ParaCoulombInter, DefineList::ParaCoulombIntra, DefineList::ParaHundCoupling, DefineList::ParaInterAll_Diagonal, SetDiagonalChemi(), SetDiagonalCoulombInter(), SetDiagonalCoulombIntra(), SetDiagonalHund(), SetDiagonalInterAll(), stdoutMPI, and TimeKeeper().

Referenced by CalcSpectrum(), and main().

87  {
88 
89  FILE *fp;
90  long unsigned int i,j;
91  long unsigned int isite1,isite2;
92  long unsigned int spin;
93  double tmp_V;
94 
95  /*[s] For InterAll*/
96  long unsigned int A_spin,B_spin;
97  /*[e] For InterAll*/
98  long unsigned int i_max=X->Check.idim_max;
99 
100  fprintf(stdoutMPI, "%s", cProStartCalcDiag);
101  TimeKeeper(X, cFileNameTimeKeep, cDiagonalCalcStart, "a");
102 
103 #pragma omp parallel for default(none) private(j) shared(list_Diagonal) firstprivate(i_max)
104  for(j = 1;j <= i_max; j++){
105  list_Diagonal[j]=0.0;
106  }
107 
108  if(X->Def.NCoulombIntra>0){
109  if(childfopenMPI(cFileNameCheckCoulombIntra, "w", &fp)!=0){
110  return -1;
111  }
112  for(i = 0; i < X->Def.NCoulombIntra; i++){
113  isite1 = X->Def.CoulombIntra[i][0]+1;
114  tmp_V = X->Def.ParaCoulombIntra[i];
115  fprintf(fp,"i=%ld isite1=%ld tmp_V=%lf \n",i,isite1,tmp_V);
116  SetDiagonalCoulombIntra(isite1, tmp_V, X);
117  }
118  fclose(fp);
119  }
120 
121  if(X->Def.EDNChemi>0){
122  if(childfopenMPI(cFileNameCheckChemi,"w", &fp)!=0){
123  return -1;
124  }
125  for(i = 0; i < X->Def.EDNChemi; i++){
126  isite1 = X->Def.EDChemi[i]+1;
127  spin = X->Def.EDSpinChemi[i];
128  tmp_V = -X->Def.EDParaChemi[i];
129  fprintf(fp,"i=%ld spin=%ld isite1=%ld tmp_V=%lf \n",i,spin,isite1,tmp_V);
130  if(SetDiagonalChemi(isite1, tmp_V,spin, X) !=0){
131  return -1;
132  }
133  }
134  fclose(fp);
135  }
136 
137  if(X->Def.NCoulombInter>0){
138  if(childfopenMPI(cFileNameCheckInterU,"w", &fp)!=0){
139  return -1;
140  }
141  for(i = 0; i < X->Def.NCoulombInter; i++){
142  isite1 = X->Def.CoulombInter[i][0]+1;
143  isite2 = X->Def.CoulombInter[i][1]+1;
144  tmp_V = X->Def.ParaCoulombInter[i];
145  fprintf(fp,"i=%ld isite1=%ld isite2=%ld tmp_V=%lf \n",i,isite1,isite2,tmp_V);
146  if(SetDiagonalCoulombInter(isite1, isite2, tmp_V, X) !=0){
147  return -1;
148  }
149  }
150  fclose(fp);
151  }
152  if(X->Def.NHundCoupling>0){
153  if(childfopenMPI(cFileNameCheckHund,"w", &fp) !=0){
154  return -1;
155  }
156  for(i = 0; i < X->Def.NHundCoupling; i++){
157  isite1 = X->Def.HundCoupling[i][0]+1;
158  isite2 = X->Def.HundCoupling[i][1]+1;
159  tmp_V = -X->Def.ParaHundCoupling[i];
160  if(SetDiagonalHund(isite1, isite2, tmp_V, X) !=0){
161  return -1;
162  }
163  fprintf(fp,"i=%ld isite1=%ld isite2=%ld tmp_V=%lf \n",i,isite1,isite2,tmp_V);
164  }
165  fclose(fp);
166  }
167 
168  if(X->Def.NInterAll_Diagonal>0){
169  if(childfopenMPI(cFileNameCheckInterAll,"w", &fp) !=0){
170  return -1;
171  }
172  for(i = 0; i < X->Def.NInterAll_Diagonal; i++){
173  isite1=X->Def.InterAll_Diagonal[i][0]+1;
174  A_spin=X->Def.InterAll_Diagonal[i][1];
175  isite2=X->Def.InterAll_Diagonal[i][2]+1;
176  B_spin=X->Def.InterAll_Diagonal[i][3];
177  tmp_V = X->Def.ParaInterAll_Diagonal[i];
178  fprintf(fp,"i=%ld isite1=%ld A_spin=%ld isite2=%ld B_spin=%ld tmp_V=%lf \n", i, isite1, A_spin, isite2, B_spin, tmp_V);
179  SetDiagonalInterAll(isite1, isite2, A_spin, B_spin, tmp_V, X);
180  }
181  fclose(fp);
182  }
183 
184  TimeKeeper(X, cFileNameTimeKeep, cDiagonalCalcFinish, "a");
185  fprintf(stdoutMPI, "%s", cProEndCalcDiag);
186  return 0;
187 }
cFileNameCheckChemi
const char * cFileNameCheckChemi
Definition: global.c:28
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
DefineList::NHundCoupling
unsigned int NHundCoupling
Number of Hund coupling.
Definition: struct.h:133
DefineList::EDSpinChemi
int * EDSpinChemi
[DefineList::Nsite]
Definition: struct.h:99
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
cDiagonalCalcStart
const char * cDiagonalCalcStart
Definition: LogMessage.c:56
DefineList::ParaInterAll_Diagonal
double * ParaInterAll_Diagonal
[DefineList::NInterAll_Diagonal] Coupling constant of diagonal inter-all term. malloc in setmem_def()...
Definition: struct.h:168
cFileNameCheckHund
const char * cFileNameCheckHund
Definition: global.c:30
childfopenMPI
int childfopenMPI(const char *_cPathChild, const char *_cmode, FILE **_fp)
Only the root process open file in output/ directory.
Definition: FileIO.c:27
SetDiagonalInterAll
int SetDiagonalInterAll(long unsigned int isite1, long unsigned int isite2, long unsigned int isigma1, long unsigned int isigma2, double dtmp_V, struct BindStruct *X)
Calculate the components for general two-body diagonal interaction, .
Definition: diagonalcalc.c:1110
cProStartCalcDiag
const char * cProStartCalcDiag
Definition: ProgressMessage.c:33
cFileNameCheckCoulombIntra
const char * cFileNameCheckCoulombIntra
Definition: global.c:27
DefineList::InterAll_Diagonal
int ** InterAll_Diagonal
[DefineList::NinterAll_Diagonal][4] Interacted quartet
Definition: struct.h:162
DefineList::HundCoupling
int ** HundCoupling
[DefineList::NHundCoupling][2] Index of Hund coupling. malloc in setmem_def().
Definition: struct.h:134
DefineList::NInterAll_Diagonal
unsigned int NInterAll_Diagonal
Number of interall term (diagonal)
Definition: struct.h:164
cProEndCalcDiag
const char * cProEndCalcDiag
Definition: ProgressMessage.c:34
DefineList::ParaCoulombInter
double * ParaCoulombInter
[DefineList::NCoulombInter]Coupling constant of off-site Coulomb interaction. malloc in setmem_def()...
Definition: struct.h:130
DefineList::CoulombInter
int ** CoulombInter
Definition: struct.h:128
DefineList::NCoulombIntra
unsigned int NCoulombIntra
Definition: struct.h:121
SetDiagonalHund
int SetDiagonalHund(long unsigned int isite1, long unsigned int isite2, double dtmp_V, struct BindStruct *X)
Calculate the components for Hund interaction, .
Definition: diagonalcalc.c:785
DefineList::EDParaChemi
double * EDParaChemi
[DefineList::Nsite] On-site potential parameter. malloc in setmem_def().
Definition: struct.h:100
cFileNameCheckInterAll
const char * cFileNameCheckInterAll
Definition: global.c:31
SetDiagonalCoulombIntra
int SetDiagonalCoulombIntra(long unsigned int isite1, double dtmp_V, struct BindStruct *X)
Calculate the components for Coulombintra interaction, .
Definition: diagonalcalc.c:261
cFileNameCheckInterU
const char * cFileNameCheckInterU
Definition: global.c:29
cFileNameTimeKeep
const char * cFileNameTimeKeep
Definition: global.c:23
SetDiagonalCoulombInter
int SetDiagonalCoulombInter(long unsigned int isite1, long unsigned int isite2, double dtmp_V, struct BindStruct *X)
Calculate the components for Coulombinter interaction, .
Definition: diagonalcalc.c:544
DefineList::CoulombIntra
int ** CoulombIntra
Definition: struct.h:122
cDiagonalCalcFinish
const char * cDiagonalCalcFinish
Definition: LogMessage.c:57
DefineList::ParaHundCoupling
double * ParaHundCoupling
[DefineList::NHundCoupling] Hund coupling constant. malloc in setmem_def().
Definition: struct.h:136
list_Diagonal
double * list_Diagonal
Definition: global.h:46
SetDiagonalChemi
int SetDiagonalChemi(long unsigned int isite1, double dtmp_V, long unsigned int spin, struct BindStruct *X)
Calculate the components for the chemical potential .
Definition: diagonalcalc.c:374
DefineList::NCoulombInter
unsigned int NCoulombInter
Number of off-site Coulomb interaction.
Definition: struct.h:127
DefineList::ParaCoulombIntra
double * ParaCoulombIntra
Definition: struct.h:124
DefineList::EDNChemi
unsigned int EDNChemi
Number of on-site term.
Definition: struct.h:97
DefineList::EDChemi
int * EDChemi
[DefineList::Nsite] Chemical potential. malloc in setmem_def().
Definition: struct.h:98
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
TimeKeeper
int TimeKeeper(struct BindStruct *X, const char *cFileName, const char *cTimeKeeper_Message, const char *cWriteType)
Functions for writing a time log.
Definition: log.c:42
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ diagonalcalcForTE()

int diagonalcalcForTE ( const int  _istep,
struct BindStruct *  X,
double complex *  tmp_v0,
double complex *  tmp_v1 
)

Definition at line 197 of file diagonalcalc.c.

References BindStruct::Def, DefineList::NTEChemi, DefineList::NTEInterAllDiagonal, DefineList::NTETransferDiagonal, DefineList::ParaTEChemi, DefineList::ParaTEInterAllDiagonal, DefineList::ParaTETransferDiagonal, SetDiagonalCoulombIntra(), SetDiagonalTEChemi(), SetDiagonalTEInterAll(), SetDiagonalTETransfer(), DefineList::SpinTEChemi, DefineList::TEChemi, DefineList::TEInterAllDiagonal, and DefineList::TETransferDiagonal.

Referenced by diagonalcalc(), and mltply().

202  {
203 
204  long unsigned int i;
205  long unsigned int isite1, isite2;
206  long unsigned int A_spin, B_spin;
207  double tmp_V;
208 
209  if (X->Def.NTETransferDiagonal[_istep] > 0) {
210  for (i = 0; i < X->Def.NTETransferDiagonal[_istep]; i++) {
211  isite1 = X->Def.TETransferDiagonal[_istep][i][0] + 1;
212  A_spin = X->Def.TETransferDiagonal[_istep][i][1];
213  tmp_V = X->Def.ParaTETransferDiagonal[_istep][i];
214  SetDiagonalTETransfer(isite1, tmp_V, A_spin, X, tmp_v0, tmp_v1);
215  }
216  }
217  else if (X->Def.NTEInterAllDiagonal[_istep] >0) {
218  for (i = 0; i < X->Def.NTEInterAllDiagonal[_istep]; i++) {
219  //Assume n_{1\sigma_1} n_{2\sigma_2}
220  isite1 = X->Def.TEInterAllDiagonal[_istep][i][0] + 1;
221  A_spin = X->Def.TEInterAllDiagonal[_istep][i][1];
222  isite2 = X->Def.TEInterAllDiagonal[_istep][i][2] + 1;
223  B_spin = X->Def.TEInterAllDiagonal[_istep][i][3];
224  tmp_V = X->Def.ParaTEInterAllDiagonal[_istep][i];
225 
226  if (SetDiagonalTEInterAll(isite1, isite2, A_spin, B_spin, tmp_V, X, tmp_v0, tmp_v1) != 0) {
227  return -1;
228  }
229  }
230 
231  if (X->Def.NTEChemi[_istep] > 0) {
232  for(i=0; i< X->Def.NTEChemi[_istep]; i++) {
233  isite1 = X->Def.TEChemi[_istep][i] + 1;
234  A_spin = X->Def.SpinTEChemi[_istep][i];
235  tmp_V = -X->Def.ParaTEChemi[_istep][i];
236  if (SetDiagonalTEChemi(isite1, A_spin, tmp_V, X, tmp_v0, tmp_v1) != 0) {
237  return -1;
238  }
239  }
240  }
241  }
242  return 0;
243 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
DefineList::SpinTEChemi
int ** SpinTEChemi
Definition: struct.h:295
SetDiagonalTEChemi
int SetDiagonalTEChemi(long unsigned int isite1, long unsigned int spin, double dtmp_V, struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1)
Update the vector by the chemical potential generated by the commutation relation in terms of the g...
Definition: diagonalcalc.c:1762
DefineList::TETransferDiagonal
int *** TETransferDiagonal
Definition: struct.h:262
SetDiagonalTEInterAll
int SetDiagonalTEInterAll(long unsigned int isite1, long unsigned int isite2, long unsigned int isigma1, long unsigned int isigma2, double dtmp_V, struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1)
Update the vector by the general two-body diagonal interaction, . (Using in Time Evolution mode)...
Definition: diagonalcalc.c:1431
DefineList::TEChemi
int ** TEChemi
Definition: struct.h:293
DefineList::ParaTEChemi
double ** ParaTEChemi
Definition: struct.h:296
DefineList::NTEChemi
unsigned int * NTEChemi
Definition: struct.h:294
DefineList::ParaTETransferDiagonal
double ** ParaTETransferDiagonal
Definition: struct.h:266
SetDiagonalTETransfer
int SetDiagonalTETransfer(long unsigned int isite1, double dtmp_V, long unsigned int spin, struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1)
Update the vector by the general one-body diagonal interaction, . (Using in Time Evolution mode)...
Definition: diagonalcalc.c:1940
DefineList::NTEInterAllDiagonal
unsigned int * NTEInterAllDiagonal
Definition: struct.h:276
DefineList::ParaTEInterAllDiagonal
double ** ParaTEInterAllDiagonal
Definition: struct.h:291
DefineList::TEInterAllDiagonal
int *** TEInterAllDiagonal
Definition: struct.h:284
DefineList::NTETransferDiagonal
unsigned int * NTETransferDiagonal
Definition: struct.h:258
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalChemi()

int SetDiagonalChemi ( long unsigned int  isite1,
double  dtmp_V,
long unsigned int  spin,
struct BindStruct *  X 
)

Calculate the components for the chemical potential \( \mu_{i \sigma_1} n_ {i \sigma_1} \).

Parameters
isite1[in] a site number
dtmp_V[in] A value of coulombintra interaction \( \mu_{i \sigma_1} \)
spin[in] Spin index for the chemical potential
X[in] Define list to get dimension number
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
0.1
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 374 of file diagonalcalc.c.

References BitCheckGeneral(), cErrNoModel, BindStruct::Check, BindStruct::Def, FALSE, DefineList::iCalcModel, CheckList::idim_max, DefineList::iFlgGeneralSpin, list_1, list_Diagonal, myrank, DefineList::Nsite, SetDiagonalCoulombInter(), DefineList::SiteToBit, stdoutMPI, and DefineList::Tpow.

Referenced by diagonalcalc(), and SetDiagonalCoulombIntra().

379  {
380  long unsigned int is1_up;
381  long unsigned int ibit1_up;
382  long unsigned int num1;
383  long unsigned int isigma1 =spin;
384  long unsigned int is1,ibit1;
385 
386  long unsigned int j;
387  long unsigned int i_max=X->Check.idim_max;
388 
389  /*
390  When isite1 is in the inter process region
391  */
392  if (isite1 > X->Def.Nsite){
393 
394  switch (X->Def.iCalcModel) {
395 
396  case HubbardGC:
397  case KondoGC:
398  case Hubbard:
399  case Kondo:
400 
401  if (spin == 0) {
402  is1 = X->Def.Tpow[2 * isite1 - 2];
403  }
404  else {
405  is1 = X->Def.Tpow[2 * isite1 - 1];
406  }
407  ibit1 = (unsigned long int)myrank & is1;
408  num1 = ibit1 / is1;
409 #pragma omp parallel for default(none) shared(list_Diagonal) \
410  firstprivate(i_max, dtmp_V, num1) private(j)
411  for (j = 1; j <= i_max; j++) list_Diagonal[j] += num1*dtmp_V;
412 
413  break;/*case HubbardGC, case KondoGC, Hubbard, Kondo:*/
414 
415  case SpinGC:
416  case Spin:
417 
418  if (X->Def.iFlgGeneralSpin == FALSE) {
419  is1_up = X->Def.Tpow[isite1 - 1];
420  ibit1_up = (((unsigned long int)myrank& is1_up) / is1_up) ^ (1 - spin);
421 #pragma omp parallel for default(none) shared(list_Diagonal) \
422 firstprivate(i_max, dtmp_V, ibit1_up) private(j)
423  for (j = 1; j <= i_max; j++) list_Diagonal[j] += dtmp_V * ibit1_up;
424  } /*if (X->Def.iFlgGeneralSpin == FALSE)*/
425  else /*if (X->Def.iFlgGeneralSpin == TRUE)*/ {
426  num1 = BitCheckGeneral((unsigned long int)myrank,
427  isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
428  if (num1 != 0) {
429 #pragma omp parallel for default(none) shared(list_Diagonal) \
430 firstprivate(i_max, dtmp_V) private(j)
431  for (j = 1; j <= i_max; j++) list_Diagonal[j] += dtmp_V;
432  }/*if (num1 != 0)*/
433  }/*if (X->Def.iFlgGeneralSpin == TRUE)*/
434  break;/*case SpinGC, Spin:*/
435 
436  default:
437  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
438  return -1;
439 
440  } /*switch (X->Def.iCalcModel)*/
441 
442  return 0;
443 
444  }/*if (isite1 >= X->Def.Nsite*/
445 
446  switch (X->Def.iCalcModel){
447  case HubbardGC:
448  if(spin==0){
449  is1 = X->Def.Tpow[2*isite1-2];
450  }else{
451  is1 = X->Def.Tpow[2*isite1-1];
452  }
453 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1) private(num1, ibit1)
454  for(j = 1;j <= i_max;j++){
455 
456  ibit1 = (j-1)&is1;
457  num1 = ibit1/is1;
458  //fprintf(stdoutMPI, "DEBUG: spin=%ld is1=%ld: isite1=%ld j=%ld num1=%ld \n",spin,is1,isite1,j,num1);
459 
460  list_Diagonal[j]+=num1*dtmp_V;
461  }
462  break;
463  case KondoGC:
464  case Hubbard:
465  case Kondo:
466  if(spin==0){
467  is1 = X->Def.Tpow[2*isite1-2];
468  }else{
469  is1 = X->Def.Tpow[2*isite1-1];
470  }
471 
472 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1) private(num1, ibit1)
473  for(j = 1;j <= i_max;j++){
474 
475  ibit1 = list_1[j]&is1;
476  num1 = ibit1/is1;
477  list_Diagonal[j]+=num1*dtmp_V;
478  }
479  break;
480 
481  case SpinGC:
482  if(X->Def.iFlgGeneralSpin==FALSE){
483  is1_up = X->Def.Tpow[isite1-1];
484 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, spin) private(num1, ibit1_up)
485  for(j = 1;j <= i_max;j++){
486  ibit1_up=(((j-1)& is1_up)/is1_up)^(1-spin);
487  list_Diagonal[j] += dtmp_V * ibit1_up;
488  }
489  }
490  else{
491 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
492  for(j = 1;j <= i_max; j++){
493  num1=BitCheckGeneral (j-1, isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
494  if(num1 != 0){
495  list_Diagonal[j] += dtmp_V;
496  }
497  }
498  }
499  break;
500 
501  case Spin:
502  if(X->Def.iFlgGeneralSpin==FALSE){
503  is1_up = X->Def.Tpow[isite1-1];
504 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, spin) private(num1, ibit1_up)
505  for(j = 1;j <= i_max;j++){
506  ibit1_up=((list_1[j]& is1_up)/is1_up)^(1-spin);
507  list_Diagonal[j] += dtmp_V * ibit1_up;
508  }
509  }
510  else{
511 #pragma omp parallel for default(none) shared(list_Diagonal, list_1) firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
512  for(j = 1;j <= i_max; j++){
513  num1=BitCheckGeneral (list_1[j], isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
514  if(num1 != 0){
515  list_Diagonal[j] += dtmp_V;
516  }
517  }
518  }
519 
520  break;
521  default:
522  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
523  return -1;
524  }
525 
526  return 0;
527 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
DefineList::iFlgGeneralSpin
int iFlgGeneralSpin
Flag for the general (Sz/=1/2) spin.
Definition: struct.h:86
FALSE
#define FALSE
Definition: global.h:25
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
DefineList::SiteToBit
long int * SiteToBit
[DefineList::NsiteMPI] Similar to DefineList::Tpow. For general spin.
Definition: struct.h:94
BitCheckGeneral
int BitCheckGeneral(const long unsigned int org_bit, const unsigned int org_isite, const unsigned int target_ispin, const long int *SiteToBit, const long unsigned int *Tpow)
bit check function for general spin
Definition: bitcalc.c:393
list_Diagonal
double * list_Diagonal
Definition: global.h:46
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalCoulombInter()

int SetDiagonalCoulombInter ( long unsigned int  isite1,
long unsigned int  isite2,
double  dtmp_V,
struct BindStruct *  X 
)

Calculate the components for Coulombinter interaction, \( V_{ij} n_ {i}n_{j} \).

Parameters
isite1[in] a site number \(i \)
isite2[in] a site number \(j \)
dtmp_V[in] A value of coulombinter interaction \( V_{ij} \)
X[in] Define list to get the operator information.
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
0.1
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 544 of file diagonalcalc.c.

References cErrNoModel, BindStruct::Check, BindStruct::Def, DefineList::iCalcModel, CheckList::idim_max, list_1, list_Diagonal, myrank, DefineList::Nsite, SetDiagonalHund(), stdoutMPI, and DefineList::Tpow.

Referenced by diagonalcalc(), and SetDiagonalChemi().

550 {
551 
552  long unsigned int is1_up, is1_down;
553  long unsigned int ibit1_up, ibit1_down;
554  long unsigned int num1;
555  long unsigned int is2_up, is2_down;
556  long unsigned int ibit2_up, ibit2_down;
557  long unsigned int num2;
558 
559  long unsigned int j;
560  long unsigned int i_max=X->Check.idim_max;
561 
562  /*
563  Force isite1 <= isite2
564  */
565  if (isite2 < isite1) {
566  j = isite2;
567  isite2 = isite1;
568  isite1 = j;
569  }/*if (isite2 < isite1)*/
570  /*
571  When isite1 & site2 are in the inter process region
572  */
573  if (/*isite2 => */ isite1 > X->Def.Nsite) {
574 
575  switch (X->Def.iCalcModel) {
576 
577  case HubbardGC:
578  case KondoGC:
579  case Hubbard:
580  case Kondo:
581 
582  is1_up = X->Def.Tpow[2 * isite1 - 2];
583  is1_down = X->Def.Tpow[2 * isite1 - 1];
584  is2_up = X->Def.Tpow[2 * isite2 - 2];
585  is2_down = X->Def.Tpow[2 * isite2 - 1];
586 
587  num1 = 0;
588  num2 = 0;
589 
590  ibit1_up = (unsigned long int)myrank&is1_up;
591  num1 += ibit1_up / is1_up;
592  ibit1_down = (unsigned long int)myrank&is1_down;
593  num1 += ibit1_down / is1_down;
594 
595  ibit2_up = (unsigned long int)myrank&is2_up;
596  num2 += ibit2_up / is2_up;
597  ibit2_down = (unsigned long int)myrank&is2_down;
598  num2 += ibit2_down / is2_down;
599 
600 #pragma omp parallel for default(none) shared(list_Diagonal) \
601  firstprivate(i_max, dtmp_V, num1, num2) private(j)
602  for (j = 1; j <= i_max; j++) list_Diagonal[j] += num1*num2*dtmp_V;
603 
604  break;/*case HubbardGC, KondoGC, Hubbard, Kondo:*/
605 
606  case Spin:
607  case SpinGC:
608 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V)
609  for (j = 1; j <= i_max; j++) {
610  list_Diagonal[j] += dtmp_V;
611  }
612  break;/*case Spin, SpinGC*/
613 
614  default:
615  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
616  return -1;
617 
618  }/*switch (X->Def.iCalcModel)*/
619 
620  return 0;
621 
622  }/*if (isite1 > X->Def.Nsite)*/
623  else if (isite2 > X->Def.Nsite /* => isite1 */) {
624 
625  switch(X->Def.iCalcModel){
626  case HubbardGC:
627  case KondoGC:
628  case Hubbard:
629  case Kondo:
630  is1_up = X->Def.Tpow[2 * isite1 - 2];
631  is1_down = X->Def.Tpow[2 * isite1 - 1];
632  is2_up = X->Def.Tpow[2 * isite2 - 2];
633  is2_down = X->Def.Tpow[2 * isite2 - 1];
634  num2 = 0;
635  ibit2_up = (unsigned long int)myrank&is2_up;
636  num2 += ibit2_up / is2_up;
637  ibit2_down = (unsigned long int)myrank&is2_down;
638  num2 += ibit2_down / is2_down;
639  break;
640 
641  case Spin:
642  case SpinGC:
643  break;
644 
645  default:
646  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
647  return -1;
648  }
649 
650  switch (X->Def.iCalcModel) {
651 
652  case HubbardGC:
653 
654 #pragma omp parallel for default(none) shared(list_Diagonal) \
655 firstprivate(i_max, dtmp_V, num2, is1_up, is1_down) \
656 private(num1, ibit1_up, ibit1_down, j)
657  for (j = 1; j <= i_max; j++) {
658  num1 = 0;
659  ibit1_up = (j - 1)&is1_up;
660  num1 += ibit1_up / is1_up;
661  ibit1_down = (j - 1)&is1_down;
662  num1 += ibit1_down / is1_down;
663 
664  list_Diagonal[j] += num1*num2*dtmp_V;
665  }
666 
667  break;/*case HubbardGC*/
668 
669  case KondoGC:
670  case Hubbard:
671  case Kondo:
672 
673 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) \
674 firstprivate(i_max, dtmp_V, is1_up, is1_down, num2) \
675 private(num1, ibit1_up, ibit1_down, j)
676  for (j = 1; j <= i_max; j++) {
677  num1 = 0;
678  ibit1_up = list_1[j] & is1_up;
679  num1 += ibit1_up / is1_up;
680  ibit1_down = list_1[j] & is1_down;
681  num1 += ibit1_down / is1_down;
682 
683  list_Diagonal[j] += num1*num2*dtmp_V;
684  }
685  break;/*case KondoGC, Hubbard, Kondo:*/
686 
687  case Spin:
688  case SpinGC:
689 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V)
690  for (j = 1; j <= i_max; j++) {
691  list_Diagonal[j] += dtmp_V;
692  }
693  break;/* case Spin, SpinGC:*/
694 
695  default:
696  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
697  return -1;
698 
699  }/*switch (X->Def.iCalcModel)*/
700 
701  return 0;
702 
703  }/*else if (isite2 > X->Def.Nsite)*/
704  else{
705  switch (X->Def.iCalcModel){
706  case HubbardGC: //list_1[j] -> j-1
707  is1_up = X->Def.Tpow[2*isite1-2];
708  is1_down = X->Def.Tpow[2*isite1-1];
709  is2_up = X->Def.Tpow[2*isite2-2];
710  is2_down = X->Def.Tpow[2*isite2-1];
711 #pragma omp parallel for default(none) shared( list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is1_down, is2_up, is2_down) private(num1, ibit1_up, ibit1_down, num2, ibit2_up, ibit2_down)
712  for(j = 1;j <= i_max;j++){
713  num1=0;
714  num2=0;
715  ibit1_up=(j-1)&is1_up;
716  num1+=ibit1_up/is1_up;
717  ibit1_down=(j-1)&is1_down;
718  num1+=ibit1_down/is1_down;
719 
720  ibit2_up=(j-1)&is2_up;
721  num2+=ibit2_up/is2_up;
722  ibit2_down=(j-1)&is2_down;
723  num2+=ibit2_down/is2_down;
724 
725  list_Diagonal[j]+=num1*num2*dtmp_V;
726  }
727  break;
728  case KondoGC:
729  case Hubbard:
730  case Kondo:
731  is1_up = X->Def.Tpow[2*isite1-2];
732  is1_down = X->Def.Tpow[2*isite1-1];
733  is2_up = X->Def.Tpow[2*isite2-2];
734  is2_down = X->Def.Tpow[2*isite2-1];
735 
736 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is1_down, is2_up, is2_down) private(num1, ibit1_up, ibit1_down, num2, ibit2_up, ibit2_down)
737  for(j = 1;j <= i_max;j++){
738  num1=0;
739  num2=0;
740  ibit1_up=list_1[j]&is1_up;
741  num1+=ibit1_up/is1_up;
742  ibit1_down=list_1[j]&is1_down;
743  num1+=ibit1_down/is1_down;
744 
745  ibit2_up=list_1[j]&is2_up;
746  num2+=ibit2_up/is2_up;
747  ibit2_down=list_1[j]&is2_down;
748  num2+=ibit2_down/is2_down;
749 
750  list_Diagonal[j]+=num1*num2*dtmp_V;
751  }
752  break;
753 
754  case Spin:
755  case SpinGC:
756 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V)
757  for(j = 1;j <= i_max; j++){
758  list_Diagonal[j] += dtmp_V;
759  }
760  break;
761  default:
762  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
763  return -1;
764  }
765  }
766 
767  return 0;
768 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
list_Diagonal
double * list_Diagonal
Definition: global.h:46
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalCoulombIntra()

int SetDiagonalCoulombIntra ( long unsigned int  isite1,
double  dtmp_V,
struct BindStruct *  X 
)

Calculate the components for Coulombintra interaction, \( U_i n_ {i \uparrow}n_{i \downarrow} \).

Parameters
isite1[in] a site number
dtmp_V[in] A value of coulombintra interaction \( U_i \)
X[in] Define list to get dimension number
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
0.1
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 261 of file diagonalcalc.c.

References cErrNoModel, BindStruct::Check, BindStruct::Def, DefineList::iCalcModel, CheckList::idim_max, list_1, list_Diagonal, myrank, DefineList::Nsite, SetDiagonalChemi(), stdoutMPI, and DefineList::Tpow.

Referenced by diagonalcalc(), and diagonalcalcForTE().

265  {
266  long unsigned int is;
267  long unsigned int ibit;
268  long unsigned int is1_up, is1_down;
269 
270  long unsigned int j;
271  long unsigned int i_max=X->Check.idim_max;
272 
273  /*
274  When isite1 is in the inter process region
275  */
276  if (isite1 > X->Def.Nsite){
277 
278  switch (X->Def.iCalcModel) {
279 
280  case HubbardGC:
281  case KondoGC:
282  case Hubbard:
283  case Kondo:
284 
285  is1_up = X->Def.Tpow[2 * isite1 - 2];
286  is1_down = X->Def.Tpow[2 * isite1 - 1];
287  is = is1_up + is1_down;
288  ibit = (unsigned long int)myrank & is;
289  if (ibit == is) {
290 #pragma omp parallel for default(none) shared(list_Diagonal) \
291  firstprivate(i_max, dtmp_V) private(j)
292  for (j = 1; j <= i_max; j++) list_Diagonal[j] += dtmp_V;
293  }
294 
295  break; /*case HubbardGC, KondoGC, Hubbard, Kondo:*/
296 
297  case Spin:
298  case SpinGC:
299  /*
300  They do not have the Coulomb term
301  */
302  break;
303 
304  default:
305  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
306  return -1;
307  //break;
308 
309  }/*switch (X->Def.iCalcModel)*/
310 
311  return 0;
312 
313  }/*if (isite1 >= X->Def.Nsite*/
314  else{
315  switch (X->Def.iCalcModel){
316  case HubbardGC:
317  is1_up = X->Def.Tpow[2*isite1-2];
318  is1_down = X->Def.Tpow[2*isite1-1];
319  is=is1_up+is1_down;
320 #pragma omp parallel for default(none) shared(list_Diagonal, list_1) firstprivate(i_max, is, dtmp_V) private(ibit)
321  for(j = 1;j <= i_max;j++){
322  ibit=(j-1)&is;
323  if(ibit==is){
324  list_Diagonal[j]+=dtmp_V;
325  }
326  }
327 
328  break;
329  case KondoGC:
330  case Hubbard:
331  case Kondo:
332  is1_up = X->Def.Tpow[2*isite1-2];
333  is1_down = X->Def.Tpow[2*isite1-1];
334  is=is1_up+is1_down;
335 #pragma omp parallel for default(none) shared(list_Diagonal, list_1) firstprivate(i_max, is, dtmp_V) private(ibit)
336  for(j = 1;j <= i_max;j++){
337  ibit=list_1[j]&is;
338  if(ibit==is){
339  list_Diagonal[j]+=dtmp_V;
340  }
341  }
342  break;
343 
344  case Spin:
345  case SpinGC:
346  break;
347 
348  default:
349  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
350  return -1;
351  //break;
352  }
353  }
354  return 0;
355 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
list_Diagonal
double * list_Diagonal
Definition: global.h:46
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalHund()

int SetDiagonalHund ( long unsigned int  isite1,
long unsigned int  isite2,
double  dtmp_V,
struct BindStruct *  X 
)

Calculate the components for Hund interaction, \( H_{ij}(n_ {i\uparrow}n_{j\uparrow}+ n_ {i\downarrow}n_{j\downarrow})\).

Parameters
isite1[in] a site number \(i \)
isite2[in] a site number \(j \)
dtmp_V[in] A value of Hund interaction \( H_{ij} \)
X[in] Define list to get the operator information.
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
0.1
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 785 of file diagonalcalc.c.

References cErrNoModel, BindStruct::Check, BindStruct::Def, DefineList::iCalcModel, CheckList::idim_max, list_1, list_Diagonal, myrank, DefineList::Nsite, SetDiagonalInterAll(), stdoutMPI, and DefineList::Tpow.

Referenced by diagonalcalc(), and SetDiagonalCoulombInter().

790  {
791 
792  long unsigned int is1_up, is1_down;
793  long unsigned int ibit1_up, ibit1_down;
794  long unsigned int num1_up, num1_down;
795  long unsigned int is2_up, is2_down;
796  long unsigned int ibit2_up, ibit2_down;
797  long unsigned int num2_up, num2_down;
798 
799  long unsigned int is_up;
800  long unsigned int ibit;
801  long unsigned int j;
802  long unsigned int i_max=X->Check.idim_max;
803  /*
804  Force isite1 <= isite2
805  */
806  if (isite2 < isite1) {
807  j = isite2;
808  isite2 = isite1;
809  isite1 = j;
810  }
811  /*
812  When isite1 & site2 are in the inter process region
813  */
814  if (/*isite2 >= */ isite1 > X->Def.Nsite){
815 
816  switch (X->Def.iCalcModel) {
817 
818  case HubbardGC:
819  case KondoGC:
820  case Hubbard:
821  case Kondo:
822 
823  is1_up = X->Def.Tpow[2 * isite1 - 2];
824  is1_down = X->Def.Tpow[2 * isite1 - 1];
825  is2_up = X->Def.Tpow[2 * isite2 - 2];
826  is2_down = X->Def.Tpow[2 * isite2 - 1];
827 
828  num1_up = 0;
829  num1_down = 0;
830  num2_up = 0;
831  num2_down = 0;
832 
833  ibit1_up = (unsigned long int)myrank &is1_up;
834  num1_up = ibit1_up / is1_up;
835  ibit1_down = (unsigned long int)myrank &is1_down;
836  num1_down = ibit1_down / is1_down;
837 
838  ibit2_up = (unsigned long int)myrank &is2_up;
839  num2_up = ibit2_up / is2_up;
840  ibit2_down = (unsigned long int)myrank &is2_down;
841  num2_down = ibit2_down / is2_down;
842 
843 #pragma omp parallel for default(none) shared(list_Diagonal) \
844  firstprivate(i_max, dtmp_V, num1_up, num1_down, num2_up, num2_down) private(j)
845  for (j = 1; j <= i_max; j++)
846  list_Diagonal[j] += dtmp_V*(num1_up*num2_up + num1_down*num2_down);
847 
848  break;/*case HubbardGC, KondoGC, Hubbard, Kondo:*/
849 
850  case SpinGC:
851  case Spin:
852 
853  is1_up = X->Def.Tpow[isite1 - 1];
854  is2_up = X->Def.Tpow[isite2 - 1];
855  is_up = is1_up + is2_up;
856  ibit = (unsigned long int)myrank & is_up;
857  if (ibit == 0 || ibit == is_up) {
858 #pragma omp parallel for default(none) shared(list_Diagonal) \
859 firstprivate(i_max, dtmp_V) private(j)
860  for (j = 1; j <= i_max; j++) list_Diagonal[j] += dtmp_V;
861  }
862  break;/*case SpinGC, Spin:*/
863 
864  default:
865  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
866  return -1;
867  }
868 
869  return 0;
870 
871  }/*if (isite1 > X->Def.Nsite)*/
872  else if (isite2 > X->Def.Nsite /* >= isite1 */) {
873 
874  switch (X->Def.iCalcModel) {
875 
876  case HubbardGC:
877 
878  is1_up = X->Def.Tpow[2 * isite1 - 2];
879  is1_down = X->Def.Tpow[2 * isite1 - 1];
880  is2_up = X->Def.Tpow[2 * isite2 - 2];
881  is2_down = X->Def.Tpow[2 * isite2 - 1];
882 
883  num2_up = 0;
884  num2_down = 0;
885 
886  ibit2_up = (unsigned long int)myrank &is2_up;
887  num2_up = ibit2_up / is2_up;
888  ibit2_down = (unsigned long int)myrank &is2_down;
889  num2_down = ibit2_down / is2_down;
890 
891 #pragma omp parallel for default(none) shared( list_Diagonal) \
892 firstprivate(i_max, dtmp_V, num2_up, num2_down, is1_up, is1_down) \
893 private(num1_up, num1_down, ibit1_up, ibit1_down, j)
894  for (j = 1; j <= i_max; j++) {
895  num1_up = 0;
896  num1_down = 0;
897 
898  ibit1_up = (j - 1)&is1_up;
899  num1_up = ibit1_up / is1_up;
900  ibit1_down = (j - 1)&is1_down;
901  num1_down = ibit1_down / is1_down;
902 
903  list_Diagonal[j] += dtmp_V*(num1_up*num2_up + num1_down*num2_down);
904  }
905  break;/*case HubbardGC:*/
906 
907  case KondoGC:
908  case Hubbard:
909  case Kondo:
910 
911  is1_up = X->Def.Tpow[2 * isite1 - 2];
912  is1_down = X->Def.Tpow[2 * isite1 - 1];
913  is2_up = X->Def.Tpow[2 * isite2 - 2];
914  is2_down = X->Def.Tpow[2 * isite2 - 1];
915 
916  num2_up = 0;
917  num2_down = 0;
918 
919  ibit2_up = (unsigned long int)myrank&is2_up;
920  num2_up = ibit2_up / is2_up;
921  ibit2_down = (unsigned long int)myrank&is2_down;
922  num2_down = ibit2_down / is2_down;
923 
924 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) \
925 firstprivate(i_max, dtmp_V, num2_up, num2_down, is1_up, is1_down) \
926 private(num1_up, num1_down, ibit1_up, ibit1_down, j)
927  for (j = 1; j <= i_max; j++) {
928  num1_up = 0;
929  num1_down = 0;
930 
931  ibit1_up = list_1[j] & is1_up;
932  num1_up = ibit1_up / is1_up;
933  ibit1_down = list_1[j] & is1_down;
934  num1_down = ibit1_down / is1_down;
935 
936  list_Diagonal[j] += dtmp_V*(num1_up*num2_up + num1_down*num2_down);
937  }
938  break;/*case KondoGC, Hubbard, Kondo:*/
939 
940  case SpinGC:
941  is1_up = X->Def.Tpow[isite1 - 1];
942  is2_up = X->Def.Tpow[isite2 - 1];
943  ibit2_up = (unsigned long int)myrank & is2_up;
944 
945  if (ibit2_up == is2_up) {
946 #pragma omp parallel for default(none) shared(list_Diagonal) \
947 firstprivate(i_max, dtmp_V, is1_up) private(j, ibit1_up)
948  for (j = 1; j <= i_max; j++) {
949  ibit1_up = (j - 1) & is1_up;
950  if (ibit1_up == is1_up) {
951  list_Diagonal[j] += dtmp_V;
952  }
953  }
954  }
955  else if(ibit2_up == 0){
956 #pragma omp parallel for default(none) shared(list_Diagonal) \
957 firstprivate(i_max, dtmp_V, is1_up) private(j, ibit1_up)
958  for (j = 1; j <= i_max; j++) {
959  ibit1_up = (j - 1) & is1_up;
960  if (ibit1_up == 0) {
961  list_Diagonal[j] += dtmp_V;
962  }
963  }
964  }
965  break;/*case SpinGC:*/
966 
967  case Spin:
968  is1_up = X->Def.Tpow[isite1 - 1];
969  is2_up = X->Def.Tpow[isite2 - 1];
970  ibit2_up = (unsigned long int)myrank & is2_up;
971 
972  if (ibit2_up == is2_up) {
973 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) \
974 firstprivate(i_max, dtmp_V, is1_up) private(j, ibit1_up)
975  for (j = 1; j <= i_max; j++) {
976  ibit1_up = list_1[j] & is1_up;
977  if (ibit1_up == is1_up) {
978  list_Diagonal[j] += dtmp_V;
979  }
980  }
981  }
982  else if (ibit2_up == 0) {
983 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) \
984 firstprivate(i_max, dtmp_V, is1_up) private(j, ibit1_up)
985  for (j = 1; j <= i_max; j++) {
986  ibit1_up = list_1[j] & is1_up;
987  if (ibit1_up == 0) {
988  list_Diagonal[j] += dtmp_V;
989  }
990  }
991  }
992  break;/*case Spin:*/
993 
994  default:
995  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
996  return -1;
997 
998  }/*switch (X->Def.iCalcModel)*/
999 
1000  return 0;
1001 
1002  }/*else if (isite2 > X->Def.Nsite)*/
1003  else{
1004  switch (X->Def.iCalcModel){
1005  case HubbardGC: // list_1[j] -> j-1
1006  is1_up = X->Def.Tpow[2*isite1-2];
1007  is1_down = X->Def.Tpow[2*isite1-1];
1008  is2_up = X->Def.Tpow[2*isite2-2];
1009  is2_down = X->Def.Tpow[2*isite2-1];
1010 
1011 #pragma omp parallel for default(none) shared( list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is1_down, is2_up, is2_down) private(num1_up, num1_down, num2_up, num2_down, ibit1_up, ibit1_down, ibit2_up, ibit2_down)
1012  for(j = 1; j <= i_max;j++){
1013  num1_up=0;
1014  num1_down=0;
1015  num2_up=0;
1016  num2_down=0;
1017 
1018  ibit1_up=(j-1)&is1_up;
1019  num1_up=ibit1_up/is1_up;
1020  ibit1_down=(j-1)&is1_down;
1021  num1_down=ibit1_down/is1_down;
1022 
1023  ibit2_up=(j-1)&is2_up;
1024  num2_up=ibit2_up/is2_up;
1025  ibit2_down=(j-1)&is2_down;
1026  num2_down=ibit2_down/is2_down;
1027 
1028  list_Diagonal[j]+=dtmp_V*(num1_up*num2_up+num1_down*num2_down);
1029  }
1030  break;
1031  case KondoGC:
1032  case Hubbard:
1033  case Kondo:
1034  is1_up = X->Def.Tpow[2*isite1-2];
1035  is1_down = X->Def.Tpow[2*isite1-1];
1036  is2_up = X->Def.Tpow[2*isite2-2];
1037  is2_down = X->Def.Tpow[2*isite2-1];
1038 
1039 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is1_down, is2_up, is2_down) private(num1_up, num1_down, num2_up, num2_down, ibit1_up, ibit1_down, ibit2_up, ibit2_down)
1040  for(j = 1; j <= i_max;j++){
1041  num1_up=0;
1042  num1_down=0;
1043  num2_up=0;
1044  num2_down=0;
1045 
1046  ibit1_up=list_1[j]&is1_up;
1047  num1_up=ibit1_up/is1_up;
1048  ibit1_down=list_1[j]&is1_down;
1049  num1_down=ibit1_down/is1_down;
1050 
1051  ibit2_up=list_1[j]&is2_up;
1052  num2_up=ibit2_up/is2_up;
1053  ibit2_down=list_1[j]&is2_down;
1054  num2_down=ibit2_down/is2_down;
1055 
1056  list_Diagonal[j]+=dtmp_V*(num1_up*num2_up+num1_down*num2_down);
1057  }
1058  break;
1059 
1060  case SpinGC:
1061  is1_up = X->Def.Tpow[isite1-1];
1062  is2_up = X->Def.Tpow[isite2-1];
1063  is_up = is1_up+is2_up;
1064 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is2_up, is_up) private(j, ibit)
1065  for(j = 1;j <= i_max;j++){
1066  ibit = (j-1) & is_up;
1067  if(ibit == 0 || ibit == is_up){
1068  list_Diagonal[j]+= dtmp_V;
1069  }
1070  }
1071  break;
1072 
1073  case Spin:
1074  is1_up = X->Def.Tpow[isite1-1];
1075  is2_up = X->Def.Tpow[isite2-1];
1076  is_up = is1_up+is2_up;
1077 #pragma omp parallel for default(none) shared(list_1, list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is2_up, is_up) private(j, ibit)
1078  for(j = 1;j <= i_max;j++){
1079  ibit = list_1[j] & is_up;
1080  if(ibit == 0 || ibit == is_up){
1081  list_Diagonal[j]+= dtmp_V;
1082  }
1083  }
1084  break;
1085  default:
1086  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1087  return -1;
1088  }
1089  }
1090  return 0;
1091 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
list_Diagonal
double * list_Diagonal
Definition: global.h:46
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalInterAll()

int SetDiagonalInterAll ( long unsigned int  isite1,
long unsigned int  isite2,
long unsigned int  isigma1,
long unsigned int  isigma2,
double  dtmp_V,
struct BindStruct *  X 
)

Calculate the components for general two-body diagonal interaction, \( H_{i\sigma_1 j\sigma_2} n_ {i\sigma_1}n_{j\sigma_2}\).

Parameters
isite1[in] a site number \(i \)
isite2[in] a site number \(j \)
isigma1[in] a spin index at \(i \) site.
isigma2[in] a spin index at \(j \) site.
dtmp_V[in] A value of general two-body diagonal interaction \( H_{i\sigma_1 j\sigma_2} \)
X[in] Define list to get the operator information.
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
0.1
Author
Takahiro Misawa (The University of Tokyo)
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 1110 of file diagonalcalc.c.

References BitCheckGeneral(), cErrNoModel, BindStruct::Check, BindStruct::Def, FALSE, DefineList::iCalcModel, CheckList::idim_max, DefineList::iFlgGeneralSpin, list_1, list_Diagonal, myrank, DefineList::Nsite, DefineList::SiteToBit, stdoutMPI, DefineList::Tpow, X_Spin_CisAis(), and X_SpinGC_CisAis().

Referenced by diagonalcalc(), and SetDiagonalHund().

1118 {
1119  long unsigned int is1_spin;
1120  long unsigned int is2_spin;
1121  long unsigned int is1_up;
1122  long unsigned int is2_up;
1123 
1124  long unsigned int ibit1_spin;
1125  long unsigned int ibit2_spin;
1126 
1127  long unsigned int num1;
1128  long unsigned int num2;
1129 
1130  long unsigned int j;
1131  long unsigned int i_max=X->Check.idim_max;
1132 
1133  /*
1134  Forse isite1 <= isite2
1135  */
1136  if (isite2 < isite1) {
1137  j = isite2;
1138  isite2 = isite1;
1139  isite1 = j;
1140  j = isigma2;
1141  isigma2 = isigma1;
1142  isigma1 = j;
1143  }
1144  /*
1145  When isite1 & site2 are in the inter process regino
1146  */
1147  if (isite1 > X->Def.Nsite) {
1148 
1149  switch (X->Def.iCalcModel) {
1150 
1151  case HubbardGC:
1152  case KondoGC:
1153  case Hubbard:
1154  case Kondo:
1155 
1156  is1_spin = X->Def.Tpow[2 * isite1 - 2 + isigma1];
1157  is2_spin = X->Def.Tpow[2 * isite2 - 2 + isigma2];
1158 
1159  num1 = 0;
1160  ibit1_spin = (unsigned long int)myrank&is1_spin;
1161  num1 += ibit1_spin / is1_spin;
1162 
1163  num2 = 0;
1164  ibit2_spin = (unsigned long int)myrank&is2_spin;
1165  num2 += ibit2_spin / is2_spin;
1166 
1167 #pragma omp parallel for default(none) shared(list_Diagonal) \
1168 firstprivate(i_max, dtmp_V, num2, num1) private(ibit1_spin, j)
1169  for (j = 1; j <= i_max; j++) list_Diagonal[j] += num1*num2*dtmp_V;
1170 
1171  break;/*case HubbardGC, KondoGC, Hubbard, Kondo:*/
1172 
1173  case SpinGC:
1174  case Spin:
1175 
1176  if (X->Def.iFlgGeneralSpin == FALSE) {
1177  is1_up = X->Def.Tpow[isite1 - 1];
1178  is2_up = X->Def.Tpow[isite2 - 1];
1179  num1 = X_SpinGC_CisAis((unsigned long int)myrank + 1, X, is1_up, isigma1);
1180  num2 = X_SpinGC_CisAis((unsigned long int)myrank + 1, X, is2_up, isigma2);
1181 
1182 #pragma omp parallel for default(none) shared(list_Diagonal) \
1183 firstprivate(i_max, dtmp_V, is1_up, isigma1, X, num1, num2) private(j)
1184  for (j = 1; j <= i_max; j++) {
1185  list_Diagonal[j] += num1*num2*dtmp_V;
1186  }
1187  }/*if (X->Def.iFlgGeneralSpin == FALSE)*/
1188  else {//start:generalspin
1189  num1 = BitCheckGeneral((unsigned long int)myrank, isite1, isigma1,
1190  X->Def.SiteToBit, X->Def.Tpow);
1191  num2 = BitCheckGeneral((unsigned long int)myrank, isite2, isigma2,
1192  X->Def.SiteToBit, X->Def.Tpow);
1193  if (num1 !=0 && num2 != 0) {
1194 #pragma omp parallel for default(none) shared(list_Diagonal) \
1195 firstprivate(i_max, dtmp_V, num1, X) private(j)
1196  for (j = 1; j <= i_max; j++) list_Diagonal[j] += dtmp_V*num1;
1197  }
1198  }/*if (X->Def.iFlgGeneralSpin == TRUE)*/
1199 
1200  break;/*case SpinGC, Spin:*/
1201 
1202  default:
1203  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1204  return -1;
1205 
1206  }/*if (isite1 > X->Def.Nsite)*/
1207 
1208  return 0;
1209 
1210  }/*if (isite1 > X->Def.Nsite)*/
1211  else if (isite2 > X->Def.Nsite) {
1212 
1213  switch (X->Def.iCalcModel) {
1214 
1215  case HubbardGC:
1216 
1217  is1_spin = X->Def.Tpow[2 * isite1 - 2 + isigma1];
1218  is2_spin = X->Def.Tpow[2 * isite2 - 2 + isigma2];
1219 
1220  num2 = 0;
1221  ibit2_spin = (unsigned long int)myrank&is2_spin;
1222  num2 += ibit2_spin / is2_spin;
1223 
1224 #pragma omp parallel for default(none) shared(list_Diagonal) \
1225 firstprivate(i_max, dtmp_V, is1_spin, num2) private(num1, ibit1_spin, j)
1226  for (j = 1; j <= i_max; j++) {
1227  num1 = 0;
1228  ibit1_spin = (j - 1)&is1_spin;
1229  num1 += ibit1_spin / is1_spin;
1230  list_Diagonal[j] += num1*num2*dtmp_V;
1231  }
1232  break;/*case HubbardGC:*/
1233 
1234  case KondoGC:
1235  case Hubbard:
1236  case Kondo:
1237 
1238  is1_spin = X->Def.Tpow[2 * isite1 - 2 + isigma1];
1239  is2_spin = X->Def.Tpow[2 * isite2 - 2 + isigma2];
1240 
1241  num2 = 0;
1242  ibit2_spin = (unsigned long int)myrank&is2_spin;
1243  num2 += ibit2_spin / is2_spin;
1244 
1245 #pragma omp parallel for default(none) shared(list_Diagonal, list_1) \
1246 firstprivate(i_max, dtmp_V, is1_spin, num2) private(num1, ibit1_spin, j)
1247  for (j = 1; j <= i_max; j++) {
1248  num1 = 0;
1249  ibit1_spin = list_1[j] & is1_spin;
1250  num1 += ibit1_spin / is1_spin;
1251  list_Diagonal[j] += num1*num2*dtmp_V;
1252  }
1253  break;/*case KondoGC, Hubbard, Kondo:*/
1254 
1255  case SpinGC:
1256 
1257  if (X->Def.iFlgGeneralSpin == FALSE) {
1258  is1_up = X->Def.Tpow[isite1 - 1];
1259  is2_up = X->Def.Tpow[isite2 - 1];
1260  num2 = X_SpinGC_CisAis((unsigned long int)myrank + 1, X, is2_up, isigma2);
1261 
1262 #pragma omp parallel for default(none) shared(list_Diagonal) \
1263 firstprivate(i_max, dtmp_V, is1_up, isigma1, X, num2) private(j, num1)
1264  for (j = 1; j <= i_max; j++) {
1265  num1 = X_SpinGC_CisAis(j, X, is1_up, isigma1);
1266  list_Diagonal[j] += num1*num2*dtmp_V;
1267  }
1268  }/* if (X->Def.iFlgGeneralSpin == FALSE)*/
1269  else {//start:generalspin
1270  num2 = BitCheckGeneral((unsigned long int)myrank, isite2, isigma2,
1271  X->Def.SiteToBit, X->Def.Tpow);
1272  if (num2 != 0) {
1273 #pragma omp parallel for default(none) shared(list_Diagonal) \
1274 firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
1275  for (j = 1; j <= i_max; j++) {
1276  num1 = BitCheckGeneral(j - 1, isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1277  list_Diagonal[j] += dtmp_V*num1;
1278  }
1279  }
1280  }/* if (X->Def.iFlgGeneralSpin == TRUE)*/
1281 
1282  break;/*case SpinGC:*/
1283 
1284  case Spin:
1285 
1286  if (X->Def.iFlgGeneralSpin == FALSE) {
1287  is1_up = X->Def.Tpow[isite1 - 1];
1288  is2_up = X->Def.Tpow[isite2 - 1];
1289  num2 = X_SpinGC_CisAis((unsigned long int)myrank + 1, X, is2_up, isigma2);
1290 
1291 #pragma omp parallel for default(none) shared(list_Diagonal) \
1292 firstprivate(i_max, dtmp_V, is1_up, isigma1, X, num2) private(j, num1)
1293  for (j = 1; j <= i_max; j++) {
1294  num1 = X_Spin_CisAis(j, X, is1_up, isigma1);
1295  list_Diagonal[j] += num1*num2*dtmp_V;
1296  }
1297  }/* if (X->Def.iFlgGeneralSpin == FALSE)*/
1298  else /* if (X->Def.iFlgGeneralSpin == TRUE)*/{
1299  num2 = BitCheckGeneral((unsigned long int)myrank, isite2, isigma2, \
1300  X->Def.SiteToBit, X->Def.Tpow);
1301  if (num2 != 0) {
1302 #pragma omp parallel for default(none) shared(list_Diagonal, list_1) \
1303 firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
1304  for (j = 1; j <= i_max; j++) {
1305  num1 = BitCheckGeneral(list_1[j], isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1306  list_Diagonal[j] += dtmp_V*num1;
1307  }
1308  }
1309  } /* if (X->Def.iFlgGeneralSpin == TRUE)*/
1310 
1311  break;/*case Spin:*/
1312 
1313  default:
1314  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1315  return -1;
1316 
1317  }/*switch (X->Def.iCalcModel)*/
1318 
1319  return 0;
1320 
1321  }/*else if (isite2 > X->Def.Nsite)*/
1322 
1323  switch (X->Def.iCalcModel){
1324  case HubbardGC: //list_1[j] -> j-1
1325  is1_spin = X->Def.Tpow[2*isite1-2+isigma1];
1326  is2_spin = X->Def.Tpow[2*isite2-2+isigma2];
1327 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V, is1_spin, is2_spin) private(num1, ibit1_spin, num2, ibit2_spin)
1328  for(j = 1;j <= i_max;j++){
1329  num1=0;
1330  num2=0;
1331  ibit1_spin=(j-1)&is1_spin;
1332  num1+=ibit1_spin/is1_spin;
1333  ibit2_spin=(j-1)&is2_spin;
1334  num2+=ibit2_spin/is2_spin;
1335  list_Diagonal[j]+=num1*num2*dtmp_V;
1336  }
1337  break;
1338  case KondoGC:
1339  case Hubbard:
1340  case Kondo:
1341  is1_spin = X->Def.Tpow[2*isite1-2+isigma1];
1342  is2_spin = X->Def.Tpow[2*isite2-2+isigma2];
1343 
1344 #pragma omp parallel for default(none) shared(list_Diagonal, list_1) firstprivate(i_max, dtmp_V, is1_spin, is2_spin) private(num1, ibit1_spin, num2, ibit2_spin)
1345  for(j = 1;j <= i_max;j++){
1346  num1=0;
1347  num2=0;
1348  ibit1_spin=list_1[j]&is1_spin;
1349  num1+=ibit1_spin/is1_spin;
1350 
1351  ibit2_spin=list_1[j]&is2_spin;
1352  num2+=ibit2_spin/is2_spin;
1353  list_Diagonal[j]+=num1*num2*dtmp_V;
1354  }
1355  break;
1356 
1357  case Spin:
1358  if(X->Def.iFlgGeneralSpin==FALSE){
1359  is1_up = X->Def.Tpow[isite1-1];
1360  is2_up = X->Def.Tpow[isite2-1];
1361 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is2_up, isigma1, isigma2, X) private(j, num1, num2)
1362  for(j = 1;j <= i_max; j++){
1363  num1=X_Spin_CisAis(j, X, is1_up, isigma1);
1364  num2=X_Spin_CisAis(j, X, is2_up, isigma2);
1365  list_Diagonal[j] += num1*num2*dtmp_V;
1366  }
1367  }
1368  else{
1369 #pragma omp parallel for default(none) shared(list_Diagonal, list_1) firstprivate(i_max, dtmp_V, isite1, isite2, isigma1, isigma2, X) private(j, num1)
1370  for(j = 1;j <= i_max; j++){
1371  num1=BitCheckGeneral (list_1[j], isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1372  if(num1 != 0){
1373  num1=BitCheckGeneral (list_1[j], isite2, isigma2, X->Def.SiteToBit, X->Def.Tpow);
1374  list_Diagonal[j] += dtmp_V*num1;
1375  }
1376  }
1377 
1378  }
1379  break;
1380 
1381  case SpinGC:
1382  if(X->Def.iFlgGeneralSpin==FALSE){
1383  is1_up = X->Def.Tpow[isite1-1];
1384  is2_up = X->Def.Tpow[isite2-1];
1385 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V, is1_up, is2_up, isigma1, isigma2, X) private(j, num1, num2)
1386  for(j = 1;j <= i_max; j++){
1387  num1=X_SpinGC_CisAis(j, X, is1_up, isigma1);
1388  num2=X_SpinGC_CisAis(j, X, is2_up, isigma2);
1389  list_Diagonal[j] += num1*num2*dtmp_V;
1390  }
1391  }
1392  else{//start:generalspin
1393 #pragma omp parallel for default(none) shared(list_Diagonal) firstprivate(i_max, dtmp_V, isite1, isite2, isigma1, isigma2, X) private(j, num1)
1394  for(j = 1;j <= i_max; j++){
1395  num1=BitCheckGeneral (j-1, isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1396  if(num1 != 0){
1397  num1=BitCheckGeneral (j-1, isite2, isigma2, X->Def.SiteToBit, X->Def.Tpow);
1398  list_Diagonal[j] += dtmp_V*num1;
1399  }
1400  }
1401  }
1402  break;
1403 
1404  default:
1405  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1406  return -1;
1407  }
1408 
1409  return 0;
1410 
1411 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
X_SpinGC_CisAis
int X_SpinGC_CisAis(long unsigned int j, struct BindStruct *X, long unsigned int is1_spin, long unsigned int sigma1)
Compute the grandcanonical spin state with bit mask is1_spin.
Definition: mltplySpinCore.c:183
X_Spin_CisAis
int X_Spin_CisAis(long unsigned int j, struct BindStruct *X, long unsigned int is1_spin, long unsigned int sigma1)
Compute the spin state with bit mask is1_spin.
Definition: mltplySpinCore.c:166
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
DefineList::iFlgGeneralSpin
int iFlgGeneralSpin
Flag for the general (Sz/=1/2) spin.
Definition: struct.h:86
FALSE
#define FALSE
Definition: global.h:25
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
DefineList::SiteToBit
long int * SiteToBit
[DefineList::NsiteMPI] Similar to DefineList::Tpow. For general spin.
Definition: struct.h:94
BitCheckGeneral
int BitCheckGeneral(const long unsigned int org_bit, const unsigned int org_isite, const unsigned int target_ispin, const long int *SiteToBit, const long unsigned int *Tpow)
bit check function for general spin
Definition: bitcalc.c:393
list_Diagonal
double * list_Diagonal
Definition: global.h:46
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalTEChemi()

int SetDiagonalTEChemi ( long unsigned int  isite1,
long unsigned int  spin,
double  dtmp_V,
struct BindStruct *  X,
double complex *  tmp_v0,
double complex *  tmp_v1 
)

Update the vector by the chemical potential \( \mu_{i \sigma_1} n_ {i \sigma_1} \)
generated by the commutation relation in terms of the general two-body interaction,
\( c_ {i \sigma_1} a_{j\sigma_2}c_ {j \sigma_2}a_ {i \sigma_1} = c_ {i \sigma_1}a_ {i \sigma_1}-c_ {i \sigma_1} a_ {i \sigma_1} c_ {j \sigma_2}a_{j\sigma_2}\) . (Using in Time Evolution mode).

Parameters
isite1[in] a site number
spin[in] a spin number
dtmp_V[in] A value of coulombintra interaction \( \mu_{i \sigma_1} \)
X[in] Define list to get dimension number
tmp_v0[in,out] Result vector
tmp_v1[in] Input produced vector
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
2.1
Author
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 1762 of file diagonalcalc.c.

References BitCheckGeneral(), cErrNoModel, BindStruct::Check, BindStruct::Def, FALSE, DefineList::iCalcModel, CheckList::idim_max, DefineList::iFlgGeneralSpin, BindStruct::Large, list_1, myrank, DefineList::Nsite, LargeList::prdct, SetDiagonalTETransfer(), DefineList::SiteToBit, stdoutMPI, SumMPI_dc(), and DefineList::Tpow.

Referenced by diagonalcalcForTE().

1769  {
1770  long unsigned int is1_up;
1771  long unsigned int num1;
1772  long unsigned int isigma1 =spin;
1773  long unsigned int is1,ibit1;
1774 
1775  long unsigned int j;
1776  long unsigned int i_max=X->Check.idim_max;
1777  double complex dam_pr=0;
1778 
1779  /*
1780  When isite1 is in the inter process region
1781  */
1782  if (isite1 > X->Def.Nsite){
1783 
1784  switch (X->Def.iCalcModel) {
1785 
1786  case HubbardGC:
1787  case KondoGC:
1788  case Hubbard:
1789  case Kondo:
1790 
1791  if (spin == 0) {
1792  is1 = X->Def.Tpow[2 * isite1 - 2];
1793  }
1794  else {
1795  is1 = X->Def.Tpow[2 * isite1 - 1];
1796  }
1797  ibit1 = (unsigned long int)myrank & is1;
1798  num1 = ibit1 / is1;
1799  break;/*case HubbardGC, case KondoGC, Hubbard, Kondo:*/
1800 
1801  case SpinGC:
1802  case Spin:
1803 
1804  if (X->Def.iFlgGeneralSpin == FALSE) {
1805  is1_up = X->Def.Tpow[isite1 - 1];
1806  num1 = (((unsigned long int)myrank& is1_up) / is1_up) ^ (1 - spin);
1807  } /*if (X->Def.iFlgGeneralSpin == FALSE)*/
1808  else /*if (X->Def.iFlgGeneralSpin == TRUE)*/ {
1809  num1 = BitCheckGeneral((unsigned long int)myrank,
1810  isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1811  }/*if (X->Def.iFlgGeneralSpin == TRUE)*/
1812  break;/*case SpinGC, Spin:*/
1813 
1814  default:
1815  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1816  return -1;
1817 
1818  } /*switch (X->Def.iCalcModel)*/
1819  if (num1 != 0) {
1820 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) \
1821 firstprivate(i_max, dtmp_V) private(j)
1822  for (j = 1; j <= i_max; j++){
1823  tmp_v0[j] += dtmp_V * tmp_v1[j];
1824  dam_pr += dtmp_V * conj(tmp_v1[j]) * tmp_v1[j];
1825  }
1826  }/*if (num1 != 0)*/
1827  dam_pr=SumMPI_dc(dam_pr);
1828  X->Large.prdct += dam_pr;
1829  return 0;
1830 
1831  }/*if (isite1 >= X->Def.Nsite*/
1832 
1833  switch (X->Def.iCalcModel){
1834  case HubbardGC:
1835  if(spin==0){
1836  is1 = X->Def.Tpow[2*isite1-2];
1837  }else{
1838  is1 = X->Def.Tpow[2*isite1-1];
1839  }
1840 
1841  #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, is1) private(num1, ibit1)
1842  for(j = 1;j <= i_max;j++){
1843  ibit1 = (j-1)&is1;
1844  num1 = ibit1/is1;
1845  tmp_v0[j] += dtmp_V * num1*tmp_v1[j];
1846  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1847  }
1848  break;
1849  case KondoGC:
1850  case Hubbard:
1851  case Kondo:
1852  if(spin==0){
1853  is1 = X->Def.Tpow[2*isite1-2];
1854  }else{
1855  is1 = X->Def.Tpow[2*isite1-1];
1856  }
1857 
1858 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, is1) private(num1, ibit1)
1859  for(j = 1;j <= i_max;j++){
1860 
1861  ibit1 = list_1[j]&is1;
1862  num1 = ibit1/is1;
1863  tmp_v0[j] += dtmp_V * num1*tmp_v1[j];
1864  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1865  }
1866  break;
1867 
1868  case SpinGC:
1869  if(X->Def.iFlgGeneralSpin==FALSE){
1870  is1_up = X->Def.Tpow[isite1-1];
1871 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, is1_up, spin) private(num1)
1872  for(j = 1;j <= i_max;j++){
1873  num1=(((j-1)& is1_up)/is1_up)^(1-spin);
1874  tmp_v0[j] += dtmp_V * num1*tmp_v1[j];
1875  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1876  }
1877  }
1878  else{
1879 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
1880  for(j = 1;j <= i_max; j++){
1881  num1=BitCheckGeneral (j-1, isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1882  if(num1 != 0){
1883  tmp_v0[j] += dtmp_V * tmp_v1[j];
1884  dam_pr += dtmp_V * conj(tmp_v1[j]) * tmp_v1[j];
1885  }
1886  }
1887  }
1888  break;
1889 
1890  case Spin:
1891  if(X->Def.iFlgGeneralSpin==FALSE){
1892  is1_up = X->Def.Tpow[isite1-1];
1893 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, is1_up, spin) private(num1)
1894  for(j = 1;j <= i_max;j++){
1895  num1=((list_1[j]& is1_up)/is1_up)^(1-spin);
1896  tmp_v0[j] += dtmp_V * num1*tmp_v1[j];
1897  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1898  }
1899  }
1900  else{
1901 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1, list_1) firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
1902  for(j = 1;j <= i_max; j++){
1903  num1=BitCheckGeneral (list_1[j], isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1904  if(num1 != 0){
1905  tmp_v0[j] += dtmp_V * tmp_v1[j];
1906  dam_pr += dtmp_V * conj(tmp_v1[j]) * tmp_v1[j];
1907 
1908  }
1909  }
1910  }
1911 
1912  break;
1913  default:
1914  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1915  return -1;
1916  }
1917  dam_pr=SumMPI_dc(dam_pr);
1918  X->Large.prdct += dam_pr;
1919  return 0;
1920 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
SumMPI_dc
double complex SumMPI_dc(double complex norm)
MPI wrapper function to obtain sum of Double complex across processes.
Definition: wrapperMPI.c:205
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
LargeList::prdct
double complex prdct
The expectation value of the energy.
Definition: struct.h:314
BindStruct::Large
struct LargeList Large
Variables for Matrix-Vector product.
Definition: struct.h:412
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
DefineList::iFlgGeneralSpin
int iFlgGeneralSpin
Flag for the general (Sz/=1/2) spin.
Definition: struct.h:86
FALSE
#define FALSE
Definition: global.h:25
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
DefineList::SiteToBit
long int * SiteToBit
[DefineList::NsiteMPI] Similar to DefineList::Tpow. For general spin.
Definition: struct.h:94
BitCheckGeneral
int BitCheckGeneral(const long unsigned int org_bit, const unsigned int org_isite, const unsigned int target_ispin, const long int *SiteToBit, const long unsigned int *Tpow)
bit check function for general spin
Definition: bitcalc.c:393
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalTEInterAll()

int SetDiagonalTEInterAll ( long unsigned int  isite1,
long unsigned int  isite2,
long unsigned int  isigma1,
long unsigned int  isigma2,
double  dtmp_V,
struct BindStruct *  X,
double complex *  tmp_v0,
double complex *  tmp_v1 
)

Update the vector by the general two-body diagonal interaction, \( H_{i\sigma_1 j\sigma_2} n_ {i\sigma_1}n_{j\sigma_2}\).
(Using in Time Evolution mode).

Parameters
isite1[in] a site number \(i \)
isite2[in] a site number \(j \)
isigma1[in] a spin index at \(i \) site.
isigma2[in] a spin index at \(j \) site.
dtmp_V[in] A value of general two-body diagonal interaction \( H_{i\sigma_1 j\sigma_2} \)
X[in] Define list to get the operator information.
tmp_v0[in,out] Result vector
tmp_v1[in] Input produced vector
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
2.1
Author
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 1431 of file diagonalcalc.c.

References BitCheckGeneral(), cErrNoModel, BindStruct::Check, BindStruct::Def, FALSE, DefineList::iCalcModel, CheckList::idim_max, DefineList::iFlgGeneralSpin, BindStruct::Large, list_1, myrank, DefineList::Nsite, LargeList::prdct, DefineList::SiteToBit, stdoutMPI, SumMPI_dc(), DefineList::Tpow, X_Spin_CisAis(), and X_SpinGC_CisAis().

Referenced by diagonalcalcForTE().

1440  {
1441  long unsigned int is1_spin;
1442  long unsigned int is2_spin;
1443  long unsigned int is1_up;
1444  long unsigned int is2_up;
1445 
1446  long unsigned int ibit1_spin;
1447  long unsigned int ibit2_spin;
1448 
1449  long unsigned int num1;
1450  long unsigned int num2;
1451 
1452  long unsigned int j;
1453  long unsigned int i_max=X->Check.idim_max;
1454  double complex dam_pr=0.0;
1455 
1456 
1457  /*
1458  Forse isite1 <= isite2
1459  */
1460  if (isite2 < isite1) {
1461  j = isite2;
1462  isite2 = isite1;
1463  isite1 = j;
1464  j = isigma2;
1465  isigma2 = isigma1;
1466  isigma1 = j;
1467  }
1468  /*
1469  When isite1 & site2 are in the inter process regino
1470  */
1471  if (isite1 > X->Def.Nsite) {
1472 
1473  switch (X->Def.iCalcModel) {
1474 
1475  case HubbardGC:
1476  case KondoGC:
1477  case Hubbard:
1478  case Kondo:
1479  is1_spin = X->Def.Tpow[2 * isite1 - 2 + isigma1];
1480  is2_spin = X->Def.Tpow[2 * isite2 - 2 + isigma2];
1481  num1 = 0;
1482  ibit1_spin = (unsigned long int)myrank&is1_spin;
1483  num1 += ibit1_spin / is1_spin;
1484  num2 = 0;
1485  ibit2_spin = (unsigned long int)myrank&is2_spin;
1486  num2 += ibit2_spin / is2_spin;
1487  break;/*case HubbardGC, KondoGC, Hubbard, Kondo:*/
1488 
1489  case SpinGC:
1490  case Spin:
1491  if (X->Def.iFlgGeneralSpin == FALSE) {
1492  is1_up = X->Def.Tpow[isite1 - 1];
1493  is2_up = X->Def.Tpow[isite2 - 1];
1494  num1 = X_SpinGC_CisAis((unsigned long int) myrank + 1, X, is1_up, isigma1);
1495  num2 = X_SpinGC_CisAis((unsigned long int) myrank + 1, X, is2_up, isigma2);
1496  }/*if (X->Def.iFlgGeneralSpin == FALSE)*/
1497  else {//start:generalspin
1498  num1 = BitCheckGeneral((unsigned long int) myrank, isite1, isigma1,
1499  X->Def.SiteToBit, X->Def.Tpow);
1500  num2 = BitCheckGeneral((unsigned long int) myrank, isite2, isigma2,
1501  X->Def.SiteToBit, X->Def.Tpow);
1502  }
1503  break;/*case SpinGC, Spin:*/
1504 
1505  default:
1506  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1507  return -1;
1508  }/*if (isite1 > X->Def.Nsite)*/
1509 
1510  if (num1 * num2 != 0) {
1511 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) \
1512 firstprivate(i_max, dtmp_V) private(j)
1513  for (j = 1; j <= i_max; j++) {
1514  tmp_v0[j] += dtmp_V * tmp_v1[j];
1515  dam_pr += dtmp_V * conj(tmp_v1[j]) * tmp_v1[j];
1516  }
1517  }
1518  dam_pr=SumMPI_dc(dam_pr);
1519  X->Large.prdct += dam_pr;
1520  return 0;
1521 
1522  }/*if (isite1 > X->Def.Nsite)*/
1523  else if (isite2 > X->Def.Nsite) {
1524 
1525  switch (X->Def.iCalcModel) {
1526 
1527  case HubbardGC:
1528 
1529  is1_spin = X->Def.Tpow[2 * isite1 - 2 + isigma1];
1530  is2_spin = X->Def.Tpow[2 * isite2 - 2 + isigma2];
1531 
1532  num2 = 0;
1533  ibit2_spin = (unsigned long int)myrank&is2_spin;
1534  num2 += ibit2_spin / is2_spin;
1535  if(num2 !=0) {
1536 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1)\
1537  firstprivate(i_max, dtmp_V, is1_spin) private(num1, ibit1_spin, j)
1538  for (j = 1; j <= i_max; j++) {
1539  num1 = 0;
1540  ibit1_spin = (j - 1) & is1_spin;
1541  num1 += ibit1_spin / is1_spin;
1542  tmp_v0[j] += dtmp_V * num1 * tmp_v1[j];
1543  dam_pr += dtmp_V * num1 * conj(tmp_v1[j]) * tmp_v1[j];
1544  }
1545  }
1546  break;/*case HubbardGC:*/
1547 
1548  case KondoGC:
1549  case Hubbard:
1550  case Kondo:
1551 
1552  is1_spin = X->Def.Tpow[2 * isite1 - 2 + isigma1];
1553  is2_spin = X->Def.Tpow[2 * isite2 - 2 + isigma2];
1554 
1555  num2 = 0;
1556  ibit2_spin = (unsigned long int)myrank&is2_spin;
1557  num2 += ibit2_spin / is2_spin;
1558  if(num2 !=0) {
1559 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1, list_1)\
1560  firstprivate(i_max, dtmp_V, is1_spin) private(num1, ibit1_spin, j)
1561  for (j = 1; j <= i_max; j++) {
1562  num1 = 0;
1563  ibit1_spin = list_1[j] & is1_spin;
1564  num1 += ibit1_spin / is1_spin;
1565  tmp_v0[j] += dtmp_V *num1*tmp_v1[j];
1566  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1567  }
1568  }
1569  break;/*case KondoGC, Hubbard, Kondo:*/
1570 
1571  case SpinGC:
1572 
1573  if (X->Def.iFlgGeneralSpin == FALSE) {
1574  is1_up = X->Def.Tpow[isite1 - 1];
1575  is2_up = X->Def.Tpow[isite2 - 1];
1576  num2 = X_SpinGC_CisAis((unsigned long int)myrank + 1, X, is2_up, isigma2);
1577 
1578  if(num2 !=0) {
1579 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1)\
1580  firstprivate(i_max, dtmp_V, is1_up, isigma1, X) private(num1, j)
1581  for (j = 1; j <= i_max; j++) {
1582  num1 = X_SpinGC_CisAis(j, X, is1_up, isigma1);
1583  tmp_v0[j] += dtmp_V * num1 * tmp_v1[j];
1584  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1585  }
1586  }
1587  }/* if (X->Def.iFlgGeneralSpin == FALSE)*/
1588  else {//start:generalspin
1589  num2 = BitCheckGeneral((unsigned long int)myrank, isite2, isigma2,
1590  X->Def.SiteToBit, X->Def.Tpow);
1591  if (num2 != 0) {
1592 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) \
1593 firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
1594  for (j = 1; j <= i_max; j++) {
1595  num1 = BitCheckGeneral(j - 1, isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1596  tmp_v0[j] += dtmp_V * num1 * tmp_v1[j];
1597  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1598  }
1599  }
1600  }/* if (X->Def.iFlgGeneralSpin == TRUE)*/
1601 
1602  break;/*case SpinGC:*/
1603 
1604  case Spin:
1605 
1606  if (X->Def.iFlgGeneralSpin == FALSE) {
1607  is1_up = X->Def.Tpow[isite1 - 1];
1608  is2_up = X->Def.Tpow[isite2 - 1];
1609  num2 = X_SpinGC_CisAis((unsigned long int)myrank + 1, X, is2_up, isigma2);
1610 
1611  if(num2 !=0) {
1612 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) \
1613 firstprivate(i_max, dtmp_V, is1_up, isigma1, X, num2) private(j, num1)
1614  for (j = 1; j <= i_max; j++) {
1615  num1 = X_Spin_CisAis(j, X, is1_up, isigma1);
1616  tmp_v0[j] += dtmp_V * num1 * tmp_v1[j];
1617  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1618  }
1619  }
1620  }/* if (X->Def.iFlgGeneralSpin == FALSE)*/
1621  else /* if (X->Def.iFlgGeneralSpin == TRUE)*/{
1622  num2 = BitCheckGeneral((unsigned long int)myrank, isite2, isigma2, \
1623  X->Def.SiteToBit, X->Def.Tpow);
1624  if (num2 != 0) {
1625 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1, list_1)\
1626 firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
1627  for (j = 1; j <= i_max; j++) {
1628  num1 = BitCheckGeneral(list_1[j], isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1629  tmp_v0[j] += dtmp_V * num1 * tmp_v1[j];
1630  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1631  }
1632  }
1633  } /* if (X->Def.iFlgGeneralSpin == TRUE)*/
1634 
1635  break;/*case Spin:*/
1636 
1637  default:
1638  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1639  return -1;
1640 
1641  }/*switch (X->Def.iCalcModel)*/
1642  dam_pr=SumMPI_dc(dam_pr);
1643  X->Large.prdct += dam_pr;
1644  return 0;
1645  }/*else if (isite2 > X->Def.Nsite)*/
1646 
1647  switch (X->Def.iCalcModel){
1648  case HubbardGC: //list_1[j] -> j-1
1649  is1_spin = X->Def.Tpow[2*isite1-2+isigma1];
1650  is2_spin = X->Def.Tpow[2*isite2-2+isigma2];
1651 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, is1_spin, is2_spin) private(num1, ibit1_spin, num2, ibit2_spin)
1652  for(j = 1;j <= i_max;j++){
1653  num1=0;
1654  num2=0;
1655  ibit1_spin=(j-1)&is1_spin;
1656  num1+=ibit1_spin/is1_spin;
1657  ibit2_spin=(j-1)&is2_spin;
1658  num2+=ibit2_spin/is2_spin;
1659  tmp_v0[j] += dtmp_V * num1*num2*tmp_v1[j];
1660  dam_pr += dtmp_V * num1*num2*conj(tmp_v1[j]) * tmp_v1[j];
1661  }
1662  break;
1663  case KondoGC:
1664  case Hubbard:
1665  case Kondo:
1666  is1_spin = X->Def.Tpow[2*isite1-2+isigma1];
1667  is2_spin = X->Def.Tpow[2*isite2-2+isigma2];
1668 
1669 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1, list_1) firstprivate(i_max, dtmp_V, is1_spin, is2_spin) private(num1, ibit1_spin, num2, ibit2_spin)
1670  for(j = 1;j <= i_max;j++){
1671  num1=0;
1672  num2=0;
1673  ibit1_spin=list_1[j]&is1_spin;
1674  num1+=ibit1_spin/is1_spin;
1675 
1676  ibit2_spin=list_1[j]&is2_spin;
1677  num2+=ibit2_spin/is2_spin;
1678  tmp_v0[j] += dtmp_V * num1*num2*tmp_v1[j];
1679  dam_pr += dtmp_V * num1*num2*conj(tmp_v1[j]) * tmp_v1[j];
1680  }
1681  break;
1682 
1683  case Spin:
1684  if(X->Def.iFlgGeneralSpin==FALSE){
1685  is1_up = X->Def.Tpow[isite1-1];
1686  is2_up = X->Def.Tpow[isite2-1];
1687 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, is1_up, is2_up, isigma1, isigma2, X) private(j, num1, num2)
1688  for(j = 1;j <= i_max; j++){
1689  num1=X_Spin_CisAis(j, X, is1_up, isigma1);
1690  num2=X_Spin_CisAis(j, X, is2_up, isigma2);
1691  tmp_v0[j] += dtmp_V * num1*num2*tmp_v1[j];
1692  dam_pr += dtmp_V * num1*num2*conj(tmp_v1[j]) * tmp_v1[j];
1693  }
1694  }
1695  else{
1696 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1, list_1) firstprivate(i_max, dtmp_V, isite1, isite2, isigma1, isigma2, X) private(j, num1)
1697  for(j = 1;j <= i_max; j++){
1698  num1=BitCheckGeneral (list_1[j], isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1699  if(num1 != 0){
1700  num1=BitCheckGeneral (list_1[j], isite2, isigma2, X->Def.SiteToBit, X->Def.Tpow);
1701  tmp_v0[j] += dtmp_V *num1*tmp_v1[j];
1702  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1703  }
1704  }
1705 
1706  }
1707  break;
1708 
1709  case SpinGC:
1710  if(X->Def.iFlgGeneralSpin==FALSE){
1711  is1_up = X->Def.Tpow[isite1-1];
1712  is2_up = X->Def.Tpow[isite2-1];
1713 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, is1_up, is2_up, isigma1, isigma2, X) private(j, num1, num2)
1714  for(j = 1;j <= i_max; j++){
1715  num1=X_SpinGC_CisAis(j, X, is1_up, isigma1);
1716  num2=X_SpinGC_CisAis(j, X, is2_up, isigma2);
1717  tmp_v0[j] += dtmp_V * num1*num2*tmp_v1[j];
1718  dam_pr += dtmp_V * num1*num2*conj(tmp_v1[j]) * tmp_v1[j];
1719  }
1720  }
1721  else{//start:generalspin
1722 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) firstprivate(i_max, dtmp_V, isite1, isite2, isigma1, isigma2, X) private(j, num1)
1723  for(j = 1;j <= i_max; j++){
1724  num1=BitCheckGeneral (j-1, isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1725  if(num1 != 0){
1726  num1=BitCheckGeneral (j-1, isite2, isigma2, X->Def.SiteToBit, X->Def.Tpow);
1727  tmp_v0[j] += dtmp_V *num1*tmp_v1[j];
1728  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
1729  }
1730  }
1731  }
1732  break;
1733 
1734  default:
1735  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1736  return -1;
1737  }
1738  dam_pr=SumMPI_dc(dam_pr);
1739  X->Large.prdct += dam_pr;
1740  return 0;
1741 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
SumMPI_dc
double complex SumMPI_dc(double complex norm)
MPI wrapper function to obtain sum of Double complex across processes.
Definition: wrapperMPI.c:205
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
LargeList::prdct
double complex prdct
The expectation value of the energy.
Definition: struct.h:314
X_SpinGC_CisAis
int X_SpinGC_CisAis(long unsigned int j, struct BindStruct *X, long unsigned int is1_spin, long unsigned int sigma1)
Compute the grandcanonical spin state with bit mask is1_spin.
Definition: mltplySpinCore.c:183
BindStruct::Large
struct LargeList Large
Variables for Matrix-Vector product.
Definition: struct.h:412
X_Spin_CisAis
int X_Spin_CisAis(long unsigned int j, struct BindStruct *X, long unsigned int is1_spin, long unsigned int sigma1)
Compute the spin state with bit mask is1_spin.
Definition: mltplySpinCore.c:166
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
DefineList::iFlgGeneralSpin
int iFlgGeneralSpin
Flag for the general (Sz/=1/2) spin.
Definition: struct.h:86
FALSE
#define FALSE
Definition: global.h:25
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
DefineList::SiteToBit
long int * SiteToBit
[DefineList::NsiteMPI] Similar to DefineList::Tpow. For general spin.
Definition: struct.h:94
BitCheckGeneral
int BitCheckGeneral(const long unsigned int org_bit, const unsigned int org_isite, const unsigned int target_ispin, const long int *SiteToBit, const long unsigned int *Tpow)
bit check function for general spin
Definition: bitcalc.c:393
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetDiagonalTETransfer()

int SetDiagonalTETransfer ( long unsigned int  isite1,
double  dtmp_V,
long unsigned int  spin,
struct BindStruct *  X,
double complex *  tmp_v0,
double complex *  tmp_v1 
)

Update the vector by the general one-body diagonal interaction, \( \mu_{i\sigma_1} n_ {i\sigma_1}\).
(Using in Time Evolution mode).

Parameters
isite1[in] a site number \(i \)
dtmp_V[in] A value of general one-body diagonal interaction \( \mu_{i\sigma_1} \)
spin[in] a spin index at \(i \) site.
X[in] Define list to get the operator information.
tmp_v0[in,out] Result vector
tmp_v1[in] Input produced vector
Return values
-1fail to calculate the diagonal component.
0succeed to calculate the diagonal component.
Version
2.1
Author
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 1940 of file diagonalcalc.c.

References BitCheckGeneral(), cErrNoModel, BindStruct::Check, BindStruct::Def, FALSE, DefineList::iCalcModel, CheckList::idim_max, DefineList::iFlgGeneralSpin, BindStruct::Large, list_1, myrank, DefineList::Nsite, LargeList::prdct, DefineList::SiteToBit, stdoutMPI, SumMPI_dc(), and DefineList::Tpow.

Referenced by diagonalcalcForTE(), and SetDiagonalTEChemi().

1947  {
1948  long unsigned int is1_up;
1949  long unsigned int ibit1_up;
1950  long unsigned int num1;
1951  long unsigned int isigma1 =spin;
1952  long unsigned int is1,ibit1;
1953  double dam_pr=0.0;
1954 
1955  long unsigned int j;
1956  long unsigned int i_max=X->Check.idim_max;
1957 
1958  /*
1959  When isite1 is in the inter process region
1960  */
1961  if (isite1 > X->Def.Nsite){
1962 
1963  switch (X->Def.iCalcModel) {
1964 
1965  case HubbardGC:
1966  case KondoGC:
1967  case Hubbard:
1968  case Kondo:
1969  if (spin == 0) {
1970  is1 = X->Def.Tpow[2 * isite1 - 2];
1971  }
1972  else {
1973  is1 = X->Def.Tpow[2 * isite1 - 1];
1974  }
1975  ibit1 = (unsigned long int)myrank & is1;
1976  num1 = ibit1 / is1;
1977  break;/*case HubbardGC, case KondoGC, Hubbard, Kondo:*/
1978 
1979  case SpinGC:
1980  case Spin:
1981  if (X->Def.iFlgGeneralSpin == FALSE) {
1982  is1_up = X->Def.Tpow[isite1 - 1];
1983  num1 = (((unsigned long int)myrank& is1_up) / is1_up) ^ (1 - spin);
1984  } /*if (X->Def.iFlgGeneralSpin == FALSE)*/
1985  else /*if (X->Def.iFlgGeneralSpin == TRUE)*/ {
1986  num1 = BitCheckGeneral((unsigned long int)myrank,
1987  isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
1988  }/*if (X->Def.iFlgGeneralSpin == TRUE)*/
1989  break;/*case SpinGC, Spin:*/
1990 
1991  default:
1992  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
1993  return -1;
1994 
1995  } /*switch (X->Def.iCalcModel)*/
1996 
1997  if(num1 !=0) {
1998 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1)\
1999  firstprivate(i_max, dtmp_V) private(j)
2000  for (j = 1; j <= i_max; j++) {
2001  tmp_v0[j] += dtmp_V * tmp_v1[j];
2002  dam_pr += dtmp_V * conj(tmp_v1[j]) * tmp_v1[j];
2003  }
2004  }
2005  }/*if (isite1 >= X->Def.Nsite*/
2006  else {//(isite1 < X->Def.Nsite)
2007  switch (X->Def.iCalcModel) {
2008  case HubbardGC:
2009  if (spin == 0) {
2010  is1 = X->Def.Tpow[2 * isite1 - 2];
2011  } else {
2012  is1 = X->Def.Tpow[2 * isite1 - 1];
2013  }
2014 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1) \
2015  firstprivate(i_max, dtmp_V, is1) private(num1, ibit1)
2016  for (j = 1; j <= i_max; j++) {
2017  ibit1 = (j - 1) & is1;
2018  num1 = ibit1 / is1;
2019  tmp_v0[j] += dtmp_V * num1*tmp_v1[j];
2020  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
2021  }
2022  break;
2023 
2024  case KondoGC:
2025  case Hubbard:
2026  case Kondo:
2027  if (spin == 0) {
2028  is1 = X->Def.Tpow[2 * isite1 - 2];
2029  } else {
2030  is1 = X->Def.Tpow[2 * isite1 - 1];
2031  }
2032 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1) \
2033  firstprivate(i_max, dtmp_V, is1) private(num1, ibit1)
2034  for (j = 1; j <= i_max; j++) {
2035  ibit1 = list_1[j] & is1;
2036  num1 = ibit1 / is1;
2037  tmp_v0[j] += dtmp_V * num1*tmp_v1[j];
2038  dam_pr += dtmp_V * num1*conj(tmp_v1[j]) * tmp_v1[j];
2039  }
2040  break;
2041 
2042  case SpinGC:
2043  if (X->Def.iFlgGeneralSpin == FALSE) {
2044  is1_up = X->Def.Tpow[isite1 - 1];
2045 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1) \
2046  firstprivate(i_max, dtmp_V, is1_up, spin) private(num1, ibit1_up)
2047  for (j = 1; j <= i_max; j++) {
2048  ibit1_up = (((j - 1) & is1_up) / is1_up) ^ (1 - spin);
2049  tmp_v0[j] += dtmp_V * ibit1_up*tmp_v1[j];
2050  dam_pr += dtmp_V * ibit1_up*conj(tmp_v1[j]) * tmp_v1[j];
2051  }
2052  } else {
2053 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(tmp_v0, tmp_v1) \
2054  firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
2055  for (j = 1; j <= i_max; j++) {
2056  num1 = BitCheckGeneral(j - 1, isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
2057  if (num1 != 0) {
2058  tmp_v0[j] += dtmp_V *tmp_v1[j];
2059  dam_pr += dtmp_V *conj(tmp_v1[j]) * tmp_v1[j];
2060  }
2061  }
2062  }
2063  break;
2064 
2065  case Spin:
2066  if (X->Def.iFlgGeneralSpin == FALSE) {
2067  is1_up = X->Def.Tpow[isite1 - 1];
2068 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1)\
2069  firstprivate(i_max, dtmp_V, is1_up, spin) private(num1, ibit1_up)
2070  for (j = 1; j <= i_max; j++) {
2071  ibit1_up = ((list_1[j] & is1_up) / is1_up) ^ (1 - spin);
2072  tmp_v0[j] += dtmp_V * ibit1_up * tmp_v1[j];
2073  dam_pr += dtmp_V * ibit1_up * conj(tmp_v1[j]) * tmp_v1[j];
2074  }
2075  } else {
2076 #pragma omp parallel for default(none) reduction(+:dam_pr) shared(list_1, tmp_v0, tmp_v1)\
2077  firstprivate(i_max, dtmp_V, isite1, isigma1, X) private(j, num1)
2078  for (j = 1; j <= i_max; j++) {
2079  num1 = BitCheckGeneral(list_1[j], isite1, isigma1, X->Def.SiteToBit, X->Def.Tpow);
2080  tmp_v0[j] += dtmp_V * num1 * tmp_v1[j];
2081  dam_pr += dtmp_V * num1 * conj(tmp_v1[j]) * tmp_v1[j];
2082  }
2083  }
2084  break;
2085 
2086  default:
2087  fprintf(stdoutMPI, cErrNoModel, X->Def.iCalcModel);
2088  return -1;
2089  }
2090  }
2091  dam_pr=SumMPI_dc(dam_pr);
2092  X->Large.prdct += dam_pr;
2093  return 0;
2094 }
BindStruct::Def
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
SumMPI_dc
double complex SumMPI_dc(double complex norm)
MPI wrapper function to obtain sum of Double complex across processes.
Definition: wrapperMPI.c:205
CheckList::idim_max
unsigned long int idim_max
The dimension of the Hilbert space of this process.
Definition: struct.h:303
LargeList::prdct
double complex prdct
The expectation value of the energy.
Definition: struct.h:314
BindStruct::Large
struct LargeList Large
Variables for Matrix-Vector product.
Definition: struct.h:412
DefineList::Nsite
unsigned int Nsite
Number of sites in the INTRA process region.
Definition: struct.h:56
DefineList::Tpow
long unsigned int * Tpow
[2 * DefineList::NsiteMPI] malloc in setmem_def().
Definition: struct.h:90
DefineList::iFlgGeneralSpin
int iFlgGeneralSpin
Flag for the general (Sz/=1/2) spin.
Definition: struct.h:86
FALSE
#define FALSE
Definition: global.h:25
list_1
long unsigned int * list_1
Definition: global.h:47
cErrNoModel
char * cErrNoModel
Error Message in diagonal calc.c.
Definition: ErrorMessage.c:99
DefineList::SiteToBit
long int * SiteToBit
[DefineList::NsiteMPI] Similar to DefineList::Tpow. For general spin.
Definition: struct.h:94
BitCheckGeneral
int BitCheckGeneral(const long unsigned int org_bit, const unsigned int org_isite, const unsigned int target_ispin, const long int *SiteToBit, const long unsigned int *Tpow)
bit check function for general spin
Definition: bitcalc.c:393
DefineList::iCalcModel
int iCalcModel
Switch for model. 0:Hubbard, 1:Spin, 2:Kondo, 3:HubbardGC, 4:SpinGC, 5:KondoGC, 6:HubbardNConserved.
Definition: struct.h:198
myrank
int myrank
Process ID, defined in InitializeMPI()
Definition: global.h:163
BindStruct::Check
struct CheckList Check
Size of the Hilbert space.
Definition: struct.h:411
stdoutMPI
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the call graph for this function:
+ Here is the caller graph for this function: