HΦ  3.1.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  fprintf(fp, "set arrow from %f, %f to %f, %f nohead ls %d\n", xi, yi, xj, yj, connect);
917  StdFace_FindSite(StdI, iW, iL, 0, diW, diL, 0, isiteUC, jsiteUC, isite, jsite, Cphase, dR);
918 
919  xi = StdI->direct[1][0] * ((double)iL + StdI->tau[isiteUC][1])
920  + StdI->direct[0][0] * ((double)iW + StdI->tau[isiteUC][0]);
921  yi = StdI->direct[1][1] * ((double)iL + StdI->tau[isiteUC][1])
922  + StdI->direct[0][1] * ((double)iW + StdI->tau[isiteUC][0]);
923 
924  xj = StdI->direct[0][0] * ((double)(iW + diW) + StdI->tau[jsiteUC][0])
925  + StdI->direct[1][0] * ((double)(iL + diL) + StdI->tau[jsiteUC][1]);
926  yj = StdI->direct[0][1] * ((double)(iW + diW) + StdI->tau[jsiteUC][0])
927  + StdI->direct[1][1] * ((double)(iL + diL) + StdI->tau[jsiteUC][1]);
928 
929  if (*isite < 10)fprintf(fp, "set label \"%1d\" at %f, %f center front\n", *isite, xi, yi);
930  else fprintf(fp, "set label \"%2d\" at %f, %f center front\n", *isite, xi, yi);
931  if (*jsite < 10)fprintf(fp, "set label \"%1d\" at %f, %f center front\n", *jsite, xj, yj);
932  else fprintf(fp, "set label \"%2d\" at %f, %f center front\n", *jsite, xj, yj);
933  fprintf(fp, "set arrow from %f, %f to %f, %f nohead ls %d\n", xi, yi, xj, yj, connect);
934 }/*void StdFace_SetLabel*/
938 void StdFace_PrintXSF(struct StdIntList *StdI) {
939  FILE *fp;
940  int ii, jj, kk, isite, iCell;
941  double vec[3];
942 
943  fp = fopen("lattice.xsf", "w");
944  fprintf(fp, "CRYSTAL\n");
945  fprintf(fp, "PRIMVEC\n");
946  for (ii = 0; ii < 3; ii++) {
947  for (jj = 0; jj < 3; jj++) {
948  vec[jj] = 0.0;
949  for (kk = 0; kk < 3; kk++)
950  vec[jj] += (double)StdI->box[ii][kk] * StdI->direct[kk][jj];
951  }
952  fprintf(fp, "%15.5f %15.5f %15.5f\n", vec[0], vec[1], vec[2]);
953  }
954  fprintf(fp, "PRIMCOORD\n");
955  fprintf(fp, "%d 1\n", StdI->NCell * StdI->NsiteUC);
956  for (iCell = 0; iCell < StdI->NCell; iCell++) {
957  for (isite = 0; isite < StdI->NsiteUC; isite++) {
958  for (jj = 0; jj < 3; jj++) {
959  vec[jj] = 0.0;
960  for (kk = 0; kk < 3; kk++)
961  vec[jj] += ((double)StdI->Cell[iCell][kk] + StdI->tau[isite][kk])
962  * StdI->direct[kk][jj];
963  }
964  fprintf(fp, "H %15.5f %15.5f %15.5f\n", vec[0], vec[1], vec[2]);
965  }
966  }
967  fflush(fp);
968  fclose(fp);
969 }/*void StdFace_PrintXSF*/
973 void StdFace_InputSpinNN(
974  struct StdIntList *StdI,
975  double J0[3][3],
976  double J0All,
977  char *J0name
978 )
979 {
980  int i1, i2, i3, i4;
981  char Jname[3][3][10];
982 
983  strcpy(Jname[0][0], "x\0");
984  strcpy(Jname[0][1], "xy\0");
985  strcpy(Jname[0][2], "xz\0");
986  strcpy(Jname[1][0], "yx\0");
987  strcpy(Jname[1][1], "y\0");
988  strcpy(Jname[1][2], "yz\0");
989  strcpy(Jname[2][0], "zx\0");
990  strcpy(Jname[2][1], "zy\0");
991  strcpy(Jname[2][2], "z\0");
992 
993  if (isnan(StdI->JAll) == 0 && isnan(J0All) == 0) {
994  fprintf(stdout, "\n ERROR! J and %s conflict !\n\n", J0name);
995  StdFace_exit(-1);
996  }
997  for (i1 = 0; i1 < 3; i1++) {
998  for (i2 = 0; i2 < 3; i2++) {
999  if (isnan(StdI->JAll) == 0 && isnan(StdI->J[i1][i2]) == 0) {
1000  fprintf(stdout, "\n ERROR! J and J%s conflict !\n\n", Jname[i1][i2]);
1001  StdFace_exit(-1);
1002  }
1003  else if (isnan(J0All) == 0 && isnan(StdI->J[i1][i2]) == 0) {
1004  fprintf(stdout, "\n ERROR! %s and J%s conflict !\n\n",
1005  J0name, Jname[i1][i2]);
1006  StdFace_exit(-1);
1007  }
1008  else if (isnan(J0All) == 0 && isnan(J0[i1][i2]) == 0) {
1009  fprintf(stdout, "\n ERROR! %s and %s%s conflict !\n\n", J0name,
1010  J0name, Jname[i1][i2]);
1011  StdFace_exit(-1);
1012  }
1013  else if (isnan(J0[i1][i2]) == 0 && isnan(StdI->JAll) == 0) {
1014  fprintf(stdout, "\n ERROR! %s%s and J conflict !\n\n",
1015  J0name, Jname[i1][i2]);
1016  StdFace_exit(-1);
1017  }
1018  }/*for (j = 0; j < 3; j++)*/
1019  }/*for (i = 0; i < 3; i++)*/
1020 
1021  for (i1 = 0; i1 < 3; i1++) {
1022  for (i2 = 0; i2 < 3; i2++) {
1023  for (i3 = 0; i3 < 3; i3++) {
1024  for (i4 = 0; i4 < 3; i4++) {
1025  if (isnan(J0[i1][i2]) == 0 && isnan(StdI->J[i3][i4]) == 0) {
1026  fprintf(stdout, "\n ERROR! %s%s and J%s conflict !\n\n",
1027  J0name, Jname[i1][i2], Jname[i3][i4]);
1028  StdFace_exit(-1);
1029  }
1030  }/*for (i4 = 0; i4 < 3; i4++)*/
1031  }/*for (i3 = 0; i3 < 3; i3++)*/
1032  }/*for (j = 0; j < 3; j++)*/
1033  }/*for (i = 0; i < 3; i++)*/
1034 
1035  for (i1 = 0; i1 < 3; i1++) {
1036  for (i2 = 0; i2 < 3; i2++) {
1037  if (isnan(J0[i1][i2]) == 0)
1038  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1039  else if (isnan(StdI->J[i1][i2]) == 0) {
1040  J0[i1][i2] = StdI->J[i1][i2];
1041  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1042  }
1043  else if (i1 == i2 && isnan(J0All) == 0) {
1044  J0[i1][i2] = J0All;
1045  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1046  }
1047  else if (i1 == i2 && isnan(StdI->JAll) == 0) {
1048  J0[i1][i2] = StdI->JAll;
1049  fprintf(stdout, " %14s%s = %-10.5f\n", J0name, Jname[i1][i2], J0[i1][i2]);
1050  }
1051  else {
1052  J0[i1][i2] = 0.0;
1053  }
1054  }/*for (i2 = 0; i2 < 3; i2++)*/
1055  }/*for (i = 0; i < 3; i++)*/
1056 }/*void StdFace_InputSpinNN*/
1060 void StdFace_InputSpin(
1061  struct StdIntList *StdI,
1062  double Jp[3][3],
1063  double JpAll,
1064  char *Jpname
1065 )
1066 {
1067  int i1, i2;
1068  char Jname[3][3][10];
1069 
1070  strcpy(Jname[0][0], "x\0");
1071  strcpy(Jname[0][1], "xy\0");
1072  strcpy(Jname[0][2], "xz\0");
1073  strcpy(Jname[1][0], "yx\0");
1074  strcpy(Jname[1][1], "y\0");
1075  strcpy(Jname[1][2], "yz\0");
1076  strcpy(Jname[2][0], "zx\0");
1077  strcpy(Jname[2][1], "zy\0");
1078  strcpy(Jname[2][2], "z\0");
1079 
1080  for (i1 = 0; i1 < 3; i1++) {
1081  for (i2 = 0; i2 < 3; i2++) {
1082  if (isnan(JpAll) == 0 && isnan(Jp[i1][i2]) == 0) {
1083  fprintf(stdout, "\n ERROR! %s and %s%s conflict !\n\n", Jpname,
1084  Jpname, Jname[i1][i2]);
1085  StdFace_exit(-1);
1086  }
1087  }/*for (j = 0; j < 3; j++)*/
1088  }/*for (i = 0; i < 3; i++)*/
1089 
1090  for (i1 = 0; i1 < 3; i1++) {
1091  for (i2 = 0; i2 < 3; i2++) {
1092  if (isnan(Jp[i1][i2]) == 0)
1093  fprintf(stdout, " %14s%s = %-10.5f\n", Jpname, Jname[i1][i2], Jp[i1][i2]);
1094  else if (i1 == i2 && isnan(JpAll) == 0) {
1095  Jp[i1][i2] = JpAll;
1096  fprintf(stdout, " %14s%s = %-10.5f\n", Jpname, Jname[i1][i2], Jp[i1][i2]);
1097  }
1098  else {
1099  Jp[i1][i2] = 0.0;
1100  }
1101  }/*for (i2 = 0; i2 < 3; i2++)*/
1102  }/*for (i = 0; i < 3; i++)*/
1103 }/*void StdFace_InputSpin*/
1109 void StdFace_InputCoulombV(
1110  struct StdIntList *StdI,
1111  double *V0,
1112  char *V0name
1113 )
1114 {
1115  if (isnan(StdI->V) == 0 && isnan(*V0) == 0) {
1116  fprintf(stdout, "\n ERROR! V and %s conflict !\n\n", V0name);
1117  StdFace_exit(-1);
1118  }
1119  else if (isnan(*V0) == 0)
1120  fprintf(stdout, " %15s = %-10.5f\n", V0name, *V0);
1121  else if (isnan(StdI->V) == 0) {
1122  *V0 = StdI->V;
1123  fprintf(stdout, " %15s = %-10.5f\n", V0name, *V0);
1124  }
1125  else {
1126  *V0 = 0.0;
1127  }
1128 }/*void StdFace_InputCoulombV*/
1134 void StdFace_InputHopp(
1135  struct StdIntList *StdI,
1136  double complex *t0,
1137  char *t0name
1138 )
1139 {
1140  if (isnan(creal(StdI->t)) == 0 && isnan(creal(*t0)) == 0) {
1141  fprintf(stdout, "\n ERROR! t and %s conflict !\n\n", t0name);
1142  StdFace_exit(-1);
1143  }
1144  else if (isnan(creal(*t0)) == 0)
1145  fprintf(stdout, " %15s = %-10.5f %-10.5f\n", t0name, creal(*t0), cimag(*t0));
1146  else if (isnan(creal(StdI->t)) == 0) {
1147  *t0 = StdI->t;
1148  fprintf(stdout, " %15s = %-10.5f %-10.5f\n", t0name, creal(*t0), cimag(*t0));
1149  }
1150  else {
1151  *t0 = 0.0;
1152  }
1153 }/*void StdFace_InputHopp*/
1157 void StdFace_PrintGeometry(struct StdIntList *StdI) {
1158  FILE *fp;
1159  int isite, iCell, ii;
1160 
1161  fp = fopen("geometry.dat", "w");
1162 
1163  for (ii = 0; ii < 3; ii++)
1164  fprintf(fp, "%25.15e %25.15e %25.15e\n",
1165  StdI->direct[ii][0], StdI->direct[ii][1], StdI->direct[ii][2]);
1166  fprintf(fp, "%25.15e %25.15e %25.15e\n",
1167  StdI->phase[0], StdI->phase[1], StdI->phase[2]);
1168  for (ii = 0; ii < 3; ii++)
1169  fprintf(fp, "%d %d %d\n",
1170  StdI->box[ii][0], StdI->box[ii][1], StdI->box[ii][2]);
1171 
1172  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1173  for (isite = 0; isite < StdI->NsiteUC; isite++) {
1174  fprintf(fp, "%d %d %d %d\n",
1175  StdI->Cell[iCell][0] - StdI->Cell[0][0],
1176  StdI->Cell[iCell][1] - StdI->Cell[0][1],
1177  StdI->Cell[iCell][2] - StdI->Cell[0][2],
1178  isite);
1179  }/*for (isite = 0; isite < StdI->NsiteUC; isite++)*/
1180  }/* for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1181  if (strcmp(StdI->model, "kondo") == 0) {
1182  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1183  for (isite = 0; isite < StdI->NsiteUC; isite++) {
1184  fprintf(fp, "%d %d %d %d\n",
1185  StdI->Cell[iCell][0] - StdI->Cell[0][0],
1186  StdI->Cell[iCell][1] - StdI->Cell[0][1],
1187  StdI->Cell[iCell][2] - StdI->Cell[0][2],
1188  isite + StdI->NsiteUC);
1189  }/*for (isite = 0; isite < StdI->NsiteUC; isite++)*/
1190  }/* for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1191  }
1192  fflush(fp);
1193  fclose(fp);
1194 }/*void StdFace_PrintGeometry()*/
1198 void StdFace_MallocInteractions(
1199  struct StdIntList *StdI,
1200  int ntransMax,
1201  int nintrMax
1202 ) {
1203  int ii;
1204 #if defined(_HPhi)
1205  int it;
1206 #endif
1207 
1210  StdI->transindx = (int **)malloc(sizeof(int*) * ntransMax);
1211  StdI->trans = (double complex *)malloc(sizeof(double complex) * ntransMax);
1212  for (ii = 0; ii < ntransMax; ii++) {
1213  StdI->transindx[ii] = (int *)malloc(sizeof(int) * 4);
1214  }
1215  StdI->ntrans = 0;
1216 #if defined(_HPhi)
1217  if (strcmp(StdI->method, "timeevolution") == 0 && StdI->PumpBody == 1) {
1218  StdI->npump = (int *)malloc(sizeof(int) * StdI->Lanczos_max);
1219  StdI->pumpindx = (int ***)malloc(sizeof(int**) * StdI->Lanczos_max);
1220  StdI->pump = (double complex **)malloc(sizeof(double complex*) * StdI->Lanczos_max);
1221  for (it = 0; it < StdI->Lanczos_max; it++) {
1222  StdI->npump[it] = 0;
1223  StdI->pumpindx[it] = (int **)malloc(sizeof(int*) * ntransMax);
1224  StdI->pump[it] = (double complex *)malloc(sizeof(double complex) * ntransMax);
1225  for (ii = 0; ii < ntransMax; ii++) {
1226  StdI->pumpindx[it][ii] = (int *)malloc(sizeof(int) * 4);
1227  }
1228  }/*for (it = 0; it < StdI->Lanczos_max;)*/
1229  }/*if (strcmp(StdI->method, "timeevolution") == 0*/
1230 #endif
1231 
1234  StdI->intrindx = (int **)malloc(sizeof(int*) * nintrMax);
1235  StdI->intr = (double complex *)malloc(sizeof(double complex) * nintrMax);
1236  for (ii = 0; ii < nintrMax; ii++) {
1237  StdI->intrindx[ii] = (int *)malloc(sizeof(int) * 8);
1238  }
1239  StdI->nintr = 0;
1243  StdI->CintraIndx = (int **)malloc(sizeof(int*) * nintrMax);
1244  StdI->Cintra = (double *)malloc(sizeof(double) * nintrMax);
1245  for (ii = 0; ii < nintrMax; ii++) {
1246  StdI->CintraIndx[ii] = (int *)malloc(sizeof(int) * 1);
1247  }
1248  StdI->NCintra = 0;
1252  StdI->CinterIndx = (int **)malloc(sizeof(int*) * nintrMax);
1253  StdI->Cinter = (double *)malloc(sizeof(double) * nintrMax);
1254  for (ii = 0; ii < nintrMax; ii++) {
1255  StdI->CinterIndx[ii] = (int *)malloc(sizeof(int) * 2);
1256  }
1257  StdI->NCinter = 0;
1261  StdI->HundIndx = (int **)malloc(sizeof(int*) * nintrMax);
1262  StdI->Hund = (double *)malloc(sizeof(double) * nintrMax);
1263  for (ii = 0; ii < nintrMax; ii++) {
1264  StdI->HundIndx[ii] = (int *)malloc(sizeof(int) * 2);
1265  }
1266  StdI->NHund = 0;
1270  StdI->ExIndx = (int **)malloc(sizeof(int*) * nintrMax);
1271  StdI->Ex = (double *)malloc(sizeof(double) * nintrMax);
1272  for (ii = 0; ii < nintrMax; ii++) {
1273  StdI->ExIndx[ii] = (int *)malloc(sizeof(int) * 2);
1274  }
1275  StdI->NEx = 0;
1279  StdI->PLIndx = (int **)malloc(sizeof(int*) * nintrMax);
1280  StdI->PairLift = (double *)malloc(sizeof(double) * nintrMax);
1281  for (ii = 0; ii < nintrMax; ii++) {
1282  StdI->PLIndx[ii] = (int *)malloc(sizeof(int) * 2);
1283  }
1284  StdI->NPairLift = 0;
1288  StdI->PHIndx = (int **)malloc(sizeof(int*) * nintrMax);
1289  StdI->PairHopp = (double *)malloc(sizeof(double) * nintrMax);
1290  for (ii = 0; ii < nintrMax; ii++) {
1291  StdI->PHIndx[ii] = (int *)malloc(sizeof(int) * 2);
1292  }
1293  StdI->NPairHopp = 0;
1294 }/*void StdFace_MallocInteractions*/
1295 #if defined(_mVMC)
1296 
1300 static void StdFace_FoldSiteSub(
1301  struct StdIntList *StdI,
1302  int iCellV[3],
1303  int nBox[3],
1304  int iCellV_fold[3]
1305 )
1306 {
1307  int ii, jj, iCellV_frac[3];
1311  for (ii = 0; ii < 3; ii++) {
1312  iCellV_frac[ii] = 0;
1313  for (jj = 0; jj < 3; jj++)iCellV_frac[ii] += StdI->rboxsub[ii][jj] * iCellV[jj];
1314  }
1318  for (ii = 0; ii < 3; ii++)
1319  nBox[ii] = (iCellV_frac[ii] + StdI->NCellsub * 1000) / StdI->NCellsub - 1000;
1323  for (ii = 0; ii < 3; ii++)
1324  iCellV_frac[ii] -= StdI->NCellsub*(nBox[ii]);
1325 
1326  for (ii = 0; ii < 3; ii++) {
1327  iCellV_fold[ii] = 0;
1328  for (jj = 0; jj < 3; jj++) iCellV_fold[ii] += StdI->boxsub[jj][ii] * iCellV_frac[jj];
1329  iCellV_fold[ii] = (iCellV_fold[ii] + StdI->NCellsub * 1000) / StdI->NCellsub - 1000;
1330  }
1331 }/*static void StdFace_FoldSiteSub*/
1336 void StdFace_Proj(struct StdIntList *StdI)
1337 {
1338  FILE *fp;
1339  int jsite, iCell, jCell, kCell;
1340  int nBox[3], iCellV[3], jCellV[3], ii;
1341  int iSym;
1342  int **Sym, **Anti;
1343 
1344  Sym = (int **)malloc(sizeof(int*) * StdI->nsite);
1345  Anti = (int **)malloc(sizeof(int*) * StdI->nsite);
1346  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1347  Sym[jsite] = (int *)malloc(sizeof(int) * StdI->nsite);
1348  Anti[jsite] = (int *)malloc(sizeof(int) * StdI->nsite);
1349  }
1353  StdI->NSym = 0;
1354  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1355 
1356  StdFace_FoldSiteSub(StdI, StdI->Cell[iCell], nBox, iCellV);
1357 
1358  StdFace_FoldSite(StdI, iCellV, nBox, iCellV);
1359 
1360  if (iCellV[0] == StdI->Cell[iCell][0] &&
1361  iCellV[1] == StdI->Cell[iCell][1] &&
1362  iCellV[2] == StdI->Cell[iCell][2]) {
1366  for (jCell = 0; jCell < StdI->NCell; jCell++) {
1367 
1368  for (ii = 0; ii < 3; ii++)jCellV[ii] = StdI->Cell[jCell][ii] + iCellV[ii];
1369  StdFace_FoldSite(StdI, jCellV, nBox, jCellV);
1370 
1371  for (kCell = 0; kCell < StdI->NCell; kCell++) {
1372  if (jCellV[0] == StdI->Cell[kCell][0] &&
1373  jCellV[1] == StdI->Cell[kCell][1] &&
1374  jCellV[2] == StdI->Cell[kCell][2])
1375  {
1376 
1377  for (jsite = 0; jsite < StdI->NsiteUC; jsite++) {
1378 
1379  Sym[StdI->NSym][jCell*StdI->NsiteUC + jsite] = kCell*StdI->NsiteUC + jsite;
1380  Anti[StdI->NSym][jCell*StdI->NsiteUC + jsite]
1381  = StdI->AntiPeriod[0] * nBox[0]
1382  + StdI->AntiPeriod[1] * nBox[1]
1383  + StdI->AntiPeriod[2] * nBox[2];
1384 
1385  if (strcmp(StdI->model, "kondo") == 0) {
1386  Sym[StdI->NSym][StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite] = StdI->nsite / 2 + kCell*StdI->NsiteUC + jsite;
1387  Anti[StdI->NSym][StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite]
1388  = StdI->AntiPeriod[0] * nBox[0]
1389  + StdI->AntiPeriod[1] * nBox[1]
1390  + StdI->AntiPeriod[2] * nBox[2];
1391  }/*if (strcmp(StdI->model, "kondo") == 0)*/
1392 
1393  }/*for (jsite = 0; jsite < StdI->NsiteUC; jsite++)*/
1394 
1395  }/*if (jWfold == StdI->Cell[kCell][0] && jLfold == StdI->Cell[kCell][1])*/
1396  }/*for (kCell = 0; kCell < StdI->NCell; kCell++)*/
1397  }/*for (jCell = 0; jCell < StdI->NCell; jCell++)*/
1398  StdI->NSym += 1;
1399  }/*if (iWfold == iW && iLfold == iL)*/
1400  }/*for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1401 
1402  fp = fopen("qptransidx.def", "w");
1403  fprintf(fp, "=============================================\n");
1404  fprintf(fp, "NQPTrans %10d\n", StdI->NSym);
1405  fprintf(fp, "=============================================\n");
1406  fprintf(fp, "======== TrIdx_TrWeight_and_TrIdx_i_xi ======\n");
1407  fprintf(fp, "=============================================\n");
1408  for (iSym = 0; iSym < StdI->NSym; iSym++) {
1409  fprintf(fp, "%d %10.5f\n", iSym, 1.0);
1410  }
1411  for (iSym = 0; iSym < StdI->NSym; iSym++) {
1412  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1413  if (Anti[iSym][jsite] % 2 == 0) Anti[iSym][jsite] = 1;
1414  else Anti[iSym][jsite] = -1;
1415  if (StdI->AntiPeriod[0] == 1 || StdI->AntiPeriod[1] == 1 || StdI->AntiPeriod[2] == 1) {
1416  fprintf(fp, "%5d %5d %5d %5d\n", iSym, jsite, Sym[iSym][jsite], Anti[iSym][jsite]);
1417  }
1418  else {
1419  fprintf(fp, "%5d %5d %5d\n", iSym, jsite, Sym[iSym][jsite]);
1420  }
1421  }
1422  }
1423  fflush(fp);
1424  fclose(fp);
1425  fprintf(stdout, " qptransidx.def is written.\n");
1426 
1427  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1428  free(Anti[jsite]);
1429  free(Sym[jsite]);
1430  }
1431  free(Sym);
1432  free(Anti);
1433 }/*void StdFace_Proj(struct StdIntList *StdI)*/
1438 static void StdFace_InitSiteSub(struct StdIntList *StdI)
1439 {
1440  int ii, jj, kk, prod;
1444  if ((StdI->Lsub != StdI->NaN_i || StdI->Wsub != StdI->NaN_i || StdI->Hsub != StdI->NaN_i)
1445  && (StdI->boxsub[0][0] != StdI->NaN_i || StdI->boxsub[0][1] != StdI->NaN_i || StdI->boxsub[0][2] != StdI->NaN_i ||
1446  StdI->boxsub[1][0] != StdI->NaN_i || StdI->boxsub[1][1] != StdI->NaN_i || StdI->boxsub[1][2] != StdI->NaN_i ||
1447  StdI->boxsub[2][0] != StdI->NaN_i || StdI->boxsub[2][1] != StdI->NaN_i || StdI->boxsub[2][2] != StdI->NaN_i))
1448  {
1449  fprintf(stdout, "\nERROR ! (Lsub, Wsub, Hsub) and (a0Wsub, ..., a2Hsub) conflict !\n\n");
1450  StdFace_exit(-1);
1451  }
1452  else if (StdI->Wsub != StdI->NaN_i || StdI->Lsub != StdI->NaN_i || StdI->Hsub != StdI->NaN_i) {
1453  StdFace_PrintVal_i("Lsub", &StdI->Lsub, 1);
1454  StdFace_PrintVal_i("Wsub", &StdI->Wsub, 1);
1455  StdFace_PrintVal_i("Hsub", &StdI->Hsub, 1);
1456  for (ii = 0; ii < 3; ii++) for (jj = 0; jj < 3; jj++)
1457  StdI->boxsub[ii][jj] = 0;
1458  StdI->boxsub[0][0] = StdI->Wsub;
1459  StdI->boxsub[1][1] = StdI->Lsub;
1460  StdI->boxsub[2][2] = StdI->Hsub;
1461  }
1462  else {
1463  StdFace_PrintVal_i("a0Wsub", &StdI->boxsub[0][0], StdI->box[0][0]);
1464  StdFace_PrintVal_i("a0Lsub", &StdI->boxsub[0][1], StdI->box[0][1]);
1465  StdFace_PrintVal_i("a0Hsub", &StdI->boxsub[0][2], StdI->box[0][2]);
1466  StdFace_PrintVal_i("a1Wsub", &StdI->boxsub[1][0], StdI->box[1][0]);
1467  StdFace_PrintVal_i("a1Lsub", &StdI->boxsub[1][1], StdI->box[1][1]);
1468  StdFace_PrintVal_i("a1Hsub", &StdI->boxsub[1][2], StdI->box[1][2]);
1469  StdFace_PrintVal_i("a2Wsub", &StdI->boxsub[2][0], StdI->box[2][0]);
1470  StdFace_PrintVal_i("a2Lsub", &StdI->boxsub[2][1], StdI->box[2][1]);
1471  StdFace_PrintVal_i("a2Hsub", &StdI->boxsub[2][2], StdI->box[2][2]);
1472  }
1476  StdI->NCellsub = 0;
1477  for (ii = 0; ii < 3; ii++) {
1478  StdI->NCellsub += StdI->boxsub[0][ii]
1479  * StdI->boxsub[1][(ii + 1) % 3]
1480  * StdI->boxsub[2][(ii + 2) % 3]
1481  - StdI->boxsub[0][ii]
1482  * StdI->boxsub[1][(ii + 2) % 3]
1483  * StdI->boxsub[2][(ii + 1) % 3];
1484  }
1485  printf(" Number of Cell in the sublattice: %d\n", abs(StdI->NCellsub));
1486  if (StdI->NCellsub == 0) {
1487  StdFace_exit(-1);
1488  }
1489 
1490  for (ii = 0; ii < 3; ii++) {
1491  for (jj = 0; jj < 3; jj++) {
1492  StdI->rboxsub[ii][jj] = StdI->boxsub[(ii + 1) % 3][(jj + 1) % 3] * StdI->boxsub[(ii + 2) % 3][(jj + 2) % 3]
1493  - StdI->boxsub[(ii + 1) % 3][(jj + 2) % 3] * StdI->boxsub[(ii + 2) % 3][(jj + 1) % 3];
1494  }
1495  }
1496  if (StdI->NCellsub < 0) {
1497  for (ii = 0; ii < 3; ii++)
1498  for (jj = 0; jj < 3; jj++)
1499  StdI->rboxsub[ii][jj] *= -1;
1500  StdI->NCellsub *= -1;
1501  }/*if (StdI->NCell < 0)*/
1505  for (ii = 0; ii < 3; ii++) {
1506  for (jj = 0; jj < 3; jj++) {
1507  prod = 0.0;
1508  for (kk = 0; kk < 3; kk++) prod += StdI->rboxsub[ii][kk] * (double)StdI->box[jj][kk];
1509  if (prod % StdI->NCellsub != 0) {
1510  printf("\n ERROR ! Sublattice is INCOMMENSURATE !\n\n");
1511  StdFace_exit(-1);
1512  }/*if (prod % StdI->NCellsub != 0)*/
1513  }
1514  }
1515 }/*void StdFace_InitSiteSub*/
1519 void StdFace_generate_orb(struct StdIntList *StdI) {
1520  int iCell, jCell, kCell, iCell2, jCell2, iOrb, isite, jsite, Anti;
1521  int nBox[3], iCellV[3], jCellV[3], dCellV[3], ii;
1522  int **CellDone;
1523 
1524  StdFace_InitSiteSub(StdI);
1525 
1526  StdI->Orb = (int **)malloc(sizeof(int*) * StdI->nsite);
1527  StdI->AntiOrb = (int **)malloc(sizeof(int*) * StdI->nsite);
1528  for (isite = 0; isite < StdI->nsite; isite++) {
1529  StdI->Orb[isite] = (int *)malloc(sizeof(int) * StdI->nsite);
1530  StdI->AntiOrb[isite] = (int *)malloc(sizeof(int) * StdI->nsite);
1531  }
1532  CellDone = (int **)malloc(sizeof(int*) * StdI->NCell);
1533  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1534  CellDone[iCell] = (int *)malloc(sizeof(int) * StdI->NCell);
1535  for (jCell = 0; jCell < StdI->NCell; jCell++) {
1536  CellDone[iCell][jCell] = 0;
1537  }
1538  }
1539 
1540  iOrb = 0;
1541  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1542 
1543  StdFace_FoldSiteSub(StdI, StdI->Cell[iCell], nBox, iCellV);
1544 
1545  StdFace_FoldSite(StdI, iCellV, nBox, iCellV);
1546 
1547  for (kCell = 0; kCell < StdI->NCell; kCell++) {
1548  if (iCellV[0] == StdI->Cell[kCell][0] &&
1549  iCellV[1] == StdI->Cell[kCell][1] &&
1550  iCellV[2] == StdI->Cell[kCell][2])
1551  {
1552  iCell2 = kCell;
1553  }
1554  }/*for (kCell = 0; kCell < StdI->NCell; kCell++)*/
1555 
1556  for (jCell = 0; jCell < StdI->NCell; jCell++) {
1557 
1558  for (ii = 0; ii < 3; ii++)
1559  jCellV[ii] = StdI->Cell[jCell][ii] + iCellV[ii] - StdI->Cell[iCell][ii];
1560 
1561  StdFace_FoldSite(StdI, jCellV, nBox, jCellV);
1562 
1563  for (kCell = 0; kCell < StdI->NCell; kCell++) {
1564  if (jCellV[0] == StdI->Cell[kCell][0] &&
1565  jCellV[1] == StdI->Cell[kCell][1] &&
1566  jCellV[2] == StdI->Cell[kCell][2])
1567  {
1568  jCell2 = kCell;
1569  }
1570  }/*for (kCell = 0; kCell < StdI->NCell; kCell++)*/
1571  /*
1572  AntiPeriodic factor
1573  */
1574  for (ii = 0; ii < 3; ii++)
1575  dCellV[ii] = StdI->Cell[jCell][ii] - StdI->Cell[iCell][ii];
1576  StdFace_FoldSite(StdI, dCellV, nBox, dCellV);
1577  Anti = 0;
1578  for (ii = 0; ii < 3; ii++)Anti += StdI->AntiPeriod[ii] * nBox[ii];
1579  if (Anti % 2 == 0) Anti = 1;
1580  else Anti = -1;
1581 
1582  for (isite = 0; isite < StdI->NsiteUC; isite++) {
1583  for (jsite = 0; jsite < StdI->NsiteUC; jsite++) {
1584 
1585  if (CellDone[iCell2][jCell2] == 0) {
1586  StdI->Orb[iCell2*StdI->NsiteUC + isite][jCell2*StdI->NsiteUC + jsite] = iOrb;
1587  StdI->AntiOrb[iCell2*StdI->NsiteUC + isite][jCell2*StdI->NsiteUC + jsite] = Anti;
1588  iOrb += 1;
1589  }
1590  StdI->Orb[iCell*StdI->NsiteUC + isite][jCell*StdI->NsiteUC + jsite]
1591  = StdI->Orb[iCell2*StdI->NsiteUC + isite][jCell2*StdI->NsiteUC + jsite];
1592  StdI->AntiOrb[iCell*StdI->NsiteUC + isite][jCell*StdI->NsiteUC + jsite] = Anti;
1593 
1594  if (strcmp(StdI->model, "kondo") == 0) {
1595  if (CellDone[iCell2][jCell2] == 0) {
1596  StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1597  [ jCell2*StdI->NsiteUC + jsite] = iOrb;
1598  StdI->AntiOrb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1599  [ jCell2*StdI->NsiteUC + jsite] = Anti;
1600  iOrb += 1;
1601  StdI->Orb[ iCell2*StdI->NsiteUC + isite]
1602  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = iOrb;
1603  StdI->AntiOrb[ iCell2*StdI->NsiteUC + isite]
1604  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = Anti;
1605  iOrb += 1;
1606  StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1607  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = iOrb;
1608  StdI->AntiOrb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1609  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite] = Anti;
1610  iOrb += 1;
1611  }
1612  StdI->Orb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1613  [ jCell*StdI->NsiteUC + jsite]
1614  = StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1615  [ jCell2*StdI->NsiteUC + jsite];
1616  StdI->AntiOrb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1617  [ jCell*StdI->NsiteUC + jsite] = Anti;
1618  StdI->Orb[ iCell*StdI->NsiteUC + isite]
1619  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite]
1620  = StdI->Orb[ iCell2*StdI->NsiteUC + isite]
1621  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite];
1622  StdI->AntiOrb[iCell*StdI->NsiteUC + isite]
1623  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite] = Anti;
1624  StdI->Orb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1625  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite]
1626  = StdI->Orb[StdI->nsite / 2 + iCell2*StdI->NsiteUC + isite]
1627  [StdI->nsite / 2 + jCell2*StdI->NsiteUC + jsite];
1628  StdI->AntiOrb[StdI->nsite / 2 + iCell*StdI->NsiteUC + isite]
1629  [StdI->nsite / 2 + jCell*StdI->NsiteUC + jsite] = Anti;
1630  }/*if (strcmp(StdI->model, "kondo") == 0)*/
1631 
1632  }/*for (jsite = 0; jsite < StdI->NsiteUC; jsite++)*/
1633  }/*for (isite = 0; isite < StdI->NsiteUC; isite++)*/
1634  CellDone[iCell2][jCell2] = 1;
1635  }/*for (jCell = 0; jCell < StdI->NCell; jCell++)*/
1636 
1637  }/*for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1638  StdI->NOrb = iOrb;
1639 
1640  for (iCell = 0; iCell < StdI->NCell; iCell++) free(CellDone[iCell]);
1641  free(CellDone);
1642 }/*void StdFace_generate_orb*/
1647 void PrintJastrow(struct StdIntList *StdI) {
1648  FILE *fp;
1649  int isite, jsite, isiteUC, jsiteUC, revarsal, isite1, jsite1, iorb;
1650  int NJastrow, iJastrow;
1651  int dCell, iCell;//, jCell, dCellv[3];
1652  int **Jastrow;
1653  double complex Cphase;
1654  double dR[3];
1655 
1656  Jastrow = (int **)malloc(sizeof(int*) * StdI->nsite);
1657  for (isite = 0; isite < StdI->nsite; isite++)
1658  Jastrow[isite] = (int *)malloc(sizeof(int) * StdI->nsite);
1659 
1660  if (abs(StdI->NMPTrans) == 1 || StdI->NMPTrans == StdI->NaN_i) {
1664  for (isite = 0; isite < StdI->nsite; isite++) {
1665  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1666  Jastrow[isite][jsite] = StdI->Orb[isite][jsite];
1667  }/*for (jsite = 0; jsite < isite; jsite++)*/
1668  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1672  for (iorb = 0; iorb < StdI->NOrb; iorb++) {
1673  for (isite = 0; isite < StdI->nsite; isite++) {
1674  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1675  if (Jastrow[isite][jsite] == iorb) {
1676  Jastrow[jsite][isite] = Jastrow[isite][jsite];
1677  }
1678  }/*for (jsite = 0; jsite < isite; jsite++)*/
1679  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1680  }/*for (iorb = 0; iorb < StdI->NOrb; iorb++)*/
1681 
1682  if (strcmp(StdI->model, "hubbard") == 0) NJastrow = 0;
1683  else NJastrow = -1;
1684  for (isite = 0; isite < StdI->nsite; isite++) {
1685  /*
1686  For Local spin
1687  */
1688  if (StdI->locspinflag[isite] != 0) {
1689  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1690  Jastrow[isite][jsite] = -1;
1691  Jastrow[jsite][isite] = -1;
1692  }
1693  continue;
1694  }
1695 
1696  for (jsite = 0; jsite < isite; jsite++) {
1697  if (Jastrow[isite][jsite] >= 0) {
1698  iJastrow = Jastrow[isite][jsite];
1699  NJastrow -= 1;
1700  for (isite1 = 0; isite1 < StdI->nsite; isite1++) {
1701  for (jsite1 = 0; jsite1 < StdI->nsite; jsite1++) {
1702  if (Jastrow[isite1][jsite1] == iJastrow)
1703  Jastrow[isite1][jsite1] = NJastrow;
1704  }/*for (jsite1 = 0; jsite1 < StdI->nsite; jsite1++)*/
1705  }/*for (isite1 = 0; isite1 < StdI->nsite; isite1++)*/
1706  }/*if (Jastrow[isite][jsite] >= 0)*/
1707  }/*for (jsite = 0; jsite < isite; jsite++)*/
1708  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1709 
1710  NJastrow = -NJastrow;
1711  for (isite = 0; isite < StdI->nsite; isite++) {
1712  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1713  Jastrow[isite][jsite] = -1 - Jastrow[isite][jsite];
1714  }/*for (jsite = 0; jsite < isite; jsite++)*/
1715  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1716  }/*if (abs(StdI->NMPTrans) == 1)*/
1717  else {
1718 
1719  if (strcmp(StdI->model, "spin") == 0) {
1720  NJastrow = 1;
1721 
1722  for (isite = 0; isite < StdI->nsite; isite++) {
1723  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1724  Jastrow[isite][jsite] = 0;
1725  }/*for (jsite = 0; jsite < StdI->nsite; jsite++)*/
1726  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1727  }/*if (strcmp(StdI->model, "spin") == 0)*/
1728  else {
1729 
1730  NJastrow = 0;
1731 
1732  if (strcmp(StdI->model, "kondo") == 0) {
1733  /*
1734  Local spin - itererant electron part
1735  */
1736  for (isite = 0; isite < StdI->nsite; isite++) {
1737  for (jsite = 0; jsite < StdI->nsite / 2; jsite++) {
1738  Jastrow[isite][jsite] = 0;
1739  Jastrow[jsite][isite] = 0;
1740  }/*for (jsite = 0; jsite < StdI->nsite; jsite++)*/
1741  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1742 
1743  NJastrow += 1;
1744  }/*if (strcmp(StdI->model, "kondo") == 0)*/
1745 
1746  for (dCell = 0; dCell < StdI->NCell; dCell++) {
1747  StdFace_FindSite(StdI,
1748  0, 0, 0,
1749  -StdI->Cell[dCell][0], -StdI->Cell[dCell][1], -StdI->Cell[dCell][2],
1750  0, 0, &isite, &jsite, &Cphase, dR);
1751  if (strcmp(StdI->model, "kondo") == 0) jsite += -StdI->NCell * StdI->NsiteUC;
1752  iCell = jsite / StdI->NsiteUC;
1753  if (iCell < dCell) {
1754  /*
1755  If -R has been already done, skip.
1756  */
1757  continue;
1758  }
1759  else if (iCell == dCell) {
1760  /*
1761  If revarsal symmetry [Fold(-R) = R], J(R,i,j) = J(R,j,i)
1762  */
1763  revarsal = 1;
1764  }
1765  else revarsal = 0;
1766 
1767  for (isiteUC = 0; isiteUC < StdI->NsiteUC; isiteUC++) {
1768  for (jsiteUC = 0; jsiteUC < StdI->NsiteUC; jsiteUC++) {
1769  if (revarsal == 1 && jsiteUC > isiteUC) continue;/*If [Fold(-R) = R]*/
1770  if (isiteUC == jsiteUC &&
1771  StdI->Cell[dCell][0] == 0 &&
1772  StdI->Cell[dCell][1] == 0 &&
1773  StdI->Cell[dCell][2] == 0) continue;/*Diagonal*/
1774 
1775  for (iCell = 0; iCell < StdI->NCell; iCell++) {
1776  StdFace_FindSite(StdI,
1777  StdI->Cell[iCell][0], StdI->Cell[iCell][1], StdI->Cell[iCell][2],
1778  StdI->Cell[dCell][0], StdI->Cell[dCell][1], StdI->Cell[dCell][2],
1779  isiteUC, jsiteUC, &isite, &jsite, &Cphase, dR);
1780 
1781  Jastrow[isite][jsite] = NJastrow;
1782  Jastrow[jsite][isite] = NJastrow;
1783 
1784  }/*for (iCell = 0; iCell < StdI->NCell; iCell++)*/
1785 
1786  NJastrow += 1;
1787 
1788  }/*for (jsiteUC = 0; jsiteUC < StdI->NsiteUC; jsiteUC++)*/
1789  }/*for (isiteUC = 0; isiteUC < StdI->NsiteUC; isiteUC++)*/
1790  }/*for (dCell = 0; dCell < StdI->NCell; dCell++)*/
1791  }/*if (strcmp(StdI->model, "spin") != 0)*/
1792  }/*if (abs(StdI->NMPTrans) != 1)*/
1793 
1794  fp = fopen("jastrowidx.def", "w");
1795  fprintf(fp, "=============================================\n");
1796  fprintf(fp, "NJastrowIdx %10d\n", NJastrow);
1797  fprintf(fp, "ComplexType %10d\n", 0);
1798  fprintf(fp, "=============================================\n");
1799  fprintf(fp, "=============================================\n");
1800 
1801  for (isite = 0; isite < StdI->nsite; isite++) {
1802  for (jsite = 0; jsite < StdI->nsite; jsite++) {
1803  if (isite == jsite) continue;
1804  fprintf(fp, "%5d %5d %5d\n", isite, jsite, Jastrow[isite][jsite]);
1805  }/*for (jsite = 0; jsite < isite; jsite++)*/
1806  }/*for (isite = 0; isite < StdI->nsite; isite++)*/
1807 
1808  for (iJastrow = 0; iJastrow < NJastrow; iJastrow++){
1809  if (strcmp(StdI->model, "hubbard") == 0 || iJastrow > 0)
1810  fprintf(fp, "%5d %5d\n", iJastrow, 1);
1811  else
1812  fprintf(fp, "%5d %5d\n", iJastrow, 0);
1813  }/*for (iJastrow = 0; iJastrow < NJastrow; iJastrow++)*/
1814  fflush(fp);
1815  fclose(fp);
1816  fprintf(stdout, " jastrowidx.def is written.\n");
1817 
1818  for (isite = 0; isite < StdI->nsite; isite++) free(Jastrow[isite]);
1819  free(Jastrow);
1820 }/*void PrintJastrow*/
1821 #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:200
StdIntList::NHund
int NHund
Number of Hund term, counted in each lattice file.
Definition: StdFace_vals.h:176
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:102
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:134
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:100
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:155
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:132
StdFace_PrintGeometry
void StdFace_PrintGeometry(struct StdIntList *StdI)
Print geometry of sites for the pos-process of correlation function.
Definition: StdFace_ModelUtil.c:1157
StdIntList::At
double ** At
[StdIntList::nt][3] Vector potential.
Definition: StdFace_vals.h:285
StdIntList::Ex
double * Ex
[StdIntList::NEx] Coefficient of exchange term, malloc in StdFace_MallocInteractions() and set in Std...
Definition: StdFace_vals.h:189
StdIntList::JpAll
double JpAll
Isotropic, diagonal spin coupling (2nd Near), input parameter Jp.
Definition: StdFace_vals.h:84
StdIntList::PLIndx
int ** PLIndx
[StdIntList::NPairLift][2] Site indices of pair-lift term, malloc in StdFace_MallocInteractions() and...
Definition: StdFace_vals.h:194
StdIntList::trans
double complex * trans
[StdIntList::ntrans] Coefficient of one-body term, malloc in StdFace_MallocInteractions() and set in ...
Definition: StdFace_vals.h:147
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:1198
StdIntList::JAll
double JAll
Isotropic, diagonal spin coupling (1st Near.), input parameter J.
Definition: StdFace_vals.h:82
StdIntList::S2
int S2
Total spin |S| of a local spin, input from file.
Definition: StdFace_vals.h:215
StdIntList::ExIndx
int ** ExIndx
[StdIntList::NEx][2] Site indices of exchange term, malloc in StdFace_MallocInteractions() and set in...
Definition: StdFace_vals.h:186
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_InputCoulombV
void StdFace_InputCoulombV(struct StdIntList *StdI, 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:1109
StdFace_InputHopp
void StdFace_InputHopp(struct StdIntList *StdI, 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:1134
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:143
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:135
StdIntList::CintraIndx
int ** CintraIndx
[StdIntList::NCintra][1] Site indices of intra-site Coulomb term, malloc in StdFace_MallocInteraction...
Definition: StdFace_vals.h:161
StdIntList::NCinter
int NCinter
Number of inter-site Coulomb interaction, counted in each lattice file.
Definition: StdFace_vals.h:167
StdIntList::Gamma
double Gamma
Transvars magnetic field, input parameter.
Definition: StdFace_vals.h:127
StdIntList::NCintra
int NCintra
Number of intra-site Coulomb interaction, counted in each lattice file.
Definition: StdFace_vals.h:158
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 (1nd Near), input parameter J0.
Definition: StdFace_vals.h:86
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:104
StdIntList::Hund
double * Hund
[StdIntList::NHund] Coefficient of Hund term, malloc in StdFace_MallocInteractions() and set in StdFa...
Definition: StdFace_vals.h:181
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:197
StdIntList::phase
double phase[3]
Boundary phase, input parameter phase0, etc.
Definition: StdFace_vals.h:133
StdIntList::pumpindx
int *** pumpindx
[StdIntList::nt][StdIntList::npump][4] Site/spin indices of one-body term, malloc in StdFace_MallocIn...
Definition: StdFace_vals.h:279
StdIntList::NPairLift
int NPairLift
Number of pair-lift term, counted in each lattice file.
Definition: StdFace_vals.h:192
StdIntList::PumpBody
int PumpBody
one- or two-body pumping, defined from StdIntList::PumpType
Definition: StdFace_vals.h:276
StdIntList::npump
int * npump
[StdIntList::nt] Number of transfer, counted in each lattice file.
Definition: StdFace_vals.h:278
StdIntList::locspinflag
int * locspinflag
[StdIntList::nsite] LocSpin in Expert mode, malloc and set in each lattice file.
Definition: StdFace_vals.h:141
StdIntList::HundIndx
int ** HundIndx
[StdIntList::NHund][2] Site indices of Hund term, malloc in StdFace_MallocInteractions() and set in S...
Definition: StdFace_vals.h:178
StdIntList::Lanczos_max
int Lanczos_max
The maxixmum number of iterations, input from file.
Definition: StdFace_vals.h:237
StdIntList::V
double V
Off-site Coulomb potential (1st), input parameter.
Definition: StdFace_vals.h:72
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:74
StdIntList::pump
double complex ** pump
[StdIntList::nt][StdIntList::npump] Coefficient of one-body term, malloc in StdFace_MallocInteraction...
Definition: StdFace_vals.h:282
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_InputSpin
void StdFace_InputSpin(struct StdIntList *StdI, double Jp[3][3], double JpAll, char *Jpname)
Input spin-spin interaction other than nearest-neighbor.
Definition: StdFace_ModelUtil.c:1060
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:184
StdIntList::method
char method[256]
The name of method, input from file.
Definition: StdFace_vals.h:232
StdIntList::intrindx
int ** intrindx
[StdIntList::nintr][8] Site/spin indices of two-body term, malloc in StdFace_MallocInteractions() and...
Definition: StdFace_vals.h:152
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
StdFace_InputSpinNN
void StdFace_InputSpinNN(struct StdIntList *StdI, double J0[3][3], double J0All, char *J0name)
Input nearest-neighbor spin-spin interaction.
Definition: StdFace_ModelUtil.c:973
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:205
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:140
StdIntList::transindx
int ** transindx
[StdIntList::ntrans][4] Site/spin indices of one-body term, malloc in StdFace_MallocInteractions() an...
Definition: StdFace_vals.h:144
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:164
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:938
StdIntList::nintr
int nintr
Number of InterAll, counted in each lattice file.
Definition: StdFace_vals.h:150
StdIntList::h
double h
Longitudinal magnetic field, input parameter.
Definition: StdFace_vals.h:126
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:202
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:173
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
StdIntList::CinterIndx
int ** CinterIndx
[StdIntList::NCinter][2] Site indices of inter-site Coulomb term, malloc in StdFace_MallocInteraction...
Definition: StdFace_vals.h:170