HΦ  3.2.0
StdFace_ModelUtil.c
Go to the documentation of this file.
1 /*
2 HPhi-mVMC-StdFace - Common input generator
3 Copyright (C) 2015 The University of Tokyo
4 
5 This program is free software: you can redistribute it and/or modify
6 it under the terms of the GNU General Public License as published by
7 the Free Software Foundation, either version 3 of the License, or
8 (at your option) any later version.
9 
10 This program is distributed in the hope that it will be useful,
11 but WITHOUT ANY WARRANTY; without even the implied warranty of
12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 GNU General Public License for more details.
14 
15 You should have received a copy of the GNU General Public License
16 along with this program. If not, see <http://www.gnu.org/licenses/>.
17 */
21 #include <stdio.h>
22 #include <stdlib.h>
23 #include <math.h>
24 #include <complex.h>
25 #include "StdFace_vals.h"
26 #include <string.h>
27 #ifdef MPI
28 #include <mpi.h>
29 #endif
30 
35 void StdFace_exit(int errorcode
36 )
37 {
38  int ierr = 0;
39  fflush(stdout);
40  fflush(stderr);
41 #ifdef MPI
42  fprintf(stdout, "\n\n ####### You DO NOT have to WORRY about the following MPI-ERROR MESSAGE. #######\n\n");
43  ierr = MPI_Abort(MPI_COMM_WORLD, errorcode);
44  ierr = MPI_Finalize();
45 #endif
46  if (ierr != 0) fprintf(stderr, "\n MPI_Finalize() = %d\n\n", ierr);
47  exit(errorcode);
48 }
55 void StdFace_trans(
56  struct StdIntList *StdI,
57  double complex trans0,
58  int isite,
59  int ispin,
60  int jsite,
61  int jspin
62 )
63 {
64  StdI->trans[StdI->ntrans] = trans0;
65  StdI->transindx[StdI->ntrans][0] = isite;
66  StdI->transindx[StdI->ntrans][1] = ispin;
67  StdI->transindx[StdI->ntrans][2] = jsite;
68  StdI->transindx[StdI->ntrans][3] = jspin;
69  StdI->ntrans = StdI->ntrans + 1;
70 }/*void StdFace_trans*/
75 void StdFace_Hopping(
76 struct StdIntList *StdI,
77  double complex trans0,
78  int isite,
79  int jsite,
80  double *dR
81 )
82 {
83  int ispin, it, ii;
84  double complex Cphase, coef;
90 #if defined(_HPhi)
91  if (strcmp(StdI->method, "timeevolution") == 0 && StdI->PumpBody == 1) {
92  for (it = 0; it < StdI->Lanczos_max; it++) {
93  Cphase = 0.0f;
94  for (ii = 0; ii < 3; ii++) Cphase += /*2.0*StdI->pi */ StdI->At[it][ii] * dR[ii];
95  coef = cos(Cphase) + I * sin(-Cphase);
96  for (ispin = 0; ispin < 2; ispin++) {
97  StdI->pump[it][StdI->npump[it]] = coef * trans0;
98  StdI->pumpindx[it][StdI->npump[it]][0] = isite;
99  StdI->pumpindx[it][StdI->npump[it]][1] = ispin;
100  StdI->pumpindx[it][StdI->npump[it]][2] = jsite;
101  StdI->pumpindx[it][StdI->npump[it]][3] = ispin;
102  StdI->npump[it] = StdI->npump[it] + 1;
103 
104  StdI->pump[it][StdI->npump[it]] = conj(coef * trans0);
105  StdI->pumpindx[it][StdI->npump[it]][0] = jsite;
106  StdI->pumpindx[it][StdI->npump[it]][1] = ispin;
107  StdI->pumpindx[it][StdI->npump[it]][2] = isite;
108  StdI->pumpindx[it][StdI->npump[it]][3] = ispin;
109  StdI->npump[it] = StdI->npump[it] + 1;
110  }/*for (ispin = 0; ispin < 2; ispin++)*/
111  }/*for (it = 0; it < StdI->Lanczos_max; it++)*/
112  }/*if (strcmp(StdI->method, "timeevolution") == 0)*/
113  else
114 #endif
115  for (ispin = 0; ispin < 2; ispin++) {
116  StdFace_trans(StdI, trans0, jsite, ispin, isite, ispin);
117  StdFace_trans(StdI, conj(trans0), isite, ispin, jsite, ispin);
118  }/*for (ispin = 0; ispin < 2; ispin++)*/
119 }/*void StdFace_Hopping*/
125 void StdFace_HubbardLocal(
126  struct StdIntList *StdI,
127  double mu0,
128  double h0,
129  double Gamma0,
130  double U0,
131  int isite
132 )
133 {
134  StdFace_trans(StdI, mu0 + 0.5 * h0, isite, 0, isite, 0);
135  StdFace_trans(StdI, mu0 - 0.5 * h0, isite, 1, isite, 1);
136  StdFace_trans(StdI, -0.5 * Gamma0, isite, 1, isite, 0);
137  StdFace_trans(StdI, -0.5 * Gamma0, isite, 0, isite, 1);
143  StdI->Cintra[StdI->NCintra] = U0;
144  StdI->CintraIndx[StdI->NCintra][0] = isite;
145  StdI->NCintra += 1;
146 }/*void StdFace_HubbardLocal*/
151 void StdFace_MagField(
152  struct StdIntList *StdI,
153  int S2,
154  double h,
155  double Gamma,
156  int isite
157 )
158 {
159  int ispin;
160  double S, Sz;
161 
162  S = (double)S2 * 0.5;
166  for (ispin = 0; ispin <= S2; ispin++){
173  Sz = (double)ispin - S;
174  StdFace_trans(StdI, -h * Sz, isite, ispin, isite, ispin);
185  if (ispin < S2) {
186  StdFace_trans(StdI, -0.5 * Gamma * sqrt(S*(S + 1.0) - Sz*(Sz + 1.0)),
187  isite, ispin + 1, isite, ispin);
188  StdFace_trans(StdI, -0.5 * Gamma * sqrt(S*(S + 1.0) - Sz*(Sz + 1.0)),
189  isite, ispin, isite, ispin + 1);
190  }/*if (ispin < S2)*/
191  }/*for (ispin = 0; ispin <= S2; ispin++)*/
192 }/*void StdFace_MagField*/
199 void StdFace_intr(
200  struct StdIntList *StdI,
201  double complex intr0,
202  int site1,
203  int spin1,
204  int site2,
205  int spin2,
206  int site3,
207  int spin3,
208  int site4,
209  int spin4
210 )
211 {
212  StdI->intr[StdI->nintr] = intr0;
213  StdI->intrindx[StdI->nintr][0] = site1; StdI->intrindx[StdI->nintr][1] = spin1;
214  StdI->intrindx[StdI->nintr][2] = site2; StdI->intrindx[StdI->nintr][3] = spin2;
215  StdI->intrindx[StdI->nintr][4] = site3; StdI->intrindx[StdI->nintr][5] = spin3;
216  StdI->intrindx[StdI->nintr][6] = site4; StdI->intrindx[StdI->nintr][7] = spin4;
217  StdI->nintr = StdI->nintr + 1;
218 }/*void StdFace_intr*/
223 void StdFace_GeneralJ(
224 struct StdIntList *StdI,
225  double J[3][3],
226  int Si2,
227  int Sj2,
228  int isite,
229  int jsite
230 )
231 {
232  int ispin, jspin, ZGeneral, ExGeneral;
233  double Si, Sj, Siz, Sjz;
234  double complex intr0;
238  ZGeneral = 1;
239  ExGeneral = 1;
240  if (Si2 == 1 || Sj2 == 1) {
247  ZGeneral = 0;
248 
249  StdI->Hund[StdI->NHund] = -0.5 * J[2][2];
250  StdI->HundIndx[StdI->NHund][0] = isite;
251  StdI->HundIndx[StdI->NHund][1] = jsite;
252  StdI->NHund += 1;
253 
254  StdI->Cinter[StdI->NCinter] = -0.25 * J[2][2];
255  StdI->CinterIndx[StdI->NCinter][0] = isite;
256  StdI->CinterIndx[StdI->NCinter][1] = jsite;
257  StdI->NCinter += 1;
258 
259  if (fabs(J[0][1]) < 0.000001 && fabs(J[1][0]) < 0.000001
260 #if defined(_mVMC)
261  && abs(J[0][0] - J[1][1]) < 0.000001
262 #endif
263  ) {
270  ExGeneral = 0;
271 
272 #if defined(_mVMC)
273  StdI->Ex[StdI->NEx] = - 0.25 * (J[0][0] + J[1][1]);
274 #else
275  if (strcmp(StdI->model, "kondo") == 0)
276  StdI->Ex[StdI->NEx] = -0.25 * (J[0][0] + J[1][1]);
277  else
278  StdI->Ex[StdI->NEx] = 0.25 * (J[0][0] + J[1][1]);
279 #endif
280  StdI->ExIndx[StdI->NEx][0] = isite;
281  StdI->ExIndx[StdI->NEx][1] = jsite;
282  StdI->NEx += 1;
283 
284  StdI->PairLift[StdI->NPairLift] = 0.25 * (J[0][0] - J[1][1]);
285  StdI->PLIndx[StdI->NPairLift][0] = isite;
286  StdI->PLIndx[StdI->NPairLift][1] = jsite;
287  StdI->NPairLift += 1;
288  }/*if (fabs(J[0][1]) < 0.000001 && fabs(J[1][0]) < 0.000001)*/
289  }
294  Si = 0.5 * (double)Si2;
295  Sj = 0.5 * (double)Sj2;
296 
297  for (ispin = 0; ispin <= Si2; ispin++) {
298  Siz = (double)ispin - Si;
299  for (jspin = 0; jspin <= Sj2; jspin++) {
300  Sjz = (double)jspin - Sj;
308  if (ZGeneral == 1) {
309  intr0 = J[2][2] * Siz * Sjz;
310  StdFace_intr(StdI, intr0,
311  isite, ispin, isite, ispin, jsite, jspin, jsite, jspin);
312  }
324  if ((ispin < Si2 && jspin < Sj2) && ExGeneral == 1) {
325  intr0 = 0.25 * (J[0][0] + J[1][1] + I*(J[0][1] - J[1][0]))
326  * sqrt(Si * (Si + 1.0) - Siz * (Siz + 1.0))
327  * sqrt(Sj * (Sj + 1.0) - Sjz * (Sjz + 1.0));
328  StdFace_intr(StdI, intr0,
329  isite, ispin + 1, isite, ispin, jsite, jspin, jsite, jspin + 1);
330  StdFace_intr(StdI, conj(intr0),
331  isite, ispin, isite, ispin + 1, jsite, jspin + 1, jsite, jspin);
332  }
344  if ((ispin < Si2 && jspin < Sj2) && ExGeneral == 1) {
345  intr0 = 0.5 * 0.5 * (J[0][0] - J[1][1] - I*(J[0][1] + J[1][0]))
346  * sqrt(Si * (Si + 1.0) - Siz * (Siz + 1.0))
347  * sqrt(Sj * (Sj + 1.0) - Sjz * (Sjz + 1.0));
348  StdFace_intr(StdI, intr0,
349  isite, ispin + 1, isite, ispin, jsite, jspin + 1, jsite, jspin);
350  StdFace_intr(StdI, conj(intr0),
351  isite, ispin, isite, ispin + 1, jsite, jspin, jsite, jspin + 1);
352  }
363  if (ispin < Si2) {
364  intr0 = 0.5 * (J[0][2] - I * J[1][2]) * sqrt(Si * (Si + 1.0) - Siz * (Siz + 1.0)) * Sjz;
365  StdFace_intr(StdI, intr0,
366  isite, ispin + 1, isite, ispin, jsite, jspin, jsite, jspin);
367  StdFace_intr(StdI, conj(intr0),
368  jsite, jspin, jsite, jspin, isite, ispin, isite, ispin + 1);
369  }/*if (ispin < Si2)*/
380  if (jspin < Sj2) {
381  intr0 = 0.5 * (J[2][0] - I * J[2][1]) * Siz * sqrt(Sj * (Sj + 1.0) - Sjz * (Sjz + 1.0));
382  StdFace_intr(StdI, intr0,
383  isite, ispin, isite, ispin, jsite, jspin + 1, jsite, jspin);
384  StdFace_intr(StdI, conj(intr0),
385  jsite, jspin, jsite, jspin + 1, isite, ispin, isite, ispin);
386  }/*if (jspin < Sj2)*/
387  }/*for (jspin = 0; jspin <= Sj2; jspin++)*/
388  }/*for (ispin = 0; ispin <= Si2; ispin++)*/
389 }/*StdFace_GeneralJ*/
396 void StdFace_Coulomb(
397 struct StdIntList *StdI,
398  double V,
399  int isite,
400  int jsite
401 )
402 {
403  StdI->Cinter[StdI->NCinter] = V;
404  StdI->CinterIndx[StdI->NCinter][0] = isite;
405  StdI->CinterIndx[StdI->NCinter][1] = jsite;
406  StdI->NCinter += 1;
407 }/*void StdFace_Coulomb*/
414 void StdFace_PrintVal_d(
415  char* valname,
416  double *val,
417  double val0
418 )
419 {
420  if (isnan(*val) == 1) {
421  *val = val0;
422  fprintf(stdout, " %15s = %-10.5f ###### DEFAULT VALUE IS USED ######\n", valname, *val);
423  }
424  else fprintf(stdout, " %15s = %-10.5f\n", valname, *val);
425 }/*void StdFace_PrintVal_d*/
432 void StdFace_PrintVal_dd(
433  char* valname,
434  double *val,
435  double val0,
436  double val1
437 )
438 {
439  if (isnan(*val) == 1) {
443  if (isnan(val0) == 1) *val = val1;
444  else *val = val0;
445  fprintf(stdout, " %15s = %-10.5f ###### DEFAULT VALUE IS USED ######\n", valname, *val);
446  }
447  else fprintf(stdout, " %15s = %-10.5f\n", valname, *val);
448 }/*void StdFace_PrintVal_dd*/
455 void StdFace_PrintVal_c(
456  char* valname,
457  double complex *val,
458  double complex val0
459 )
460 {
461  if (isnan(creal(*val)) == 1) {
462  *val = val0;
463  fprintf(stdout, " %15s = %-10.5f %-10.5f ###### DEFAULT VALUE IS USED ######\n", valname, creal(*val), cimag(*val));
464  }
465  else fprintf(stdout, " %15s = %-10.5f %-10.5f\n", valname, creal(*val), cimag(*val));
466 }/*void StdFace_PrintVal_c*/
473 void StdFace_PrintVal_i(
474  char* valname,
475  int *val,
476  int val0
477 )
478 {
479  int NaN_i = 2147483647;/*The upper limt of Int*/
480 
481  if (*val == NaN_i) {
482  *val = val0;
483  fprintf(stdout, " %15s = %-10d ###### DEFAULT VALUE IS USED ######\n", valname, *val);
484  }
485  else fprintf(stdout, " %15s = %-10d\n", valname, *val);
486 }/*void StdFace_PrintVal_i*/
492 void StdFace_NotUsed_d(
493  char* valname,
494  double val
495 )
496 {
497  if (isnan(val) == 0) {
498  fprintf(stdout, "\n Check ! %s is SPECIFIED but will NOT be USED. \n", valname);
499  fprintf(stdout, " Please COMMENT-OUT this line \n");
500  fprintf(stdout, " or check this input is REALLY APPROPRIATE for your purpose ! \n\n");
501  StdFace_exit(-1);
502  }
503 }/*void StdFace_NotUsed_d*/
509 void StdFace_NotUsed_c(
510  char* valname,
511  double complex val
512 )
513 {
514  if (isnan(creal(val)) == 0) {
515  fprintf(stdout, "\n Check ! %s is SPECIFIED but will NOT be USED. \n", valname);
516  fprintf(stdout, " Please COMMENT-OUT this line \n");
517  fprintf(stdout, " or check this input is REALLY APPROPRIATE for your purpose ! \n\n");
518  StdFace_exit(-1);
519  }
520 }/*void StdFace_NotUsed_c*/
526 void StdFace_NotUsed_J(
527  char* valname,
528  double JAll,
529  double J[3][3]
530 )
531 {
532  int i1, i2;
533  char Jname[3][3][10];
534 
535  sprintf(Jname[0][0], "%sx", valname);
536  sprintf(Jname[0][1], "%sxy", valname);
537  sprintf(Jname[0][2], "%sxz", valname);
538  sprintf(Jname[1][0], "%syx", valname);
539  sprintf(Jname[1][1], "%sy", valname);
540  sprintf(Jname[1][2], "%syz", valname);
541  sprintf(Jname[2][0], "%szx", valname);
542  sprintf(Jname[2][1], "%szy", valname);
543  sprintf(Jname[2][2], "%sz", valname);
544 
545  StdFace_NotUsed_d(valname, JAll);
546 
547  for (i1 = 0; i1 < 3; i1++) {
548  for (i2 = 0; i2 < 3; i2++) {
549  StdFace_NotUsed_d(Jname[i1][i2], J[i1][i2]);
550  }/*for (j = 0; j < 3; j++)*/
551  }/*for (i = 0; i < 3; i++)*/
552 }/*void StdFace_NotUsed_J*/
558 void StdFace_NotUsed_i(
559  char* valname,
560  int val
561 )
562 {
563  int NaN_i = 2147483647;
564 
565  if (val != NaN_i) {
566  fprintf(stdout, "\n Check ! %s is SPECIFIED but will NOT be USED. \n", valname);
567  fprintf(stdout, " Please COMMENT-OUT this line \n");
568  fprintf(stdout, " or check this input is REALLY APPROPRIATE for your purpose ! \n\n");
569  StdFace_exit(-1);
570  }
571 }/*void StdFace_NotUsed_i*/
577 void StdFace_RequiredVal_i(
578  char* valname,
579  int val
580 )
581 {
582  int NaN_i = 2147483647;
583 
584  if (val == NaN_i){
585  fprintf(stdout, "ERROR ! %s is NOT specified !\n", valname);
586  StdFace_exit(-1);
587  }
588  else fprintf(stdout, " %15s = %-3d\n", valname, val);
589 }/*void StdFace_RequiredVal_i*/
595 static void StdFace_FoldSite(
596  struct StdIntList *StdI,
597  int iCellV[3],
598  int nBox[3],
599  int iCellV_fold[3]
601 )
602 {
603  int ii, jj, iCellV_frac[3];
607  for (ii = 0; ii < 3; ii++) {
608  iCellV_frac[ii] = 0;
609  for (jj = 0; jj < 3; jj++)iCellV_frac[ii] += StdI->rbox[ii][jj] * iCellV[jj];
610  }
614  for (ii = 0; ii < 3; ii++)
615  nBox[ii] = (iCellV_frac[ii] + StdI->NCell * 1000) / StdI->NCell - 1000;
619  for (ii = 0; ii < 3; ii++)
620  iCellV_frac[ii] -= StdI->NCell*(nBox[ii]);
621 
622  for (ii = 0; ii < 3; ii++) {
623  iCellV_fold[ii] = 0;
624  for (jj = 0; jj < 3; jj++) iCellV_fold[ii] += StdI->box[jj][ii] * iCellV_frac[jj];
625  iCellV_fold[ii] = (iCellV_fold[ii] + StdI->NCell * 1000) / StdI->NCell - 1000;
626  }/*for (ii = 0; ii < 3; ii++)*/
627 }/*static void StdFace_FoldSite*/
632 void StdFace_InitSite(
633  struct StdIntList *StdI,
634  FILE *fp,
635  int dim
636 )
637 {
638  int bound[3][2], edge, ii, jj;
639  int ipos;
640  int nBox[3], iCellV_fold[3], iCellV[3];
641  double pos[4][2], xmin, xmax/*, offset[2], scale*/;
642 
643  fprintf(stdout, "\n @ Super-Lattice setting\n\n");
647  if (
648  (StdI->L != StdI->NaN_i || StdI->W != StdI->NaN_i || StdI->Height != StdI->NaN_i)
649  &&
650  (StdI->box[0][0] != StdI->NaN_i || StdI->box[0][1] != StdI->NaN_i || StdI->box[0][2] != StdI->NaN_i ||
651  StdI->box[1][0] != StdI->NaN_i || StdI->box[1][1] != StdI->NaN_i || StdI->box[1][2] != StdI->NaN_i ||
652  StdI->box[2][0] != StdI->NaN_i || StdI->box[2][1] != StdI->NaN_i || StdI->box[2][2] != StdI->NaN_i)
653  )
654  {
655  fprintf(stdout, "\nERROR ! (L, W, Height) and (a0W, ..., a2H) conflict !\n\n");
656  StdFace_exit(-1);
657  }
658  else if (StdI->L != StdI->NaN_i || StdI->W != StdI->NaN_i || StdI->Height != StdI->NaN_i)
659  {
660  StdFace_PrintVal_i("L", &StdI->L, 1);
661  StdFace_PrintVal_i("W", &StdI->W, 1);
662  StdFace_PrintVal_i("Height", &StdI->Height, 1);
663  for (ii = 0; ii < 3; ii++) for (jj = 0; jj < 3; jj++)
664  StdI->box[ii][jj] = 0;
665  StdI->box[0][0] = StdI->W;
666  StdI->box[1][1] = StdI->L;
667  StdI->box[2][2] = StdI->Height;
668  }
669  else
670  {
671  StdFace_PrintVal_i("a0W", &StdI->box[0][0], 1);
672  StdFace_PrintVal_i("a0L", &StdI->box[0][1], 0);
673  StdFace_PrintVal_i("a0H", &StdI->box[0][2], 0);
674  StdFace_PrintVal_i("a1W", &StdI->box[1][0], 0);
675  StdFace_PrintVal_i("a1L", &StdI->box[1][1], 1);
676  StdFace_PrintVal_i("a1H", &StdI->box[1][2], 0);
677  StdFace_PrintVal_i("a2W", &StdI->box[2][0], 0);
678  StdFace_PrintVal_i("a2L", &StdI->box[2][1], 0);
679  StdFace_PrintVal_i("a2H", &StdI->box[2][2], 1);
680  }
681  /*
682  Parameters for the 3D system will not used.
683  */
684  if (dim == 2) {
685  StdI->direct[0][2] = 0.0;
686  StdI->direct[1][2] = 0.0;
687  StdI->direct[2][0] = 0.0;
688  StdI->direct[2][1] = 0.0;
689  StdI->direct[2][2] = 1.0;
690  }
695  if (dim == 2) StdI->phase[2] = 0.0;
696  for (ii = 0; ii < 3; ii++) {
697  StdI->ExpPhase[ii] = cos(StdI->pi180 * StdI->phase[ii]) + I*sin(StdI->pi180 * StdI->phase[ii]);
698  if (cabs(StdI->ExpPhase[ii] + 1.0) < 0.000001) StdI->AntiPeriod[ii] = 1;
699  else StdI->AntiPeriod[ii] = 0;
700  }
704  StdI->tau = (double **)malloc(sizeof(double*) * StdI->NsiteUC);
705  for (ii = 0; ii < StdI->NsiteUC; ii++) {
706  StdI->tau[ii] = (double *)malloc(sizeof(double) * 3);
707  }
712  StdI->NCell = 0;
713  for (ii = 0; ii < 3; ii++) {
714  StdI->NCell += StdI->box[0][ii]
715  * StdI->box[1][(ii + 1) % 3]
716  * StdI->box[2][(ii + 2) % 3]
717  - StdI->box[0][ii]
718  * StdI->box[1][(ii + 2) % 3]
719  * StdI->box[2][(ii + 1) % 3];
720  }
721  printf(" Number of Cell = %d\n", abs(StdI->NCell));
722  if (StdI->NCell == 0) {
723  StdFace_exit(-1);
724  }
725 
726  for (ii = 0; ii < 3; ii++) {
727  for (jj = 0; jj < 3; jj++) {
728  StdI->rbox[ii][jj] = StdI->box[(ii + 1) % 3][(jj + 1) % 3] * StdI->box[(ii + 2) % 3][(jj + 2) % 3]
729  - StdI->box[(ii + 1) % 3][(jj + 2) % 3] * StdI->box[(ii + 2) % 3][(jj + 1) % 3];
730  }
731  }
732  if (StdI->NCell < 0) {
733  for (ii = 0; ii < 3; ii++)
734  for (jj = 0; jj < 3; jj++)
735  StdI->rbox[ii][jj] *= -1;
736  StdI->NCell *= -1;
737  }/*if (StdI->NCell < 0)*/
742  for (ii = 0; ii < 3; ii++) {
743  bound[ii][0] = 0;
744  bound[ii][1] = 0;
745  for (nBox[2] = 0; nBox[2] < 2; nBox[2]++) {
746  for (nBox[1] = 0; nBox[1] < 2; nBox[1]++) {
747  for (nBox[0] = 0; nBox[0] < 2; nBox[0]++) {
748  edge = 0;
749  for (jj = 0; jj < 3; jj++) edge += nBox[jj] * StdI->box[jj][ii];
750  if (edge < bound[ii][0]) bound[ii][0] = edge;
751  if (edge > bound[ii][1]) bound[ii][1] = edge;
752  }
753  }
754  }
755  }
759  StdI->Cell = (int **)malloc(sizeof(int*) * StdI->NCell);
760  for (ii = 0; ii < StdI->NCell; ii++) {
761  StdI->Cell[ii] = (int *)malloc(sizeof(int) * 3);
762  }/*for (ii = 0; ii < StdI->NCell; ii++)*/
763  jj = 0;
764  for (iCellV[2] = bound[2][0]; iCellV[2] <= bound[2][1]; iCellV[2]++) {
765  for (iCellV[1] = bound[1][0]; iCellV[1] <= bound[1][1]; iCellV[1]++) {
766  for (iCellV[0] = bound[0][0]; iCellV[0] <= bound[0][1]; iCellV[0]++) {
767  StdFace_FoldSite(StdI, iCellV, nBox, iCellV_fold);
768  if (nBox[0] == 0 && nBox[1] == 0 && nBox[2] == 0) {
769  for (ii = 0; ii < 3; ii++)
770  StdI->Cell[jj][ii] = iCellV[ii];
771  jj += 1;
772  }/*if (lUC == 1)*/
773  }/*for (iCellV[0] = bound[0][0]; iCellV[0] <= bound[0][1]; iCellV[0]++*/
774  }/*for (iCellV[1] = bound[1][0]; iCellV[1] <= bound[1][1]; iCellV[1]++)*/
775  }/*for (iCellV[2] = bound[2][0]; iCellV[2] <= bound[2][1]; iCellV[2]++)*/
779  if (dim == 2) {
780  pos[0][0] = 0.0;
781  pos[0][1] = 0.0;
782  pos[1][0] = StdI->direct[0][0] * (double)StdI->box[0][0] + StdI->direct[1][0] * (double)StdI->box[0][1];
783  pos[1][1] = StdI->direct[0][1] * (double)StdI->box[0][0] + StdI->direct[1][1] * (double)StdI->box[0][1];
784  pos[2][0] = StdI->direct[0][0] * (double)StdI->box[1][0] + StdI->direct[1][0] * (double)StdI->box[1][1];
785  pos[2][1] = StdI->direct[0][1] * (double)StdI->box[1][0] + StdI->direct[1][1] * (double)StdI->box[1][1];
786  pos[3][0] = pos[1][0] + pos[2][0];
787  pos[3][1] = pos[1][1] + pos[2][1];
788 
789  xmin = 0.0;
790  xmax = 0.0;
791  for (ipos = 0; ipos < 4; ipos++) {
792  if (pos[ipos][0] < xmin) xmin = pos[ipos][0];
793  if (pos[ipos][0] > xmax) xmax = pos[ipos][0];
794  if (pos[ipos][1] < xmin) xmin = pos[ipos][1];
795  if (pos[ipos][1] > xmax) xmax = pos[ipos][1];
796  }
797  xmin -= 2.0;
798  xmax += 2.0;
799 
800  fprintf(fp, "#set terminal pdf color enhanced \\\n");
801  fprintf(fp, "#dashed dl 1.0 size 20.0cm, 20.0cm \n");
802  fprintf(fp, "#set output \"lattice.pdf\"\n");
803  fprintf(fp, "set xrange [%f: %f]\n", xmin, xmax);
804  fprintf(fp, "set yrange [%f: %f]\n", xmin, xmax);
805  fprintf(fp, "set size square\n");
806  fprintf(fp, "unset key\n");
807  fprintf(fp, "unset tics\n");
808  fprintf(fp, "unset border\n");
809 
810  fprintf(fp, "set style line 1 lc 1 lt 1\n");
811  fprintf(fp, "set style line 2 lc 5 lt 1\n");
812  fprintf(fp, "set style line 3 lc 0 lt 1\n");
813 
814  fprintf(fp, "set arrow from %f, %f to %f, %f nohead front ls 3\n", pos[0][0], pos[0][1], pos[1][0], pos[1][1]);
815  fprintf(fp, "set arrow from %f, %f to %f, %f nohead front ls 3\n", pos[1][0], pos[1][1], pos[3][0], pos[3][1]);
816  fprintf(fp, "set arrow from %f, %f to %f, %f nohead front ls 3\n", pos[3][0], pos[3][1], pos[2][0], pos[2][1]);
817  fprintf(fp, "set arrow from %f, %f to %f, %f nohead front ls 3\n", pos[2][0], pos[2][1], pos[0][0], pos[0][1]);
818  }/*if (dim == 2)*/
819 }/*void StdFace_InitSite2D*/
823 void StdFace_FindSite(
824  struct StdIntList *StdI,
825  int iW,
826  int iL,
827  int iH,
828  int diW,
829  int diL,
830  int diH,
831  int isiteUC,
832  int jsiteUC,
833  int *isite,
834  int *jsite,
835  double complex *Cphase,
836  double *dR
837 )
838 {
839  int iCell, jCell, kCell, ii;
840  int nBox[3], jCellV[3];
841 
842  dR[0] = - (double)diW + StdI->tau[isiteUC][0] - StdI->tau[jsiteUC][0];
843  dR[1] = - (double)diL + StdI->tau[isiteUC][1] - StdI->tau[jsiteUC][1];
844  dR[2] = - (double)diH + StdI->tau[isiteUC][2] - StdI->tau[jsiteUC][2];
845 
846  jCellV[0] = iW + diW;
847  jCellV[1] = iL + diL;
848  jCellV[2] = iH + diH;
849  StdFace_FoldSite(StdI, jCellV, nBox, jCellV);
850  *Cphase = 1.0;
851  for (ii = 0; ii < 3; ii++) *Cphase *= cpow(StdI->ExpPhase[ii], (double)nBox[ii]);
852 
853  for (kCell = 0; kCell < StdI->NCell; kCell++) {
854  if (jCellV[0] == StdI->Cell[kCell][0] &&
855  jCellV[1] == StdI->Cell[kCell][1] &&
856  jCellV[2] == StdI->Cell[kCell][2])
857  {
858  jCell = kCell;
859  }
860  if (iW == StdI->Cell[kCell][0] &&
861  iL == StdI->Cell[kCell][1] &&
862  iH == StdI->Cell[kCell][2])
863  {
864  iCell = kCell;
865  }
866  }/*for (iCell = 0; iCell < StdI->NCell; iCell++)*/
867  *isite = iCell * StdI->NsiteUC + isiteUC;
868  *jsite = jCell * StdI->NsiteUC + jsiteUC;
869  if (strcmp(StdI->model, "kondo") == 0) {
870  *isite += StdI->NCell * StdI->NsiteUC;
871  *jsite += StdI->NCell * StdI->NsiteUC;
872  }
873 }/*void StdFace_FindSite*/
877 void StdFace_SetLabel(
878  struct StdIntList *StdI,
879  FILE *fp,
880  int iW,
881  int iL,
882  int diW,
883  int diL,
884  int isiteUC,
885  int jsiteUC,
886  int *isite,
887  int *jsite,
888  int connect,
889  double complex *Cphase,
890  double *dR
891 )
892 {
893  double xi, yi, xj, yj;
897  StdFace_FindSite(StdI, iW, iL, 0, -diW, -diL, 0, jsiteUC, isiteUC, isite, jsite, Cphase, dR);
898 
899  xi = StdI->direct[0][0] * ((double)iW + StdI->tau[jsiteUC][0])
900  + StdI->direct[1][0] * ((double)iL + StdI->tau[jsiteUC][1]);
901  yi = StdI->direct[0][1] * ((double)iW + StdI->tau[jsiteUC][0])
902  + StdI->direct[1][1] * ((double)iL + StdI->tau[jsiteUC][1]);
903 
904  xj = StdI->direct[0][0] * ((double)(iW - diW) + StdI->tau[isiteUC][0])
905  + StdI->direct[1][0] * ((double)(iL - diL) + StdI->tau[isiteUC][1]);
906  yj = StdI->direct[0][1] * ((double)(iW - diW) + StdI->tau[isiteUC][0])
907  + StdI->direct[1][1] * ((double)(iL - diL) + StdI->tau[isiteUC][1]);
908 
909  if (*isite < 10)fprintf(fp, "set label \"%1d\" at %f, %f center front\n", *isite, xi, yi);
910  else fprintf(fp, "set label \"%2d\" at %f, %f center front\n", *isite, xi, yi);
911  if (*jsite < 10)fprintf(fp, "set label \"%1d\" at %f, %f center front\n", *jsite, xj, yj);
912  else fprintf(fp, "set label \"%2d\" at %f, %f center front\n", *jsite, xj, yj);
913  if (connect < 3)
914  fprintf(fp, "set arrow from %f, %f to %f, %f nohead ls %d\n", xi, yi, xj, yj, connect);
918  StdFace_FindSite(StdI, iW, iL, 0, diW, diL, 0, isiteUC, jsiteUC, isite, jsite, Cphase, dR);
919 
920  xi = StdI->direct[1][0] * ((double)iL + StdI->tau[isiteUC][1])
921  + StdI->direct[0][0] * ((double)iW + StdI->tau[isiteUC][0]);
922  yi = StdI->direct[1][1] * ((double)iL + StdI->tau[isiteUC][1])
923  + StdI->direct[0][1] * ((double)iW + StdI->tau[isiteUC][0]);
924 
925  xj = StdI->direct[0][0] * ((double)(iW + diW) + StdI->tau[jsiteUC][0])
926  + StdI->direct[1][0] * ((double)(iL + diL) + StdI->tau[jsiteUC][1]);
927  yj = StdI->direct[0][1] * ((double)(iW + diW) + StdI->tau[jsiteUC][0])
928  + StdI->direct[1][1] * ((double)(iL + diL) + StdI->tau[jsiteUC][1]);
929 
930  if (*isite < 10)fprintf(fp, "set label \"%1d\" at %f, %f center front\n", *isite, xi, yi);
931  else fprintf(fp, "set label \"%2d\" at %f, %f center front\n", *isite, xi, yi);
932  if (*jsite < 10)fprintf(fp, "set label \"%1d\" at %f, %f center front\n", *jsite, xj, yj);
933  else fprintf(fp, "set label \"%2d\" at %f, %f center front\n", *jsite, xj, yj);
934  if (connect < 3)
935  fprintf(fp, "set arrow from %f, %f to %f, %f nohead ls %d\n", xi, yi, xj, yj, connect);
936 }/*void StdFace_SetLabel*/
940 void StdFace_PrintXSF(struct StdIntList *StdI) {
941  FILE *fp;
942  int ii, jj, kk, isite, iCell;
943  double vec[3];
944 
945  fp = fopen("lattice.xsf", "w");
946  fprintf(fp, "CRYSTAL\n");
947  fprintf(fp, "PRIMVEC\n");
948  for (ii = 0; ii < 3; ii++) {
949  for (jj = 0; jj < 3; jj++) {
950  vec[jj] = 0.0;
951  for (kk = 0; kk < 3; kk++)
952  vec[jj] += (double)StdI->box[ii][kk] * StdI->direct[kk][jj];
953  }
954  fprintf(fp, "%15.5f %15.5f %15.5f\n", vec[0], vec[1], vec[2]);
955  }
956  fprintf(fp, "PRIMCOORD\n");
957  fprintf(fp, "%d 1\n", StdI->NCell * StdI->NsiteUC);
958  for (iCell = 0; iCell < StdI->NCell; iCell++) {
959  for (isite = 0; isite < StdI->NsiteUC; isite++) {
960  for (jj = 0; jj < 3; jj++) {
961  vec[jj] = 0.0;
962  for (kk = 0; kk < 3; kk++)
963  vec[jj] += ((double)StdI->Cell[iCell][kk] + StdI->tau[isite][kk])
964  * StdI->direct[kk][jj];
965  }
966  fprintf(fp, "H %15.5f %15.5f %15.5f\n", vec[0], vec[1], vec[2]);
967  }
968  }
969  fflush(fp);
970  fclose(fp);
971 }/*void StdFace_PrintXSF*/
975 void StdFace_InputSpinNN(
976  double J[3][3],
977  double JAll,
978  double J0[3][3],
979  double J0All,
980  char *J0name
981 )
982 {
983  int i1, i2, i3, i4;
984  char Jname[3][3][10];
985 
986  strcpy(Jname[0][0], "x\0");
987  strcpy(Jname[0][1], "xy\0");
988  strcpy(Jname[0][2], "xz\0");
989  strcpy(Jname[1][0], "yx\0");
990  strcpy(Jname[1][1], "y\0");
991  strcpy(Jname[1][2], "yz\0");
992  strcpy(Jname[2][0], "zx\0");
993  strcpy(Jname[2][1], "zy\0");
994  strcpy(Jname[2][2], "z\0");
995 
996  if (isnan(JAll) == 0 && isnan(J0All) == 0) {
997  fprintf(stdout, "\n ERROR! %s conflict !\n\n", J0name);
998  StdFace_exit(-1);
999  }
1000  for (i1 = 0; i1 < 3; i1++) {
1001  for (i2 = 0; i2 < 3; i2++) {
1002  if (isnan(JAll) == 0 && isnan(J[i1][i2]) == 0) {
1003  fprintf(stdout, "\n ERROR! J%s conflict !\n\n", Jname[i1][i2]);
1004  StdFace_exit(-1);
1005  }
1006  else if (isnan(J0All) == 0 && isnan(J[i1][i2]) == 0) {
1007  fprintf(stdout, "\n ERROR! %s and J%s conflict !\n\n",
1008  J0name, Jname[i1][i2]);
1009  StdFace_exit(-1);
1010  }
1011  else if (isnan(J0All) == 0 && isnan(J0[i1][i2]) == 0) {
1012  fprintf(stdout, "\n ERROR! %s and %s%s conflict !\n\n", J0name,
1013  J0name, Jname[i1][i2]);
1014  StdFace_exit(-1);
1015  }
1016  else if (isnan(J0[i1][i2]) == 0 && isnan(JAll) == 0) {
1017  fprintf(stdout, "\n ERROR! %s%s conflict !\n\n",
1018  J0name, Jname[i1][i2]);
1019  StdFace_exit(-1);
1020  }
1021  }/*for (j = 0; j < 3; j++)*/
1022  }/*for (i = 0; i < 3; i++)*/
1023 
1024  for (i1 = 0; i1 < 3; i1++) {
1025  for (i2 = 0; i2 < 3; i2++) {
1026  for (i3 = 0; i3 < 3; i3++) {
1027  for (i4 = 0; i4 < 3; i4++) {
1028  if (isnan(J0[i1][i2]) == 0 && isnan(J[i3][i4]) == 0) {
1029  fprintf(stdout, "\n ERROR! %s%s and J%s conflict !\n\n",
1030  J0name, Jname[i1][i2], Jname[i3][i4]);
1031  StdFace_exit(-1);
1032  }
1033  }/*for (i4 = 0; i4 < 3; i4++)*/
1034  }/*for (i3 = 0; i3 < 3; i3++)*/
1035  }/*for (j = 0; j < 3; j++)*/
1036  }/*for (i = 0; i < 3; i++)*/
1037 
1038  for (i1 = 0; i1 < 3; i1++) {
1039  for (i2 = 0; i2 < 3; i2++) {
1040  if (isnan(J0[i1][i2]) == 0)
1041  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1042  else if (isnan(J[i1][i2]) == 0) {
1043  J0[i1][i2] = J[i1][i2];
1044  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1045  }
1046  else if (i1 == i2 && isnan(J0All) == 0) {
1047  J0[i1][i2] = J0All;
1048  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1049  }
1050  else if (i1 == i2 && isnan(JAll) == 0) {
1051  J0[i1][i2] = JAll;
1052  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1053  }
1054  else {
1055  J0[i1][i2] = 0.0;
1056  }
1057  }/*for (i2 = 0; i2 < 3; i2++)*/
1058  }/*for (i = 0; i < 3; i++)*/
1059 }/*void StdFace_InputSpinNN*/
1063 void StdFace_InputSpin(
1064  double Jp[3][3],
1065  double JpAll,
1066  char *Jpname
1067 )
1068 {
1069  int i1, i2;
1070  char Jname[3][3][10];
1071 
1072  strcpy(Jname[0][0], "x\0");
1073  strcpy(Jname[0][1], "xy\0");
1074  strcpy(Jname[0][2], "xz\0");
1075  strcpy(Jname[1][0], "yx\0");
1076  strcpy(Jname[1][1], "y\0");
1077  strcpy(Jname[1][2], "yz\0");
1078  strcpy(Jname[2][0], "zx\0");
1079  strcpy(Jname[2][1], "zy\0");
1080  strcpy(Jname[2][2], "z\0");
1081 
1082  for (i1 = 0; i1 < 3; i1++) {
1083  for (i2 = 0; i2 < 3; i2++) {
1084  if (isnan(JpAll) == 0 && isnan(Jp[i1][i2]) == 0) {
1085  fprintf(stdout, "\n ERROR! %s and %s%s conflict !\n\n", Jpname,
1086  Jpname, Jname[i1][i2]);
1087  StdFace_exit(-1);
1088  }
1089  }/*for (j = 0; j < 3; j++)*/
1090  }/*for (i = 0; i < 3; i++)*/
1091 
1092  for (i1 = 0; i1 < 3; i1++) {
1093  for (i2 = 0; i2 < 3; i2++) {
1094  if (isnan(Jp[i1][i2]) == 0)
1095  fprintf(stdout, " %14s%s = %-10.5f\n", Jpname, Jname[i1][i2], Jp[i1][i2]);
1096  else if (i1 == i2 && isnan(JpAll) == 0) {
1097  Jp[i1][i2] = JpAll;
1098  fprintf(stdout, " %14s%s = %-10.5f\n", Jpname, Jname[i1][i2], Jp[i1][i2]);
1099  }
1100  else {
1101  Jp[i1][i2] = 0.0;
1102  }
1103  }/*for (i2 = 0; i2 < 3; i2++)*/
1104  }/*for (i = 0; i < 3; i++)*/
1105 }/*void StdFace_InputSpin*/
1111 void StdFace_InputCoulombV(
1112  double V,
1113  double *V0,
1114  char *V0name
1115 )
1116 {
1117  if (isnan(V) == 0 && isnan(*V0) == 0) {
1118  fprintf(stdout, "\n ERROR! %s conflicts !\n\n", V0name);
1119  StdFace_exit(-1);
1120  }
1121  else if (isnan(*V0) == 0)
1122  fprintf(stdout, " %15s = %-10.5f\n", V0name, *V0);
1123  else if (isnan(V) == 0) {
1124  *V0 = V;
1125  fprintf(stdout, " %15s = %-10.5f\n", V0name, *V0);
1126  }
1127  else {
1128  *V0 = 0.0;
1129  }
1130 }/*void StdFace_InputCoulombV*/
1136 void StdFace_InputHopp(
1137  double complex t,
1138  double complex *t0,
1139  char *t0name
1140 )
1141 {
1142  if (isnan(creal(t)) == 0 && isnan(creal(*t0)) == 0) {
1143  fprintf(stdout, "\n ERROR! %s conflicts !\n\n", t0name);
1144  StdFace_exit(-1);
1145  }
1146  else if (isnan(creal(*t0)) == 0)
1147  fprintf(stdout, " %15s = %-10.5f %-10.5f\n", t0name, creal(*t0), cimag(*t0));
1148  else if (isnan(creal(t)) == 0) {
1149  *t0 = t;
1150  fprintf(stdout, " %15s = %-10.5f %-10.5f\n", t0name, creal(*t0), cimag(*t0));
1151  }
1152  else {
1153  *t0 = 0.0;
1154  }
1155 }/*void StdFace_InputHopp*/
1159 void StdFace_PrintGeometry(struct StdIntList *StdI) {
1160  FILE *fp;
1161  int isite, iCell, ii;
1162 
1163  fp = fopen("geometry.dat", "w");
1164 
1165  for (ii = 0; ii < 3; ii++)
1166  fprintf(fp, "%25.15e %25.15e %25.15e\n",
1167  StdI->direct[ii][0], StdI->direct[ii][1], StdI->direct[ii][2]);
1168  fprintf(fp, "%25.15e %25.15e %25.15e\n",
1169  StdI->phase[0], StdI->phase[1], StdI->phase[2]);
1170  for (ii = 0; ii < 3; ii++)
1171  fprintf(fp, "%d %d %d\n",
1172  StdI->box[ii][0], StdI->box[ii][1], StdI->box[ii][2]);
1173 
1174  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1175  for (isite = 0; isite < StdI->NsiteUC; isite++) {
1176  fprintf(fp, "%d %d %d %d\n",
1177  StdI->Cell[iCell][0] - StdI->Cell[0][0],
1178  StdI->Cell[iCell][1] - StdI->Cell[0][1],
1179  StdI->Cell[iCell][2] - StdI->Cell[0][2],
1180  isite);
1181  }/*for (isite = 0; isite < StdI->NsiteUC; isite++)*/
1182  }/* for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1183  if (strcmp(StdI->model, "kondo") == 0) {
1184  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1185  for (isite = 0; isite < StdI->NsiteUC; isite++) {
1186  fprintf(fp, "%d %d %d %d\n",
1187  StdI->Cell[iCell][0] - StdI->Cell[0][0],
1188  StdI->Cell[iCell][1] - StdI->Cell[0][1],
1189  StdI->Cell[iCell][2] - StdI->Cell[0][2],
1190  isite + StdI->NsiteUC);
1191  }/*for (isite = 0; isite < StdI->NsiteUC; isite++)*/
1192  }/* for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1193  }
1194  fflush(fp);
1195  fclose(fp);
1196 }/*void StdFace_PrintGeometry()*/
1200 void StdFace_MallocInteractions(
1201  struct StdIntList *StdI,
1202  int ntransMax,
1203  int nintrMax
1204 ) {
1205  int ii;
1206 #if defined(_HPhi)
1207  int it;
1208 #endif
1209 
1212  StdI->transindx = (int **)malloc(sizeof(int*) * ntransMax);
1213  StdI->trans = (double complex *)malloc(sizeof(double complex) * ntransMax);
1214  for (ii = 0; ii < ntransMax; ii++) {
1215  StdI->transindx[ii] = (int *)malloc(sizeof(int) * 4);
1216  }
1217  StdI->ntrans = 0;
1218 #if defined(_HPhi)
1219  if (strcmp(StdI->method, "timeevolution") == 0 && StdI->PumpBody == 1) {
1220  StdI->npump = (int *)malloc(sizeof(int) * StdI->Lanczos_max);
1221  StdI->pumpindx = (int ***)malloc(sizeof(int**) * StdI->Lanczos_max);
1222  StdI->pump = (double complex **)malloc(sizeof(double complex*) * StdI->Lanczos_max);
1223  for (it = 0; it < StdI->Lanczos_max; it++) {
1224  StdI->npump[it] = 0;
1225  StdI->pumpindx[it] = (int **)malloc(sizeof(int*) * ntransMax);
1226  StdI->pump[it] = (double complex *)malloc(sizeof(double complex) * ntransMax);
1227  for (ii = 0; ii < ntransMax; ii++) {
1228  StdI->pumpindx[it][ii] = (int *)malloc(sizeof(int) * 4);
1229  }
1230  }/*for (it = 0; it < StdI->Lanczos_max;)*/
1231  }/*if (strcmp(StdI->method, "timeevolution") == 0*/
1232 #endif
1233 
1236  StdI->intrindx = (int **)malloc(sizeof(int*) * nintrMax);
1237  StdI->intr = (double complex *)malloc(sizeof(double complex) * nintrMax);
1238  for (ii = 0; ii < nintrMax; ii++) {
1239  StdI->intrindx[ii] = (int *)malloc(sizeof(int) * 8);
1240  }
1241  StdI->nintr = 0;
1245  StdI->CintraIndx = (int **)malloc(sizeof(int*) * nintrMax);
1246  StdI->Cintra = (double *)malloc(sizeof(double) * nintrMax);
1247  for (ii = 0; ii < nintrMax; ii++) {
1248  StdI->CintraIndx[ii] = (int *)malloc(sizeof(int) * 1);
1249  }
1250  StdI->NCintra = 0;
1254  StdI->CinterIndx = (int **)malloc(sizeof(int*) * nintrMax);
1255  StdI->Cinter = (double *)malloc(sizeof(double) * nintrMax);
1256  for (ii = 0; ii < nintrMax; ii++) {
1257  StdI->CinterIndx[ii] = (int *)malloc(sizeof(int) * 2);
1258  }
1259  StdI->NCinter = 0;
1263  StdI->HundIndx = (int **)malloc(sizeof(int*) * nintrMax);
1264  StdI->Hund = (double *)malloc(sizeof(double) * nintrMax);
1265  for (ii = 0; ii < nintrMax; ii++) {
1266  StdI->HundIndx[ii] = (int *)malloc(sizeof(int) * 2);
1267  }
1268  StdI->NHund = 0;
1272  StdI->ExIndx = (int **)malloc(sizeof(int*) * nintrMax);
1273  StdI->Ex = (double *)malloc(sizeof(double) * nintrMax);
1274  for (ii = 0; ii < nintrMax; ii++) {
1275  StdI->ExIndx[ii] = (int *)malloc(sizeof(int) * 2);
1276  }
1277  StdI->NEx = 0;
1281  StdI->PLIndx = (int **)malloc(sizeof(int*) * nintrMax);
1282  StdI->PairLift = (double *)malloc(sizeof(double) * nintrMax);
1283  for (ii = 0; ii < nintrMax; ii++) {
1284  StdI->PLIndx[ii] = (int *)malloc(sizeof(int) * 2);
1285  }
1286  StdI->NPairLift = 0;
1290  StdI->PHIndx = (int **)malloc(sizeof(int*) * nintrMax);
1291  StdI->PairHopp = (double *)malloc(sizeof(double) * nintrMax);
1292  for (ii = 0; ii < nintrMax; ii++) {
1293  StdI->PHIndx[ii] = (int *)malloc(sizeof(int) * 2);
1294  }
1295  StdI->NPairHopp = 0;
1296 }/*void StdFace_MallocInteractions*/
1297 #if defined(_mVMC)
1298 
1302 static void StdFace_FoldSiteSub(
1303  struct StdIntList *StdI,
1304  int iCellV[3],
1305  int nBox[3],
1306  int iCellV_fold[3]
1307 )
1308 {
1309  int ii, jj, iCellV_frac[3];
1313  for (ii = 0; ii < 3; ii++) {
1314  iCellV_frac[ii] = 0;
1315  for (jj = 0; jj < 3; jj++)iCellV_frac[ii] += StdI->rboxsub[ii][jj] * iCellV[jj];
1316  }
1320  for (ii = 0; ii < 3; ii++)
1321  nBox[ii] = (iCellV_frac[ii] + StdI->NCellsub * 1000) / StdI->NCellsub - 1000;
1325  for (ii = 0; ii < 3; ii++)
1326  iCellV_frac[ii] -= StdI->NCellsub*(nBox[ii]);
1327 
1328  for (ii = 0; ii < 3; ii++) {
1329  iCellV_fold[ii] = 0;
1330  for (jj = 0; jj < 3; jj++) iCellV_fold[ii] += StdI->boxsub[jj][ii] * iCellV_frac[jj];
1331  iCellV_fold[ii] = (iCellV_fold[ii] + StdI->NCellsub * 1000) / StdI->NCellsub - 1000;
1332  }
1333 }/*static void StdFace_FoldSiteSub*/
1338 void StdFace_Proj(struct StdIntList *StdI)
1339 {
1340  FILE *fp;
1341  int jsite, iCell, jCell, kCell;
1342  int nBox[3], iCellV[3], jCellV[3], ii;
1343  int iSym;
1344  int **Sym, **Anti;
1345 
1346  Sym = (int **)malloc(sizeof(int*) * StdI->nsite);
1347  Anti = (int **)malloc(sizeof(int*) * StdI->nsite);
1348  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1349  Sym[jsite] = (int *)malloc(sizeof(int) * StdI->nsite);
1350  Anti[jsite] = (int *)malloc(sizeof(int) * StdI->nsite);
1351  }
1355  StdI->NSym = 0;
1356  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1357 
1358  StdFace_FoldSiteSub(StdI, StdI->Cell[iCell], nBox, iCellV);
1359 
1360  StdFace_FoldSite(StdI, iCellV, nBox, iCellV);
1361 
1362  if (iCellV[0] == StdI->Cell[iCell][0] &&
1363  iCellV[1] == StdI->Cell[iCell][1] &&
1364  iCellV[2] == StdI->Cell[iCell][2]) {
1368  for (jCell = 0; jCell < StdI->NCell; jCell++) {
1369 
1370  for (ii = 0; ii < 3; ii++)jCellV[ii] = StdI->Cell[jCell][ii] + iCellV[ii];
1371  StdFace_FoldSite(StdI, jCellV, nBox, jCellV);
1372 
1373  for (kCell = 0; kCell < StdI->NCell; kCell++) {
1374  if (jCellV[0] == StdI->Cell[kCell][0] &&
1375  jCellV[1] == StdI->Cell[kCell][1] &&
1376  jCellV[2] == StdI->Cell[kCell][2])
1377  {
1378 
1379  for (jsite = 0; jsite < StdI->NsiteUC; jsite++) {
1380 
1381  Sym[StdI->NSym][jCell*StdI->NsiteUC + jsite] = kCell*StdI->NsiteUC + jsite;
1382  Anti[StdI->NSym][jCell*StdI->NsiteUC + jsite]
1383  = StdI->AntiPeriod[0] * nBox[0]
1384  + StdI->AntiPeriod[1] * nBox[1]
1385  + StdI->AntiPeriod[2] * nBox[2];
1386 
1387  if (strcmp(StdI->model, "kondo") == 0) {
1388  Sym[StdI->NSym][StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite] = StdI->nsite / 2 + kCell*StdI->NsiteUC + jsite;
1389  Anti[StdI->NSym][StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite]
1390  = StdI->AntiPeriod[0] * nBox[0]
1391  + StdI->AntiPeriod[1] * nBox[1]
1392  + StdI->AntiPeriod[2] * nBox[2];
1393  }/*if (strcmp(StdI->model, "kondo") == 0)*/
1394 
1395  }/*for (jsite = 0; jsite < StdI->NsiteUC; jsite++)*/
1396 
1397  }/*if (jWfold == StdI->Cell[kCell][0] && jLfold == StdI->Cell[kCell][1])*/
1398  }/*for (kCell = 0; kCell < StdI->NCell; kCell++)*/
1399  }/*for (jCell = 0; jCell < StdI->NCell; jCell++)*/
1400  StdI->NSym += 1;
1401  }/*if (iWfold == iW && iLfold == iL)*/
1402  }/*for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1403 
1404  fp = fopen("qptransidx.def", "w");
1405  fprintf(fp, "=============================================\n");
1406  fprintf(fp, "NQPTrans %10d\n", StdI->NSym);
1407  fprintf(fp, "=============================================\n");
1408  fprintf(fp, "======== TrIdx_TrWeight_and_TrIdx_i_xi ======\n");
1409  fprintf(fp, "=============================================\n");
1410  for (iSym = 0; iSym < StdI->NSym; iSym++) {
1411  fprintf(fp, "%d %10.5f\n", iSym, 1.0);
1412  }
1413  for (iSym = 0; iSym < StdI->NSym; iSym++) {
1414  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1415  if (Anti[iSym][jsite] % 2 == 0) Anti[iSym][jsite] = 1;
1416  else Anti[iSym][jsite] = -1;
1417  if (StdI->AntiPeriod[0] == 1 || StdI->AntiPeriod[1] == 1 || StdI->AntiPeriod[2] == 1) {
1418  fprintf(fp, "%5d %5d %5d %5d\n", iSym, jsite, Sym[iSym][jsite], Anti[iSym][jsite]);
1419  }
1420  else {
1421  fprintf(fp, "%5d %5d %5d\n", iSym, jsite, Sym[iSym][jsite]);
1422  }
1423  }
1424  }
1425  fflush(fp);
1426  fclose(fp);
1427  fprintf(stdout, " qptransidx.def is written.\n");
1428 
1429  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1430  free(Anti[jsite]);
1431  free(Sym[jsite]);
1432  }
1433  free(Sym);
1434  free(Anti);
1435 }/*void StdFace_Proj(struct StdIntList *StdI)*/
1440 static void StdFace_InitSiteSub(struct StdIntList *StdI)
1441 {
1442  int ii, jj, kk, prod;
1446  if ((StdI->Lsub != StdI->NaN_i || StdI->Wsub != StdI->NaN_i || StdI->Hsub != StdI->NaN_i)
1447  && (StdI->boxsub[0][0] != StdI->NaN_i || StdI->boxsub[0][1] != StdI->NaN_i || StdI->boxsub[0][2] != StdI->NaN_i ||
1448  StdI->boxsub[1][0] != StdI->NaN_i || StdI->boxsub[1][1] != StdI->NaN_i || StdI->boxsub[1][2] != StdI->NaN_i ||
1449  StdI->boxsub[2][0] != StdI->NaN_i || StdI->boxsub[2][1] != StdI->NaN_i || StdI->boxsub[2][2] != StdI->NaN_i))
1450  {
1451  fprintf(stdout, "\nERROR ! (Lsub, Wsub, Hsub) and (a0Wsub, ..., a2Hsub) conflict !\n\n");
1452  StdFace_exit(-1);
1453  }
1454  else if (StdI->Wsub != StdI->NaN_i || StdI->Lsub != StdI->NaN_i || StdI->Hsub != StdI->NaN_i) {
1455  StdFace_PrintVal_i("Lsub", &StdI->Lsub, 1);
1456  StdFace_PrintVal_i("Wsub", &StdI->Wsub, 1);
1457  StdFace_PrintVal_i("Hsub", &StdI->Hsub, 1);
1458  for (ii = 0; ii < 3; ii++) for (jj = 0; jj < 3; jj++)
1459  StdI->boxsub[ii][jj] = 0;
1460  StdI->boxsub[0][0] = StdI->Wsub;
1461  StdI->boxsub[1][1] = StdI->Lsub;
1462  StdI->boxsub[2][2] = StdI->Hsub;
1463  }
1464  else {
1465  StdFace_PrintVal_i("a0Wsub", &StdI->boxsub[0][0], StdI->box[0][0]);
1466  StdFace_PrintVal_i("a0Lsub", &StdI->boxsub[0][1], StdI->box[0][1]);
1467  StdFace_PrintVal_i("a0Hsub", &StdI->boxsub[0][2], StdI->box[0][2]);
1468  StdFace_PrintVal_i("a1Wsub", &StdI->boxsub[1][0], StdI->box[1][0]);
1469  StdFace_PrintVal_i("a1Lsub", &StdI->boxsub[1][1], StdI->box[1][1]);
1470  StdFace_PrintVal_i("a1Hsub", &StdI->boxsub[1][2], StdI->box[1][2]);
1471  StdFace_PrintVal_i("a2Wsub", &StdI->boxsub[2][0], StdI->box[2][0]);
1472  StdFace_PrintVal_i("a2Lsub", &StdI->boxsub[2][1], StdI->box[2][1]);
1473  StdFace_PrintVal_i("a2Hsub", &StdI->boxsub[2][2], StdI->box[2][2]);
1474  }
1478  StdI->NCellsub = 0;
1479  for (ii = 0; ii < 3; ii++) {
1480  StdI->NCellsub += StdI->boxsub[0][ii]
1481  * StdI->boxsub[1][(ii + 1) % 3]
1482  * StdI->boxsub[2][(ii + 2) % 3]
1483  - StdI->boxsub[0][ii]
1484  * StdI->boxsub[1][(ii + 2) % 3]
1485  * StdI->boxsub[2][(ii + 1) % 3];
1486  }
1487  printf(" Number of Cell in the sublattice: %d\n", abs(StdI->NCellsub));
1488  if (StdI->NCellsub == 0) {
1489  StdFace_exit(-1);
1490  }
1491 
1492  for (ii = 0; ii < 3; ii++) {
1493  for (jj = 0; jj < 3; jj++) {
1494  StdI->rboxsub[ii][jj] = StdI->boxsub[(ii + 1) % 3][(jj + 1) % 3] * StdI->boxsub[(ii + 2) % 3][(jj + 2) % 3]
1495  - StdI->boxsub[(ii + 1) % 3][(jj + 2) % 3] * StdI->boxsub[(ii + 2) % 3][(jj + 1) % 3];
1496  }
1497  }
1498  if (StdI->NCellsub < 0) {
1499  for (ii = 0; ii < 3; ii++)
1500  for (jj = 0; jj < 3; jj++)
1501  StdI->rboxsub[ii][jj] *= -1;
1502  StdI->NCellsub *= -1;
1503  }/*if (StdI->NCell < 0)*/
1507  for (ii = 0; ii < 3; ii++) {
1508  for (jj = 0; jj < 3; jj++) {
1509  prod = 0.0;
1510  for (kk = 0; kk < 3; kk++) prod += StdI->rboxsub[ii][kk] * (double)StdI->box[jj][kk];
1511  if (prod % StdI->NCellsub != 0) {
1512  printf("\n ERROR ! Sublattice is INCOMMENSURATE !\n\n");
1513  StdFace_exit(-1);
1514  }/*if (prod % StdI->NCellsub != 0)*/
1515  }
1516  }
1517 }/*void StdFace_InitSiteSub*/
1521 void StdFace_generate_orb(struct StdIntList *StdI) {
1522  int iCell, jCell, kCell, iCell2, jCell2, iOrb, isite, jsite, Anti;
1523  int nBox[3], iCellV[3], jCellV[3], dCellV[3], ii;
1524  int **CellDone;
1525 
1526  StdFace_InitSiteSub(StdI);
1527 
1528  StdI->Orb = (int **)malloc(sizeof(int*) * StdI->nsite);
1529  StdI->AntiOrb = (int **)malloc(sizeof(int*) * StdI->nsite);
1530  for (isite = 0; isite < StdI->nsite; isite++) {
1531  StdI->Orb[isite] = (int *)malloc(sizeof(int) * StdI->nsite);
1532  StdI->AntiOrb[isite] = (int *)malloc(sizeof(int) * StdI->nsite);
1533  }
1534  CellDone = (int **)malloc(sizeof(int*) * StdI->NCell);
1535  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1536  CellDone[iCell] = (int *)malloc(sizeof(int) * StdI->NCell);
1537  for (jCell = 0; jCell < StdI->NCell; jCell++) {
1538  CellDone[iCell][jCell] = 0;
1539  }
1540  }
1541 
1542  iOrb = 0;
1543  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1544 
1545  StdFace_FoldSiteSub(StdI, StdI->Cell[iCell], nBox, iCellV);
1546 
1547  StdFace_FoldSite(StdI, iCellV, nBox, iCellV);
1548 
1549  for (kCell = 0; kCell < StdI->NCell; kCell++) {
1550  if (iCellV[0] == StdI->Cell[kCell][0] &&
1551  iCellV[1] == StdI->Cell[kCell][1] &&
1552  iCellV[2] == StdI->Cell[kCell][2])
1553  {
1554  iCell2 = kCell;
1555  }
1556  }/*for (kCell = 0; kCell < StdI->NCell; kCell++)*/
1557 
1558  for (jCell = 0; jCell < StdI->NCell; jCell++) {
1559 
1560  for (ii = 0; ii < 3; ii++)
1561  jCellV[ii] = StdI->Cell[jCell][ii] + iCellV[ii] - StdI->Cell[iCell][ii];
1562 
1563  StdFace_FoldSite(StdI, jCellV, nBox, jCellV);
1564 
1565  for (kCell = 0; kCell < StdI->NCell; kCell++) {
1566  if (jCellV[0] == StdI->Cell[kCell][0] &&
1567  jCellV[1] == StdI->Cell[kCell][1] &&
1568  jCellV[2] == StdI->Cell[kCell][2])
1569  {
1570  jCell2 = kCell;
1571  }
1572  }/*for (kCell = 0; kCell < StdI->NCell; kCell++)*/
1573  /*
1574  AntiPeriodic factor
1575  */
1576  for (ii = 0; ii < 3; ii++)
1577  dCellV[ii] = StdI->Cell[jCell][ii] - StdI->Cell[iCell][ii];
1578  StdFace_FoldSite(StdI, dCellV, nBox, dCellV);
1579  Anti = 0;
1580  for (ii = 0; ii < 3; ii++)Anti += StdI->AntiPeriod[ii] * nBox[ii];
1581  if (Anti % 2 == 0) Anti = 1;
1582  else Anti = -1;
1583 
1584  for (isite = 0; isite < StdI->NsiteUC; isite++) {
1585  for (jsite = 0; jsite < StdI->NsiteUC; jsite++) {
1586 
1587  if (CellDone[iCell2][jCell2] == 0) {
1588  StdI->Orb[iCell2*StdI->NsiteUC + isite][jCell2*StdI->NsiteUC + jsite] = iOrb;
1589  StdI->AntiOrb[iCell2*StdI->NsiteUC + isite][jCell2*StdI->NsiteUC + jsite] = Anti;
1590  iOrb += 1;
1591  }
1592  StdI->Orb[iCell*StdI->NsiteUC + isite][jCell*StdI->NsiteUC + jsite]
1593  = StdI->Orb[iCell2*StdI->NsiteUC + isite][jCell2*StdI->NsiteUC + jsite];
1594  StdI->AntiOrb[iCell*StdI->NsiteUC + isite][jCell*StdI->NsiteUC + jsite] = Anti;
1595 
1596  if (strcmp(StdI->model, "kondo") == 0) {
1597  if (CellDone[iCell2][jCell2] == 0) {
1598  StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1599  [ jCell2*StdI->NsiteUC + jsite] = iOrb;
1600  StdI->AntiOrb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1601  [ jCell2*StdI->NsiteUC + jsite] = Anti;
1602  iOrb += 1;
1603  StdI->Orb[ iCell2*StdI->NsiteUC + isite]
1604  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = iOrb;
1605  StdI->AntiOrb[ iCell2*StdI->NsiteUC + isite]
1606  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = Anti;
1607  iOrb += 1;
1608  StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1609  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = iOrb;
1610  StdI->AntiOrb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1611  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = Anti;
1612  iOrb += 1;
1613  }
1614  StdI->Orb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1615  [ jCell*StdI->NsiteUC + jsite]
1616  = StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1617  [ jCell2*StdI->NsiteUC + jsite];
1618  StdI->AntiOrb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1619  [ jCell*StdI->NsiteUC + jsite] = Anti;
1620  StdI->Orb[ iCell*StdI->NsiteUC + isite]
1621  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite]
1622  = StdI->Orb[ iCell2*StdI->NsiteUC + isite]
1623  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite];
1624  StdI->AntiOrb[iCell*StdI->NsiteUC + isite]
1625  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite] = Anti;
1626  StdI->Orb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1627  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite]
1628  = StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1629  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite];
1630  StdI->AntiOrb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1631  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite] = Anti;
1632  }/*if (strcmp(StdI->model, "kondo") == 0)*/
1633 
1634  }/*for (jsite = 0; jsite < StdI->NsiteUC; jsite++)*/
1635  }/*for (isite = 0; isite < StdI->NsiteUC; isite++)*/
1636  CellDone[iCell2][jCell2] = 1;
1637  }/*for (jCell = 0; jCell < StdI->NCell; jCell++)*/
1638 
1639  }/*for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1640  StdI->NOrb = iOrb;
1641 
1642  for (iCell = 0; iCell < StdI->NCell; iCell++) free(CellDone[iCell]);
1643  free(CellDone);
1644 }/*void StdFace_generate_orb*/
1649 void PrintJastrow(struct StdIntList *StdI) {
1650  FILE *fp;
1651  int isite, jsite, isiteUC, jsiteUC, revarsal, isite1, jsite1, iorb;
1652  int NJastrow, iJastrow;
1653  int dCell, iCell;//, jCell, dCellv[3];
1654  int **Jastrow;
1655  double complex Cphase;
1656  double dR[3];
1657 
1658  Jastrow = (int **)malloc(sizeof(int*) * StdI->nsite);
1659  for (isite = 0; isite < StdI->nsite; isite++)
1660  Jastrow[isite] = (int *)malloc(sizeof(int) * StdI->nsite);
1661 
1662  if (abs(StdI->NMPTrans) == 1 || StdI->NMPTrans == StdI->NaN_i) {
1666  for (isite = 0; isite < StdI->nsite; isite++) {
1667  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1668  Jastrow[isite][jsite] = StdI->Orb[isite][jsite];
1669  }/*for (jsite = 0; jsite < isite; jsite++)*/
1670  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1674  for (iorb = 0; iorb < StdI->NOrb; iorb++) {
1675  for (isite = 0; isite < StdI->nsite; isite++) {
1676  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1677  if (Jastrow[isite][jsite] == iorb) {
1678  Jastrow[jsite][isite] = Jastrow[isite][jsite];
1679  }
1680  }/*for (jsite = 0; jsite < isite; jsite++)*/
1681  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1682  }/*for (iorb = 0; iorb < StdI->NOrb; iorb++)*/
1683 
1684  if (strcmp(StdI->model, "hubbard") == 0) NJastrow = 0;
1685  else NJastrow = -1;
1686  for (isite = 0; isite < StdI->nsite; isite++) {
1687  /*
1688  For Local spin
1689  */
1690  if (StdI->locspinflag[isite] != 0) {
1691  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1692  Jastrow[isite][jsite] = -1;
1693  Jastrow[jsite][isite] = -1;
1694  }
1695  continue;
1696  }
1697 
1698  for (jsite = 0; jsite < isite; jsite++) {
1699  if (Jastrow[isite][jsite] >= 0) {
1700  iJastrow = Jastrow[isite][jsite];
1701  NJastrow -= 1;
1702  for (isite1 = 0; isite1 < StdI->nsite; isite1++) {
1703  for (jsite1 = 0; jsite1 < StdI->nsite; jsite1++) {
1704  if (Jastrow[isite1][jsite1] == iJastrow)
1705  Jastrow[isite1][jsite1] = NJastrow;
1706  }/*for (jsite1 = 0; jsite1 < StdI->nsite; jsite1++)*/
1707  }/*for (isite1 = 0; isite1 < StdI->nsite; isite1++)*/
1708  }/*if (Jastrow[isite][jsite] >= 0)*/
1709  }/*for (jsite = 0; jsite < isite; jsite++)*/
1710  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1711 
1712  NJastrow = -NJastrow;
1713  for (isite = 0; isite < StdI->nsite; isite++) {
1714  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1715  Jastrow[isite][jsite] = -1 - Jastrow[isite][jsite];
1716  }/*for (jsite = 0; jsite < isite; jsite++)*/
1717  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1718  }/*if (abs(StdI->NMPTrans) == 1)*/
1719  else {
1720 
1721  if (strcmp(StdI->model, "spin") == 0) {
1722  NJastrow = 1;
1723 
1724  for (isite = 0; isite < StdI->nsite; isite++) {
1725  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1726  Jastrow[isite][jsite] = 0;
1727  }/*for (jsite = 0; jsite < StdI->nsite; jsite++)*/
1728  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1729  }/*if (strcmp(StdI->model, "spin") == 0)*/
1730  else {
1731 
1732  NJastrow = 0;
1733 
1734  if (strcmp(StdI->model, "kondo") == 0) {
1735  /*
1736  Local spin - itererant electron part
1737  */
1738  for (isite = 0; isite < StdI->nsite; isite++) {
1739  for (jsite = 0; jsite < StdI->nsite / 2; jsite++) {
1740  Jastrow[isite][jsite] = 0;
1741  Jastrow[jsite][isite] = 0;
1742  }/*for (jsite = 0; jsite < StdI->nsite; jsite++)*/
1743  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1744 
1745  NJastrow += 1;
1746  }/*if (strcmp(StdI->model, "kondo") == 0)*/
1747 
1748  for (dCell = 0; dCell < StdI->NCell; dCell++) {
1749  StdFace_FindSite(StdI,
1750  0, 0, 0,
1751  -StdI->Cell[dCell][0], -StdI->Cell[dCell][1], -StdI->Cell[dCell][2],
1752  0, 0, &isite, &jsite, &Cphase, dR);
1753  if (strcmp(StdI->model, "kondo") == 0) jsite += -StdI->NCell * StdI->NsiteUC;
1754  iCell = jsite / StdI->NsiteUC;
1755  if (iCell < dCell) {
1756  /*
1757  If -R has been already done, skip.
1758  */
1759  continue;
1760  }
1761  else if (iCell == dCell) {
1762  /*
1763  If revarsal symmetry [Fold(-R) = R], J(R,i,j) = J(R,j,i)
1764  */
1765  revarsal = 1;
1766  }
1767  else revarsal = 0;
1768 
1769  for (isiteUC = 0; isiteUC < StdI->NsiteUC; isiteUC++) {
1770  for (jsiteUC = 0; jsiteUC < StdI->NsiteUC; jsiteUC++) {
1771  if (revarsal == 1 && jsiteUC > isiteUC) continue;/*If [Fold(-R) = R]*/
1772  if (isiteUC == jsiteUC &&
1773  StdI->Cell[dCell][0] == 0 &&
1774  StdI->Cell[dCell][1] == 0 &&
1775  StdI->Cell[dCell][2] == 0) continue;/*Diagonal*/
1776 
1777  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1778  StdFace_FindSite(StdI,
1779  StdI->Cell[iCell][0], StdI->Cell[iCell][1], StdI->Cell[iCell][2],
1780  StdI->Cell[dCell][0], StdI->Cell[dCell][1], StdI->Cell[dCell][2],
1781  isiteUC, jsiteUC, &isite, &jsite, &Cphase, dR);
1782 
1783  Jastrow[isite][jsite] = NJastrow;
1784  Jastrow[jsite][isite] = NJastrow;
1785 
1786  }/*for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1787 
1788  NJastrow += 1;
1789 
1790  }/*for (jsiteUC = 0; jsiteUC < StdI->NsiteUC; jsiteUC++)*/
1791  }/*for (isiteUC = 0; isiteUC < StdI->NsiteUC; isiteUC++)*/
1792  }/*for (dCell = 0; dCell < StdI->NCell; dCell++)*/
1793  }/*if (strcmp(StdI->model, "spin") != 0)*/
1794  }/*if (abs(StdI->NMPTrans) != 1)*/
1795 
1796  fp = fopen("jastrowidx.def", "w");
1797  fprintf(fp, "=============================================\n");
1798  fprintf(fp, "NJastrowIdx %10d\n", NJastrow);
1799  fprintf(fp, "ComplexType %10d\n", 0);
1800  fprintf(fp, "=============================================\n");
1801  fprintf(fp, "=============================================\n");
1802 
1803  for (isite = 0; isite < StdI->nsite; isite++) {
1804  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1805  if (isite == jsite) continue;
1806  fprintf(fp, "%5d %5d %5d\n", isite, jsite, Jastrow[isite][jsite]);
1807  }/*for (jsite = 0; jsite < isite; jsite++)*/
1808  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1809 
1810  for (iJastrow = 0; iJastrow < NJastrow; iJastrow++){
1811  if (strcmp(StdI->model, "hubbard") == 0 || iJastrow > 0)
1812  fprintf(fp, "%5d %5d\n", iJastrow, 1);
1813  else
1814  fprintf(fp, "%5d %5d\n", iJastrow, 0);
1815  }/*for (iJastrow = 0; iJastrow < NJastrow; iJastrow++)*/
1816  fflush(fp);
1817  fclose(fp);
1818  fprintf(stdout, " jastrowidx.def is written.\n");
1819 
1820  for (isite = 0; isite < StdI->nsite; isite++) free(Jastrow[isite]);
1821  free(Jastrow);
1822 }/*void PrintJastrow*/
1823 #endif
StdFace_trans
void StdFace_trans(struct StdIntList *StdI, double complex trans0, int isite, int ispin, int jsite, int jspin)
Add transfer to the list set StdIntList::trans and StdIntList::transindx and increment StdIntList::nt...
Definition: StdFace_ModelUtil.c:55
StdFace_PrintVal_i
void StdFace_PrintVal_i(char *valname, int *val, int val0)
Print a valiable (integer) read from the input file if it is not specified in the input file (=214748...
Definition: StdFace_ModelUtil.c:473
StdIntList::NPairHopp
int NPairHopp
Number of pair-hopping term, counted in each lattice file.
Definition: StdFace_vals.h:221
StdIntList::NHund
int NHund
Number of Hund term, counted in each lattice file.
Definition: StdFace_vals.h:197
StdIntList::Jp
double Jp[3][3]
Isotropic, diagonal/off-diagonal spin coupling (2nd Near.), input parameter J&#39;x, J&#39;y, J&#39;z, J&#39;xy, etc.
Definition: StdFace_vals.h:114
StdIntList::box
int box[3][3]
The shape of the super-cell. Input parameter a0W, a0L, a0H, etc. or defined from StdIntList::W, etc. in StdFace_InitSite().
Definition: StdFace_vals.h:44
StdIntList::ExpPhase
double complex ExpPhase[3]
.
Definition: StdFace_vals.h:155
StdIntList::J
double J[3][3]
Isotropic, diagonal/off-diagonal spin coupling (1st Near.), input parameter Jx, Jy, Jz, Jxy, etc.
Definition: StdFace_vals.h:112
StdFace_GeneralJ
void StdFace_GeneralJ(struct StdIntList *StdI, double J[3][3], int Si2, int Sj2, int isite, int jsite)
Treat J as a 3*3 matrix [(6S + 1)*(6S&#39; + 1) interactions].
Definition: StdFace_ModelUtil.c:223
StdIntList::L
int L
Number of sites along the 2nd axis, input parameter.
Definition: StdFace_vals.h:40
StdIntList::intr
double complex * intr
[StdIntList::nintr] Coefficient of general two-body term, malloc in StdFace_MallocInteractions() and ...
Definition: StdFace_vals.h:176
StdFace_HubbardLocal
void StdFace_HubbardLocal(struct StdIntList *StdI, double mu0, double h0, double Gamma0, double U0, int isite)
Add intra-Coulomb, magnetic field, chemical potential for the itenerant electron. ...
Definition: StdFace_ModelUtil.c:125
StdIntList::pi180
double pi180
, set in StdFace_ResetVals().
Definition: StdFace_vals.h:153
StdFace_PrintGeometry
void StdFace_PrintGeometry(struct StdIntList *StdI)
Print geometry of sites for the pos-process of correlation function.
Definition: StdFace_ModelUtil.c:1159
StdIntList::At
double ** At
[StdIntList::nt][3] Vector potential.
Definition: StdFace_vals.h:314
StdIntList::Ex
double * Ex
[StdIntList::NEx] Coefficient of exchange term, malloc in StdFace_MallocInteractions() and set in Std...
Definition: StdFace_vals.h:210
StdIntList::JpAll
double JpAll
Isotropic, diagonal spin coupling (2nd Near), input parameter Jp.
Definition: StdFace_vals.h:90
StdFace_InputHopp
void StdFace_InputHopp(double complex t, double complex *t0, char *t0name)
Input hopping integral from the input file, if it is not specified, use the default value(0 or the is...
Definition: StdFace_ModelUtil.c:1136
StdIntList::PLIndx
int ** PLIndx
[StdIntList::NPairLift][2] Site indices of pair-lift term, malloc in StdFace_MallocInteractions() and...
Definition: StdFace_vals.h:215
StdIntList::trans
double complex * trans
[StdIntList::ntrans] Coefficient of one-body term, malloc in StdFace_MallocInteractions() and set in ...
Definition: StdFace_vals.h:168
StdFace_Hopping
void StdFace_Hopping(struct StdIntList *StdI, double complex trans0, int isite, int jsite, double *dR)
Add Hopping for the both spin.
Definition: StdFace_ModelUtil.c:75
StdIntList::t
double complex t
Nearest-neighbor hopping, input parameter.
Definition: StdFace_vals.h:62
StdFace_MallocInteractions
void StdFace_MallocInteractions(struct StdIntList *StdI, int ntransMax, int nintrMax)
Malloc Arrays for interactions.
Definition: StdFace_ModelUtil.c:1200
StdIntList::JAll
double JAll
Isotropic, diagonal spin coupling (1st Near.), input parameter J.
Definition: StdFace_vals.h:88
StdIntList::S2
int S2
Total spin |S| of a local spin, input from file.
Definition: StdFace_vals.h:236
StdIntList::ExIndx
int ** ExIndx
[StdIntList::NEx][2] Site indices of exchange term, malloc in StdFace_MallocInteractions() and set in...
Definition: StdFace_vals.h:207
vec
double complex ** vec
Definition: global.h:45
StdIntList::NsiteUC
int NsiteUC
Number of sites in the unit cell. Defined in the beginning of each lattice function.
Definition: StdFace_vals.h:53
StdFace_InputSpin
void StdFace_InputSpin(double Jp[3][3], double JpAll, char *Jpname)
Input spin-spin interaction other than nearest-neighbor.
Definition: StdFace_ModelUtil.c:1063
StdFace_InitSite
void StdFace_InitSite(struct StdIntList *StdI, FILE *fp, int dim)
Initialize the super-cell where simulation is performed.
Definition: StdFace_ModelUtil.c:632
StdIntList::W
int W
Number of sites along the 1st axis, input parameter.
Definition: StdFace_vals.h:39
StdIntList::ntrans
int ntrans
Number of transfer, counted in each lattice file.
Definition: StdFace_vals.h:164
StdIntList::model
char model[256]
Name of model, input parameter.
Definition: StdFace_vals.h:60
StdFace_NotUsed_J
void StdFace_NotUsed_J(char *valname, double JAll, double J[3][3])
Stop HPhi if variables (real) not used is specified in the input file (!=NaN).
Definition: StdFace_ModelUtil.c:526
StdFace_PrintVal_dd
void StdFace_PrintVal_dd(char *valname, double *val, double val0, double val1)
Print a valiable (real) read from the input file if it is not specified in the input file (=NaN)...
Definition: StdFace_ModelUtil.c:432
StdIntList::AntiPeriod
int AntiPeriod[3]
If corresponding StdIntList::phase = 180, it becomes 1.
Definition: StdFace_vals.h:156
StdIntList::CintraIndx
int ** CintraIndx
[StdIntList::NCintra][1] Site indices of intra-site Coulomb term, malloc in StdFace_MallocInteraction...
Definition: StdFace_vals.h:182
StdIntList::NCinter
int NCinter
Number of inter-site Coulomb interaction, counted in each lattice file.
Definition: StdFace_vals.h:188
StdIntList::Gamma
double Gamma
Transvars magnetic field, input parameter.
Definition: StdFace_vals.h:148
StdIntList::NCintra
int NCintra
Number of intra-site Coulomb interaction, counted in each lattice file.
Definition: StdFace_vals.h:179
StdIntList::rbox
int rbox[3][3]
The inversion of StdIntList::box. Set in StdFace_InitSite().
Definition: StdFace_vals.h:47
StdIntList::J0All
double J0All
Anisotropic, diagonal spin coupling (1st Near), input parameter J0.
Definition: StdFace_vals.h:92
StdIntList::J0
double J0[3][3]
Isotropic, diagonal/off-diagonal spin coupling (1st Near.), input parameter J0x, J0y, J0z, J0xy, etc. or set in StdFace_InputSpinNN().
Definition: StdFace_vals.h:116
StdIntList::Hund
double * Hund
[StdIntList::NHund] Coefficient of Hund term, malloc in StdFace_MallocInteractions() and set in StdFa...
Definition: StdFace_vals.h:202
StdIntList::Cell
int ** Cell
[StdIntList][3] The cell position in the fractional coordinate. Malloc and Set in StdFace_InitSite()...
Definition: StdFace_vals.h:51
StdIntList::PairLift
double * PairLift
[StdIntList::NPairLift] Coefficient of pair-lift term, malloc in StdFace_MallocInteractions() and set...
Definition: StdFace_vals.h:218
StdIntList::phase
double phase[3]
Boundary phase, input parameter phase0, etc.
Definition: StdFace_vals.h:154
StdIntList::pumpindx
int *** pumpindx
[StdIntList::nt][StdIntList::npump][4] Site/spin indices of one-body term, malloc in StdFace_MallocIn...
Definition: StdFace_vals.h:308
StdIntList::NPairLift
int NPairLift
Number of pair-lift term, counted in each lattice file.
Definition: StdFace_vals.h:213
StdIntList::PumpBody
int PumpBody
one- or two-body pumping, defined from StdIntList::PumpType
Definition: StdFace_vals.h:305
StdIntList::npump
int * npump
[StdIntList::nt] Number of transfer, counted in each lattice file.
Definition: StdFace_vals.h:307
StdIntList::locspinflag
int * locspinflag
[StdIntList::nsite] LocSpin in Expert mode, malloc and set in each lattice file.
Definition: StdFace_vals.h:162
StdIntList::HundIndx
int ** HundIndx
[StdIntList::NHund][2] Site indices of Hund term, malloc in StdFace_MallocInteractions() and set in S...
Definition: StdFace_vals.h:199
StdIntList::Lanczos_max
int Lanczos_max
The maxixmum number of iterations, input from file.
Definition: StdFace_vals.h:265
StdIntList::V
double V
Off-site Coulomb potential (1st), input parameter.
Definition: StdFace_vals.h:75
StdIntList::t0
double complex t0
Anisotropic hopping (1st), input parameter.
Definition: StdFace_vals.h:64
StdIntList::direct
double direct[3][3]
The unit direct lattice vector. Set in StdFace_InitSite().
Definition: StdFace_vals.h:42
StdFace_MagField
void StdFace_MagField(struct StdIntList *StdI, int S2, double h, double Gamma, int isite)
Add longitudinal and transvars magnetic field to the list.
Definition: StdFace_ModelUtil.c:151
StdFace_NotUsed_d
void StdFace_NotUsed_d(char *valname, double val)
Stop HPhi if a variable (real) not used is specified in the input file (!=NaN).
Definition: StdFace_ModelUtil.c:492
StdIntList::V0
double V0
Anisotropic Coulomb potential (1st), input parameter.
Definition: StdFace_vals.h:77
StdIntList::pump
double complex ** pump
[StdIntList::nt][StdIntList::npump] Coefficient of one-body term, malloc in StdFace_MallocInteraction...
Definition: StdFace_vals.h:311
StdFace_Coulomb
void StdFace_Coulomb(struct StdIntList *StdI, double V, int isite, int jsite)
Add onsite/offsite Coulomb term to the list StdIntList::Cinter and StdIntList::CinterIndx, and increase the number of them (StdIntList::NCinter).
Definition: StdFace_ModelUtil.c:396
StdFace_InputCoulombV
void StdFace_InputCoulombV(double V, double *V0, char *V0name)
Input off-site Coulomb interaction from the input file, if it is not specified, use the default value...
Definition: StdFace_ModelUtil.c:1111
StdFace_SetLabel
void StdFace_SetLabel(struct StdIntList *StdI, FILE *fp, int iW, int iL, int diW, int diL, int isiteUC, int jsiteUC, int *isite, int *jsite, int connect, double complex *Cphase, double *dR)
Set Label in the gnuplot display (Only used in 2D system)
Definition: StdFace_ModelUtil.c:877
StdIntList::NEx
int NEx
Number of exchange term, counted in each lattice file.
Definition: StdFace_vals.h:205
StdIntList::method
char method[256]
The name of method, input from file.
Definition: StdFace_vals.h:259
StdIntList::intrindx
int ** intrindx
[StdIntList::nintr][8] Site/spin indices of two-body term, malloc in StdFace_MallocInteractions() and...
Definition: StdFace_vals.h:173
StdFace_PrintVal_d
void StdFace_PrintVal_d(char *valname, double *val, double val0)
Print a valiable (real) read from the input file if it is not specified in the input file (=NaN)...
Definition: StdFace_ModelUtil.c:414
StdIntList::NCell
int NCell
The number of the unit cell in the super-cell (determinant of StdIntList::box). Set in StdFace_InitSi...
Definition: StdFace_vals.h:49
StdFace_RequiredVal_i
void StdFace_RequiredVal_i(char *valname, int val)
Stop HPhi if a variable (integer) which must be specified is absent in the input file (=2147483647...
Definition: StdFace_ModelUtil.c:577
StdFace_FoldSite
static void StdFace_FoldSite(struct StdIntList *StdI, int iCellV[3], int nBox[3], int iCellV_fold[3])
Move a site into the original supercell if it is outside the original supercell.
Definition: StdFace_ModelUtil.c:595
StdIntList::PairHopp
double * PairHopp
[StdIntList::NPairLift] Coefficient of pair-hopping term, malloc in StdFace_MallocInteractions() and ...
Definition: StdFace_vals.h:226
StdFace_NotUsed_c
void StdFace_NotUsed_c(char *valname, double complex val)
Stop HPhi if a variable (complex) not used is specified in the input file (!=NaN).
Definition: StdFace_ModelUtil.c:509
StdFace_PrintVal_c
void StdFace_PrintVal_c(char *valname, double complex *val, double complex val0)
Print a valiable (complex) read from the input file if it is not specified in the input file (=NaN)...
Definition: StdFace_ModelUtil.c:455
StdIntList::nsite
int nsite
Number of sites, set in the each lattice file.
Definition: StdFace_vals.h:161
StdIntList::transindx
int ** transindx
[StdIntList::ntrans][4] Site/spin indices of one-body term, malloc in StdFace_MallocInteractions() an...
Definition: StdFace_vals.h:165
StdIntList::Height
int Height
Number of sites along the 3rd axis, input parameter.
Definition: StdFace_vals.h:41
StdFace_vals.h
Variables used in the Standard mode. These variables are passed as a pointer of the structure(StdIntL...
StdIntList::Cintra
double * Cintra
[StdIntList::NCintra] Coefficient of intra-site Coulomb term, malloc in StdFace_MallocInteractions() ...
Definition: StdFace_vals.h:185
StdFace_NotUsed_i
void StdFace_NotUsed_i(char *valname, int val)
Stop HPhi if a variable (integer) not used is specified in the input file (!=2147483647, the upper limt of Int).
Definition: StdFace_ModelUtil.c:558
StdFace_PrintXSF
void StdFace_PrintXSF(struct StdIntList *StdI)
Print lattice.xsf (XCrysDen format)
Definition: StdFace_ModelUtil.c:940
StdIntList::nintr
int nintr
Number of InterAll, counted in each lattice file.
Definition: StdFace_vals.h:171
StdIntList::h
double h
Longitudinal magnetic field, input parameter.
Definition: StdFace_vals.h:147
StdFace_FindSite
void StdFace_FindSite(struct StdIntList *StdI, int iW, int iL, int iH, int diW, int diL, int diH, int isiteUC, int jsiteUC, int *isite, int *jsite, double complex *Cphase, double *dR)
Find the index of transfer and interaction.
Definition: StdFace_ModelUtil.c:823
StdIntList::PHIndx
int ** PHIndx
[StdIntList::NPairLift][2] Site indices of pair-hopping term, malloc in StdFace_MallocInteractions() ...
Definition: StdFace_vals.h:223
StdFace_intr
void StdFace_intr(struct StdIntList *StdI, double complex intr0, int site1, int spin1, int site2, int spin2, int site3, int spin3, int site4, int spin4)
Add interaction (InterAll) to the list Set StdIntList::intr and StdIntList::intrindx and increase the...
Definition: StdFace_ModelUtil.c:199
StdIntList::NaN_i
int NaN_i
It is used for initializing input parameter. This means that a parameter wich is not specified in inp...
Definition: StdFace_vals.h:28
StdIntList::Cinter
double * Cinter
[StdIntList::NCinter] Coefficient of inter-site Coulomb term, malloc in StdFace_MallocInteractions() ...
Definition: StdFace_vals.h:194
StdIntList::tau
double ** tau
Cell-internal site position in the fractional coordinate. Defined in the beginning of each lattice fu...
Definition: StdFace_vals.h:55
StdFace_exit
void StdFace_exit(int errorcode)
MPI Abortation wrapper.
Definition: StdFace_ModelUtil.c:35
StdFace_InputSpinNN
void StdFace_InputSpinNN(double J[3][3], double JAll, double J0[3][3], double J0All, char *J0name)
Input nearest-neighbor spin-spin interaction.
Definition: StdFace_ModelUtil.c:975
StdIntList::CinterIndx
int ** CinterIndx
[StdIntList::NCinter][2] Site indices of inter-site Coulomb term, malloc in StdFace_MallocInteraction...
Definition: StdFace_vals.h:191