diff --git a/src/mVMC/average.c b/src/mVMC/average.c index 4cb0749e..f4954818 100644 --- a/src/mVMC/average.c +++ b/src/mVMC/average.c @@ -230,7 +230,7 @@ void WeightAverageGreenFunc(MPI_Comm comm) { } if(NLanczosMode>1){ /* QCisAjsQ and QCisAjsCktAltQ */ - n = NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*NCisAjsCktAltDC; + n = NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*(NCisAjsCktAltDC+NCisAjsCktAlt); if(AllComplexFlag==0){ vec_real=QCisAjsQ_real; weightAverageReduce_real(n,vec_real,comm); diff --git a/src/mVMC/calgrn.c b/src/mVMC/calgrn.c index 0ccd0dc9..708a186d 100644 --- a/src/mVMC/calgrn.c +++ b/src/mVMC/calgrn.c @@ -84,7 +84,7 @@ void CalculateGreenFunc(const double w, const double complex ip, int *eleIdx, in tmp = GreenFunc2(ri,rj,rk,rl,s,t,ip,myEleIdx,eleCfg,myEleNum,eleProjCnt, myProjCntNew,rbmCnt,myRBMCntNew,myBuffer); - PhysCisAjsCktAltDC[idx] += w*tmp; + LocalCisAjsCktAltDC[idx] = tmp; } #pragma omp master @@ -95,6 +95,11 @@ void CalculateGreenFunc(const double w, const double complex ip, int *eleIdx, in PhysCisAjs[idx] += w*LocalCisAjs[idx]; } + #pragma omp for private(idx) nowait + for (idx=0;idx, */ @@ -308,10 +308,13 @@ double *LSLQ_real; /* [NLSHam][NLSHam]*/ //TBC double complex *QCisAjsQ; /* QCisAjsQ[NLSHam][NLSHam][NCisAjs]*/ //TBC double complex *QCisAjsCktAltQ; /* QCisAjsCktAltQ[NLSHam][NLSHam][NCisAjsCktAlt]*/ //TBC +double complex *QCisAjsCktAltQDC; /* QCisAjsCktAltQ[NLSHam][NLSHam][NCisAjsCktAlt] + DC Lanczos Calculation */ double complex *LSLCisAjs; /* [NLSHam][NCisAjs]*/ //TBC double *QCisAjsQ_real; /* QCisAjsQ[NLSHam][NLSHam][NCisAjs]*/ //TBC double *QCisAjsCktAltQ_real; /* QCisAjsCktAltQ[NLSHam][NLSHam][NCisAjsCktAlt]*/ //TBC +double *QCisAjsCktAltQDC_real; /* QCisAjsCktAltQ[NLSHam][NLSHam][NCisAjsCktAlt]*/ double *LSLCisAjs_real; /* [NLSHam][NCisAjs]*/ //TBC /***** Output File *****/ @@ -327,8 +330,10 @@ FILE *FileLS; FILE *FileLSQQQQ; FILE *FileLSQCisAjsQ; FILE *FileLSQCisAjsCktAltQ; +FILE *FileLSQCisAjsCktAltQ; FILE *FileLSCisAjs; FILE *FileLSCisAjsCktAlt; +FILE *FileLSCisAjsCktAltDC; /* FILE *FileTimerList; */ diff --git a/src/mVMC/include/physcal_lanczos.h b/src/mVMC/include/physcal_lanczos.h index 9dd02cce..27dbf7ff 100644 --- a/src/mVMC/include/physcal_lanczos.h +++ b/src/mVMC/include/physcal_lanczos.h @@ -7,6 +7,7 @@ int PhysCalLanczos_real( double *_QQQQ_real, double *_QCisAjsQ_real, double *_QCisAjsCktAltQ_real, + double *_QCisAjsCktAltQDC_real, const int _nLSHam, const int _Ns, const int _nCisAjs, @@ -14,20 +15,23 @@ int PhysCalLanczos_real( int **_iOneBodyGIdx, int **_CisAjsLzIdx, const int _nCisAjsCktAlt, - int **_CisAjsCktAlt, + const int _nCisAjsCktAltDC, + int **_CisAjsCktAltDC, const int _NLanczosmode, FILE *_FileLS, FILE *_FileLSQQQQ, FILE *_FileLSQCisAjsQ, FILE *_FileLSQCisAjsCktAltQ, FILE *_FileLSCisAjs, - FILE *_FileLSCisAjsCktAlt + FILE *_FileLSCisAjsCktAlt, + FILE *_FileLSCisAjsCktAltDC ); int PhysCalLanczos_fcmp( double complex* _QQQQ, double complex* _QCisAjsQ, double complex* _QCisAjsCktAltQ, + double complex* _QCisAjsCktAltQDC, const int _nLSHam, const int _Ns, const int _nCisAjs, @@ -35,14 +39,16 @@ int PhysCalLanczos_fcmp( int **_iOneBodyGIdx, int **_CisAjsLzIdx, const int _nCisAjsCktAlt, - int **_CisAjsCktAlt, + const int _nCisAjsCktAltDC, + int **_CisAjsCktAltDC, const int _NLanczosmode, FILE *_FileLS, FILE *_FileLSQQQQ, FILE *_FileLSQCisAjsQ, FILE *_FileLSQCisAjsCktAltQ, FILE *_FileLSCisAjs, - FILE *_FileLSCisAjsCktAlt + FILE *_FileLSCisAjsCktAlt, + FILE *_FileLSCisAjsCktAltDC ); #endif diff --git a/src/mVMC/initfile.c b/src/mVMC/initfile.c index 815f971f..80a04f38 100644 --- a/src/mVMC/initfile.c +++ b/src/mVMC/initfile.c @@ -126,9 +126,15 @@ void InitFilePhysCal(int i, int rank) { CDataFileHead, idx); FileLSCisAjs = fopen(fileName, "w"); - sprintf(fileName, "%s_ls_cisajscktalt_%03d.dat", + // CACA + sprintf(fileName, "%s_ls_cisajscktaltex_%03d.dat", CDataFileHead, idx); FileLSCisAjsCktAlt = fopen(fileName, "w"); + + // CACADC + sprintf(fileName, "%s_ls_cisajscktalt_%03d.dat", + CDataFileHead, idx); + FileLSCisAjsCktAltDC = fopen(fileName, "w"); } } @@ -176,6 +182,7 @@ void CloseFilePhysCal(int rank) { #endif fclose(FileLSCisAjs); fclose(FileLSCisAjsCktAlt); + fclose(FileLSCisAjsCktAltDC); } } diff --git a/src/mVMC/physcal_lanczos.c b/src/mVMC/physcal_lanczos.c index 26b9782a..13640256 100644 --- a/src/mVMC/physcal_lanczos.c +++ b/src/mVMC/physcal_lanczos.c @@ -29,6 +29,7 @@ int PhysCalLanczos_real double *_QQQQ_real, double *_QCisAjsQ_real, double *_QCisAjsCktAltQ_real, + double *_QCisAjsCktAltQDC_real, const int _nLSHam, const int _Ns, const int _nCisAjs, @@ -36,14 +37,16 @@ int PhysCalLanczos_real int **_iOneBodyGIdx, int **_CisAjsLzIdx, const int _nCisAjsCktAlt, - int **_CisAjsCktAlt, + const int _nCisAjsCktAltDC, + int **_CisAjsCktAltDC, const int _NLanczosmode, FILE *_FileLS, FILE *_FileLSQQQQ, FILE *_FileLSQCisAjsQ, FILE *_FileLSQCisAjsCktAltQ, FILE *_FileLSCisAjs, - FILE *_FileLSCisAjsCktAlt + FILE *_FileLSCisAjsCktAlt, + FILE *_FileLSCisAjsCktAltDC ) { int i=0; @@ -53,8 +56,10 @@ int PhysCalLanczos_real double alpha_m, ene_m, ene_vm; double *LS_CisAjs_real; double *LS_CisAjsCktAlt_real; + double *LS_CisAjsCktAltDC_real; LS_CisAjs_real = (double*)malloc(sizeof(double)*_nCisAjs); LS_CisAjsCktAlt_real = (double*)malloc(sizeof(double)*_nCisAjsCktAlt); + LS_CisAjsCktAltDC_real = (double*)malloc(sizeof(double)*_nCisAjsCktAltDC); CalculateEne(_QQQQ_real[2], _QQQQ_real[3], _QQQQ_real[10], _QQQQ_real[11], _QQQQ_real[15], &alpha_p, &ene_p, &ene_vp, &alpha_m, &ene_m, &ene_vm); @@ -114,21 +119,30 @@ int PhysCalLanczos_real fprintf(_FileLSQCisAjsCktAltQ, "\n"); #endif + CalculatePhysVal_real(_QQQQ_real[2], _QQQQ_real[3], + alpha, _QCisAjsCktAltQDC_real, _nCisAjsCktAltDC, + _nLSHam, LS_CisAjsCktAltDC_real); + /* zvo_ls_cisajscktalt.dat */ + for (i = 0; i < _nCisAjsCktAltDC; i++) { + fprintf(_FileLSCisAjsCktAltDC, "%d %d %d %d %d %d %d %d % .18e 0.0\n", + _CisAjsCktAltDC[i][0], _CisAjsCktAltDC[i][1], _CisAjsCktAltDC[i][2], _CisAjsCktAltDC[i][3], + _CisAjsCktAltDC[i][4], _CisAjsCktAltDC[i][5], _CisAjsCktAltDC[i][6], _CisAjsCktAltDC[i][7], + LS_CisAjsCktAltDC_real[i]); + } + fprintf(_FileLSCisAjsCktAltDC, "\n"); + CalculatePhysVal_real(_QQQQ_real[2], _QQQQ_real[3], alpha, _QCisAjsCktAltQ_real, _nCisAjsCktAlt, _nLSHam, LS_CisAjsCktAlt_real); - /* zvo_ls_cisajscktalt.dat */ + /* zvo_ls_cisajscktaltex.dat */ for (i = 0; i < _nCisAjsCktAlt; i++) { - fprintf(_FileLSCisAjsCktAlt, "%d %d %d %d %d %d %d %d % .18e 0.0\n", - _CisAjsCktAlt[i][0], _CisAjsCktAlt[i][1], _CisAjsCktAlt[i][2], _CisAjsCktAlt[i][3], - _CisAjsCktAlt[i][4], _CisAjsCktAlt[i][5], _CisAjsCktAlt[i][6], _CisAjsCktAlt[i][7], - LS_CisAjsCktAlt_real[i]); - + fprintf(_FileLSCisAjsCktAlt, "% .18e 0.0 ", LS_CisAjsCktAlt_real[i]); // LS_CisAjsCktAlt_real } fprintf(_FileLSCisAjsCktAlt, "\n"); } free(LS_CisAjs_real); free(LS_CisAjsCktAlt_real); + free(LS_CisAjsCktAltDC_real); return 0; } @@ -136,6 +150,7 @@ int PhysCalLanczos_fcmp( double complex* _QQQQ, double complex* _QCisAjsQ, double complex* _QCisAjsCktAltQ, + double complex* _QCisAjsCktAltQDC, const int _nLSHam, const int _Ns, const int _nCisAjs, @@ -143,14 +158,16 @@ int PhysCalLanczos_fcmp( int **_iOneBodyGIdx, int **_CisAjsLzIdx, const int _nCisAjsCktAlt, - int **_CisAjsCktAlt, + const int _nCisAjsCktAltDC, + int **_CisAjsCktAltDC, const int _NLanczosmode, FILE *_FileLS, FILE *_FileLSQQQQ, FILE *_FileLSQCisAjsQ, FILE *_FileLSQCisAjsCktAltQ, FILE *_FileLSCisAjs, - FILE *_FileLSCisAjsCktAlt + FILE *_FileLSCisAjsCktAlt, + FILE *_FileLSCisAjsCktAltDC ) { int i=0; @@ -160,8 +177,10 @@ int PhysCalLanczos_fcmp( double alpha_m, ene_m, ene_vm; double complex*LS_CisAjs; double complex*LS_CisAjsCktAlt; + double complex*LS_CisAjsCktAltDC; LS_CisAjs = (double complex*)malloc(sizeof(double complex)*_nCisAjs); LS_CisAjsCktAlt = (double complex*)malloc(sizeof(double complex)*_nCisAjsCktAlt); + LS_CisAjsCktAltDC = (double complex*)malloc(sizeof(double complex)*_nCisAjsCktAltDC); /* zvo_ls.dat */ if(!CalculateEne(creal(_QQQQ[2]),creal(_QQQQ[3]), @@ -221,20 +240,30 @@ int PhysCalLanczos_fcmp( } fprintf(_FileLSQCisAjsCktAltQ, "\n"); #endif + CalculatePhysVal_fcmp(_QQQQ[2], _QQQQ[3], + alpha, _QCisAjsCktAltQDC, _nCisAjsCktAltDC, + _nLSHam, LS_CisAjsCktAltDC); + /* zvo_ls_cisajscktalt.dat */ + for (i = 0; i < _nCisAjsCktAltDC; i++) { + fprintf(_FileLSCisAjsCktAltDC, "%d %d %d %d %d %d %d %d % .18e % .18e\n", + _CisAjsCktAltDC[i][0], _CisAjsCktAltDC[i][1], _CisAjsCktAltDC[i][2], _CisAjsCktAltDC[i][3], + _CisAjsCktAltDC[i][4], _CisAjsCktAltDC[i][5], _CisAjsCktAltDC[i][6], _CisAjsCktAltDC[i][7], + creal(LS_CisAjsCktAltDC[i]), cimag(LS_CisAjsCktAltDC[i])); + } + fprintf(_FileLSCisAjsCktAltDC, "\n"); + CalculatePhysVal_fcmp(_QQQQ[2], _QQQQ[3], alpha, _QCisAjsCktAltQ, _nCisAjsCktAlt, _nLSHam, LS_CisAjsCktAlt); - /* zvo_ls_cisajs.dat */ + /* zvo_ls_cisajscktaltex.dat */ for (i = 0; i < _nCisAjsCktAlt; i++) { - fprintf(_FileLSCisAjsCktAlt, "%d %d %d %d %d %d %d %d % .18e % .18e\n", - _CisAjsCktAlt[i][0], _CisAjsCktAlt[i][1], _CisAjsCktAlt[i][2], _CisAjsCktAlt[i][3], - _CisAjsCktAlt[i][4], _CisAjsCktAlt[i][5], _CisAjsCktAlt[i][6], _CisAjsCktAlt[i][7], - creal(LS_CisAjsCktAlt[i]), cimag(LS_CisAjsCktAlt[i])); + fprintf(_FileLSCisAjsCktAlt, "% .18e % .18e ", creal(LS_CisAjsCktAlt[i]), cimag(LS_CisAjsCktAlt[i])); // LS_CisAjsCktAlt_real } fprintf(_FileLSCisAjsCktAlt, "\n"); } free(LS_CisAjs); free(LS_CisAjsCktAlt); + free(LS_CisAjsCktAltDC); return 0; } diff --git a/src/mVMC/readdef.c b/src/mVMC/readdef.c index 9e2b0818..d19c96e6 100644 --- a/src/mVMC/readdef.c +++ b/src/mVMC/readdef.c @@ -88,10 +88,10 @@ int GetInfoInterAll(FILE *fp, int **ArrayIdx, double complex *ArrayValue, int GetInfoOptTrans(FILE *fp, int **Array, double *ArrayPara, int *ArrayOpt, int **ArraySgn, int _iFlagOptTrans, int *iOptCount, int _fidx, int _APFlag, int Nsite, int NArray, char *defname); -int GetInfoTwoBodyG(FILE *fp, int **ArrayIdx, int **ArrayIdxTwoBodyGLz, int **ArrayToIdx, int **ArrayIdxOneBodyG, - int _NLanczosMode, int Nsite, int NArray, char *defname); +int GetInfoTwoBodyG(FILE *fp, int **ArrayIdx, int Nsite, int NArray, char *defname); -int GetInfoTwoBodyGEx(FILE *fp, int **ArrayIdx, int Nsite, int NArray, char *defname); +int GetInfoTwoBodyGEx(FILE *fp, int **ArrayIdx, int **ArrayToIdx, int **ArrayIdxOneBodyG, + int Nsite, int NArray, char *defname); int GetInfoOrbitalGeneral(FILE *fp, int **Array, int *ArrayOpt, int **ArraySgn, int *iOptCount, int _fidx, int _iComplexFlag, int _iFlagOrbitalGeneral, int _APFlag, int Nsite, int NArray, @@ -191,7 +191,7 @@ int JudgeOrbitalMode(int *_iFlgOrbitalGeneral, const int _iFlgOrbitalAP, const i return iret; } -int ReadGreen(char *xNameListFile, int Nca, int **caIdx, int Ncacadc, int **cacaDCIdx, int Ns) { +int ReadGreen(char *xNameListFile, int Nca, int **caIdx, int Ncacadc, int **cacaIdx, int Ns) { FILE *fp; char defname[D_FileNameMax]; char ctmp[D_FileNameMax]; @@ -227,9 +227,10 @@ int ReadGreen(char *xNameListFile, int Nca, int **caIdx, int Ncacadc, int **caca /*cisajs.def----------------------------------------*/ if (GetInfoOneBodyG(fp, caIdx, iOneBodyGIdx, 0, Ns, Nca, defname) != 0) return (-1); break; - case KWTwoBodyG: - /*cisajscktaltdc.def--------------------------------*/ - if (GetInfoTwoBodyG(fp, cacaDCIdx, cacaDCIdx, iOneBodyGIdx, caIdx, 0, Ns, Ncacadc, defname) != 0) return (-1); + case KWTwoBodyGEx: + /*cisajscktalt.def----------------------------------*/ + /*load as if it's DC for index rearranging----------*/ + if (GetInfoTwoBodyG(fp, cacaIdx, Ns, Ncacadc, defname) != 0) return (-1); break; default: break; @@ -245,27 +246,34 @@ int ReadGreen(char *xNameListFile, int Nca, int **caIdx, int Ncacadc, int **caca /// \param Ncacadc Number of CisAjsCktAltDC /// \param Ns Number of sites /// \return Number of calculation target -int CountOneBodyGForLanczos(char *xNameListFile, int Nca, int Ncacadc, int Ns, int **caIdx, int **iFlgOneBodyG) { +int CountOneBodyGForLanczos(char *xNameListFile, int Nca, int Ncacadc, int Ns, int **iFlgOneBodyG) { int info = 0; int i, j, isite1, isite2; int icount = 0; - int **cacaDCIdx; + int **cacaIdx; + int **caIdx; - cacaDCIdx = malloc(sizeof(int *) * Ncacadc); - //pInt=cacaDCIdx[0]; - for (i = 0; i < Ncacadc; i++) { - cacaDCIdx[i] = malloc(sizeof(int) * 8); - } + cacaIdx = malloc(sizeof(int *) * Ncacadc); + for (i = 0; i < Ncacadc; i++) + cacaIdx[i] = malloc(sizeof(int) * 8); + caIdx = malloc(sizeof(int *) * Nca); + for (i = 0; i < Nca; i++) + caIdx[i] = malloc(sizeof(int) * 4); for (i = 0; i < 2 * Ns; i++) { for (j = 0; j < 2 * Ns; j++) { iFlgOneBodyG[i][j] = -1; } } - info = ReadGreen(xNameListFile, Nca, caIdx, Ncacadc, cacaDCIdx, Ns); + info = ReadGreen(xNameListFile, Nca, caIdx, Ncacadc, cacaIdx, Ns); if (info != 0) { - free(cacaDCIdx); + for (i = 0; i < Ncacadc; i++) + free(cacaIdx[i]); + free(cacaIdx); + for (i = 0; i < Nca; i++) + free(caIdx[i]); + free(caIdx); return (info); } @@ -279,8 +287,8 @@ int CountOneBodyGForLanczos(char *xNameListFile, int Nca, int Ncacadc, int Ns, i } //cisajscktalt -> cisajs, cltakt (Note: indecies of the latter Green's function are modified) for (i = 0; i < Ncacadc; i++) { - isite1 = cacaDCIdx[i][0] + cacaDCIdx[i][1] * Ns; - isite2 = cacaDCIdx[i][2] + cacaDCIdx[i][3] * Ns; + isite1 = cacaIdx[i][0] + cacaIdx[i][1] * Ns; + isite2 = cacaIdx[i][2] + cacaIdx[i][3] * Ns; if (iFlgOneBodyG[isite1][isite2] == -1) { iFlgOneBodyG[isite1][isite2] = icount; icount++; @@ -290,15 +298,20 @@ int CountOneBodyGForLanczos(char *xNameListFile, int Nca, int Ncacadc, int Ns, i isite1 = cacaDCIdx[i][4] + cacaDCIdx[i][5] * Ns; isite2 = cacaDCIdx[i][6] + cacaDCIdx[i][7] * Ns; */ - isite1 = cacaDCIdx[i][6] + cacaDCIdx[i][7] * Ns; - isite2 = cacaDCIdx[i][4] + cacaDCIdx[i][5] * Ns; + isite1 = cacaIdx[i][6] + cacaIdx[i][7] * Ns; + isite2 = cacaIdx[i][4] + cacaIdx[i][5] * Ns; if (iFlgOneBodyG[isite1][isite2] == -1) { iFlgOneBodyG[isite1][isite2] = icount; icount++; } } - free(cacaDCIdx); + for (i = 0; i < Ncacadc; i++) + free(cacaIdx[i]); + free(cacaIdx); + for (i = 0; i < Nca; i++) + free(caIdx[i]); + free(caIdx); return icount; } @@ -559,23 +572,17 @@ int ReadDefFileNInt(char *xNameListFile, MPI_Comm comm) { //TODO: LanczosMode is not supported for Sz not conserved mode. - //For Lanczos mode: Calculation of Green's function - if (bufInt[IdxLanczosMode] > 1) { - //Get info of CisAjs and CisAjsCktAltDC + //For indirect calculation of Green's function + if (bufInt[IdxNTwoBodyGEx] > 0 || bufInt[IdxLanczosMode] > 1) { + //Get info of CisAjs and CisAjsCktAlt(GreenTwoEx as if it's DC) int i; - NCisAjsLz = bufInt[IdxNOneBodyG]; - //bufInt[IdxNTwoBodyGEx] = bufInt[IdxNTwoBodyG]; - CisAjsLzIdx = malloc(sizeof(int *) * NCisAjsLz); - for (i = 0; i < NCisAjsLz; i++) { - CisAjsLzIdx[i] = malloc(sizeof(int) * 4); - } iOneBodyGIdx = malloc(sizeof(int *) * (2 * bufInt[IdxNsite])); //For spin - //pInt=iFlgOneBodyG[0]; for (i = 0; i < 2 * bufInt[IdxNsite]; i++) { iOneBodyGIdx[i] = malloc(sizeof(int) * (2 * bufInt[IdxNsite])); } - bufInt[IdxNOneBodyG] = CountOneBodyGForLanczos(xNameListFile, NCisAjsLz, bufInt[IdxNTwoBodyG], - bufInt[IdxNsite], CisAjsLzIdx, iOneBodyGIdx); + bufInt[IdxNOneBodyG] = CountOneBodyGForLanczos(xNameListFile, + bufInt[IdxNOneBodyG], bufInt[IdxNTwoBodyGEx], + bufInt[IdxNsite], iOneBodyGIdx); } //CalcNCond @@ -789,7 +796,7 @@ int ReadDefFileNInt(char *xNameListFile, MPI_Comm comm) { + Nsite * NQPTrans /* QPTransInv */ + Nsite * NQPTrans /* QPTransSgn */ + 4 * NCisAjs /* CisAjs */ - + 8 * NCisAjsCktAlt /* CisAjsCktAlt */ + + 2 * NCisAjsCktAlt /* CisAjsCktAlt */ + 8 * NCisAjsCktAltDC /* CisAjsCktAltDC */ + 8 * NInterAll /* InterAll */ + Nsite * NQPOptTrans /* QPOptTrans */ @@ -849,8 +856,6 @@ int ReadDefFileIdxPara(char *xNameListFile, MPI_Comm comm) { int x0, x1; int rank; - int iNOneBodyG; - MPI_Comm_rank(comm, &rank); if (rank == 0) { @@ -1038,19 +1043,17 @@ int ReadDefFileIdxPara(char *xNameListFile, MPI_Comm comm) { case KWOneBodyG: /*cisajs.def----------------------------------------*/ - iNOneBodyG = (NLanczosMode < 2) ? NCisAjs : NCisAjsLz; - if (GetInfoOneBodyG(fp, CisAjsIdx, iOneBodyGIdx, NLanczosMode, Nsite, iNOneBodyG, defname) != 0) info = 1; + if (GetInfoOneBodyG(fp, CisAjsIdx, iOneBodyGIdx, NCisAjsCktAlt>0, Nsite, NCisAjs, defname) != 0) info = 1; break; case KWTwoBodyGEx: /*cisajscktalt.def----------------------------------*/ - if (GetInfoTwoBodyGEx(fp, CisAjsCktAltIdx, Nsite, NCisAjsCktAlt, defname) != 0) info = 1; + if (GetInfoTwoBodyGEx(fp, CisAjsCktAltIdx, iOneBodyGIdx, CisAjsIdx, Nsite, NCisAjsCktAlt, defname) != 0) info = 1; break; case KWTwoBodyG: /*cisajscktaltdc.def--------------------------------*/ - if (GetInfoTwoBodyG(fp, CisAjsCktAltDCIdx, CisAjsCktAltLzIdx, iOneBodyGIdx, CisAjsIdx, NLanczosMode, Nsite, - NCisAjsCktAltDC, defname) != 0) + if (GetInfoTwoBodyG(fp, CisAjsCktAltDCIdx, Nsite, NCisAjsCktAltDC, defname) != 0) info = 1; break; @@ -1157,9 +1160,6 @@ int ReadDefFileIdxPara(char *xNameListFile, MPI_Comm comm) { */ #ifdef _mpi_use SafeMpiBcastInt(LocSpn, NTotalDefInt, comm); - if (NLanczosMode > 1) { - SafeMpiBcastInt(CisAjsCktAltLzIdx[0], NCisAjsCktAltDC * 2, comm); - } SafeMpiBcast_fcmp(ParaTransfer, NTransfer + NInterAll, comm); SafeMpiBcast(ParaCoulombIntra, NTotalDefDouble, comm); SafeMpiBcast_fcmp(ParaQPTrans, NQPTrans, comm); @@ -2265,7 +2265,7 @@ int GetInfoTransSym(FILE *fp, int **Array, int **ArraySgn, int **ArrayInv, doubl } int -GetInfoOneBodyG(FILE *fp, int **ArrayIdx, int **ArrayToIdx, int _NLanczosMode, int Nsite, int NArray, char *defname) { +GetInfoOneBodyG(FILE *fp, int **ArrayIdx, int **ArrayToIdx, int IndirectGFOn, int Nsite, int NArray, char *defname) { char ctmp2[256]; int idx = 0, info = 0; int x0 = 0, x1 = 0, x2 = 0, x3 = 0; @@ -2274,12 +2274,12 @@ GetInfoOneBodyG(FILE *fp, int **ArrayIdx, int **ArrayToIdx, int _NLanczosMode, i if (NArray == 0) return 0; while (fgets(ctmp2, sizeof(ctmp2) / sizeof(char), fp) != NULL) { sscanf(ctmp2, "%d %d %d %d\n", &x0, &x1, &x2, &x3); - if (_NLanczosMode < 2) { // Normal + if (!IndirectGFOn) { // Normal ArrayIdx[idx][0] = x0; ArrayIdx[idx][1] = x1; ArrayIdx[idx][2] = x2; ArrayIdx[idx][3] = x3; - } else { //For Calc Green func by Lanczos mode + } else { //For Calc Green func indirectly isite1 = x0 + x1 * Nsite; isite2 = x2 + x3 * Nsite; idxLanczos = ArrayToIdx[isite1][isite2]; @@ -2296,28 +2296,50 @@ GetInfoOneBodyG(FILE *fp, int **ArrayIdx, int **ArrayToIdx, int _NLanczosMode, i idx++; } - if (idx != NArray) info = ReadDefFileError(defname); + if (!IndirectGFOn && idx != NArray) + info = ReadDefFileError(defname); return info; } -int GetInfoTwoBodyGEx(FILE *fp, int **ArrayIdx, int Nsite, int NArray, char *defname) { +// Formerly CisAjsCktAlt +int GetInfoTwoBodyGEx(FILE *fp, int **ArrayIdx, int **ArrayToIdx, int **ArrayIdxOneBodyG, + int Nsite, int NArray, char *defname) { char ctmp2[256]; int idx = 0, info = 0; int x0 = 0, x1 = 0, x2 = 0, x3 = 0; int x4 = 0, x5 = 0, x6 = 0, x7 = 0; - //Debug + int isite1 = 0, isite2 = 0; + int idxLanczos = 0; + if (NArray == 0) return 0; while (fgets(ctmp2, sizeof(ctmp2) / sizeof(char), fp) != NULL) { sscanf(ctmp2, "%d %d %d %d %d %d %d %d\n", &x0, &x1, &x2, &x3, &x4, &x5, &x6, &x7); - ArrayIdx[idx][0] = x0; // Index to OneBodyG1 - ArrayIdx[idx][1] = x1; // Index to OneBodyG2 - ArrayIdx[idx][2] = x2; // G1:site i - ArrayIdx[idx][3] = x3; // G1:site j - ArrayIdx[idx][4] = x4; // G1:sigma1 - ArrayIdx[idx][5] = x5; // G2:site l - ArrayIdx[idx][6] = x6; // G2:site k - ArrayIdx[idx][7] = x7; // G2:sigma2 - if (CheckQuadSite(x2, x3, x5, x6, Nsite) != 0) { + + isite1 = x0 + x1 * Nsite; + isite2 = x2 + x3 * Nsite; + idxLanczos = ArrayToIdx[isite1][isite2]; + ArrayIdxOneBodyG[idxLanczos][0] = x0; + ArrayIdxOneBodyG[idxLanczos][1] = x1; + ArrayIdxOneBodyG[idxLanczos][2] = x2; + ArrayIdxOneBodyG[idxLanczos][3] = x3; + ArrayIdx[idx][0] = idxLanczos; + + isite1 = x6 + x7 * Nsite; + isite2 = x4 + x5 * Nsite; + idxLanczos = ArrayToIdx[isite1][isite2]; + /* + ArrayIdxOneBodyG[idxLanczos][0] = x4; + ArrayIdxOneBodyG[idxLanczos][1] = x5; + ArrayIdxOneBodyG[idxLanczos][2] = x6; + ArrayIdxOneBodyG[idxLanczos][3] = x7; + */ + ArrayIdxOneBodyG[idxLanczos][0] = x6; + ArrayIdxOneBodyG[idxLanczos][1] = x7; + ArrayIdxOneBodyG[idxLanczos][2] = x4; + ArrayIdxOneBodyG[idxLanczos][3] = x5; + ArrayIdx[idx][1] = idxLanczos; + + if (CheckQuadSite(x0, x2, x4, x6, Nsite) != 0) { fprintf(stderr, "Error: Site index is incorrect. \n"); info = 1; break; @@ -2328,8 +2350,8 @@ int GetInfoTwoBodyGEx(FILE *fp, int **ArrayIdx, int Nsite, int NArray, char *def return info; } -int GetInfoTwoBodyG(FILE *fp, int **ArrayIdx, int **ArrayIdxTwoBodyGLz, int **ArrayToIdx, int **ArrayIdxOneBodyG, - int _NLanczosMode, int Nsite, int NArray, char *defname) { +// Formerly CisAjsCktAltDC +int GetInfoTwoBodyG(FILE *fp, int **ArrayIdx, int Nsite, int NArray, char *defname) { char ctmp2[256]; int idx = 0, info = 0; int x0 = 0, x1 = 0, x2 = 0, x3 = 0; @@ -2354,44 +2376,6 @@ int GetInfoTwoBodyG(FILE *fp, int **ArrayIdx, int **ArrayIdxTwoBodyGLz, int **Ar break; } - if (_NLanczosMode > 1) { //Calc TwoBodyG by Lanczos method - - isite1 = x0 + x1 * Nsite; - isite2 = x2 + x3 * Nsite; - idxLanczos = ArrayToIdx[isite1][isite2]; - - ArrayIdxOneBodyG[idxLanczos][0] = x0; - ArrayIdxOneBodyG[idxLanczos][1] = x1; - ArrayIdxOneBodyG[idxLanczos][2] = x2; - ArrayIdxOneBodyG[idxLanczos][3] = x3; - /* - ArrayIdxOneBodyG[idxLanczos][0] = x2; - ArrayIdxOneBodyG[idxLanczos][1] = x3; - ArrayIdxOneBodyG[idxLanczos][2] = x0; - ArrayIdxOneBodyG[idxLanczos][3] = x1; - */ - ArrayIdxTwoBodyGLz[idx][0] = idxLanczos; - - /* - isite1 = x4 + x5 * Nsite; - isite2 = x6 + x7 * Nsite; - */ - isite1 = x6 + x7 * Nsite; - isite2 = x4 + x5 * Nsite; - idxLanczos = ArrayToIdx[isite1][isite2]; - ArrayIdxOneBodyG[idxLanczos][0] = x6; - ArrayIdxOneBodyG[idxLanczos][1] = x7; - ArrayIdxOneBodyG[idxLanczos][2] = x4; - ArrayIdxOneBodyG[idxLanczos][3] = x5; - - /* - ArrayIdxOneBodyG[idxLanczos][0] = x4; - ArrayIdxOneBodyG[idxLanczos][1] = x5; - ArrayIdxOneBodyG[idxLanczos][2] = x6; - ArrayIdxOneBodyG[idxLanczos][3] = x7; - */ - ArrayIdxTwoBodyGLz[idx][1] = idxLanczos; - } idx++; } if (idx != NArray) info = ReadDefFileError(defname); diff --git a/src/mVMC/setmemory.c b/src/mVMC/setmemory.c index 57735653..21652546 100644 --- a/src/mVMC/setmemory.c +++ b/src/mVMC/setmemory.c @@ -197,7 +197,7 @@ void SetMemoryDef() { CisAjsCktAltIdx = (int**)malloc(sizeof(int*)*NCisAjsCktAlt); for(i=0;i1){ - CisAjsCktAltLzIdx = malloc(sizeof(int*)*NCisAjsCktAltDC); - for(i=0;i0) { + for(i=0;i<2*Nsite;i++) + free(iOneBodyGIdx[i]); + free(iOneBodyGIdx); + } free(QPOptTransSgn); free(QPOptTrans); @@ -277,6 +275,7 @@ void FreeMemoryDef() { free(DoublonHolon4siteIdx); free(DoublonHolon2siteIdx); free(JastrowIdx); + free(ParaTransfer); free(ExchangeCoupling); free(PairHopping); free(HundCoupling); @@ -428,10 +427,11 @@ void SetMemory() { /***** Physical Quantity *****/ if(NVMCCalMode==1){ PhysCisAjs = (double complex*)malloc(sizeof(double complex) - *(NCisAjs+NCisAjsCktAlt+NCisAjsCktAltDC+NCisAjs)); + *(NCisAjs+NCisAjsCktAlt+NCisAjsCktAltDC+NCisAjs+NCisAjsCktAltDC)); PhysCisAjsCktAlt = PhysCisAjs + NCisAjs; PhysCisAjsCktAltDC = PhysCisAjsCktAlt + NCisAjsCktAlt; LocalCisAjs = PhysCisAjsCktAltDC + NCisAjsCktAltDC; + LocalCisAjsCktAltDC = LocalCisAjs + NCisAjs; if(NLanczosMode>0){ QQQQ = (double complex*)malloc(sizeof(double complex) @@ -444,14 +444,16 @@ void SetMemory() { if(NLanczosMode>1){ QCisAjsQ = (double complex*)malloc(sizeof(double complex) - *(NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*NCisAjsCktAltDC + NLSHam*NCisAjs) ); + *(NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*(NCisAjsCktAltDC+NCisAjsCktAlt) + NLSHam*NCisAjs) ); QCisAjsCktAltQ = QCisAjsQ + NLSHam*NLSHam*NCisAjs; - LSLCisAjs = QCisAjsCktAltQ + NLSHam*NLSHam*NCisAjsCktAltDC; + QCisAjsCktAltQDC = QCisAjsCktAltQ + NLSHam*NLSHam*NCisAjsCktAlt; + LSLCisAjs = QCisAjsCktAltQDC + NLSHam*NLSHam*NCisAjsCktAltDC; //for real QCisAjsQ_real = (double *)malloc(sizeof(double ) - *(NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*NCisAjsCktAltDC + NLSHam*NCisAjs) ); + *(NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*(NCisAjsCktAltDC+NCisAjsCktAlt) + NLSHam*NCisAjs) ); QCisAjsCktAltQ_real = QCisAjsQ_real + NLSHam*NLSHam*NCisAjs; - LSLCisAjs_real = QCisAjsCktAltQ_real + NLSHam*NLSHam*NCisAjsCktAltDC; + QCisAjsCktAltQDC_real = QCisAjsCktAltQ_real + NLSHam*NLSHam*NCisAjsCktAlt; + LSLCisAjs_real = QCisAjsCktAltQDC_real + NLSHam*NLSHam*NCisAjsCktAltDC; } } diff --git a/src/mVMC/vmccal.c b/src/mVMC/vmccal.c index 195e070a..50a3fb39 100644 --- a/src/mVMC/vmccal.c +++ b/src/mVMC/vmccal.c @@ -69,6 +69,16 @@ void calculateQCACAQ_real(double *qcacaq, const double *lslca, const double w, const int nLSHam, const int nCA, const int nCACA, int **cacaIdx); +void calculateQCACAQDC_real(double *qcacaq, const double *lslq, const double w, + const int nLSHam, const int nCA, const int nCACA, + int *eleIdx, int *eleCfg, int *eleNum, int *eleProjCnt, + const double h1, const double ip); + +void calculateQCACAQDC(double complex *qcacaq, const double complex *lslq, const double w, + const int nLSHam, const int nCA, const int nCACA, + int *eleIdx, int *eleCfg, int *eleNum, int *eleProjCnt, + const double complex h1, const double complex ip,double complex *rbmCnt); + void VMCMainCal(MPI_Comm comm) { int *eleIdx,*eleCfg,*eleNum,*eleProjCnt; double complex e,ip; @@ -276,14 +286,18 @@ void VMCMainCal(MPI_Comm comm) { LSLocalCisAjs_real(creal(e),creal(ip),eleIdx,eleCfg,eleNum,eleProjCnt); calculateQCAQ_real(QCisAjsQ_real,LSLCisAjs_real,LSLQ_real,w,NLSHam,NCisAjs); calculateQCACAQ_real(QCisAjsCktAltQ_real,LSLCisAjs_real,w,NLSHam,NCisAjs, - NCisAjsCktAltDC, CisAjsCktAltLzIdx); + NCisAjsCktAlt, CisAjsCktAltIdx); + calculateQCACAQDC_real(QCisAjsCktAltQDC_real,LSLQ_real,w,NLSHam,NCisAjs, + NCisAjsCktAltDC,eleIdx,eleCfg,eleNum,eleProjCnt,creal(e),creal(ip)); } else{ LSLocalCisAjs(e,ip,eleIdx,eleCfg,eleNum,eleProjCnt,rbmCnt); calculateQCAQ(QCisAjsQ,LSLCisAjs,LSLQ,w,NLSHam,NCisAjs); calculateQCACAQ(QCisAjsCktAltQ,LSLCisAjs,w,NLSHam,NCisAjs, - NCisAjsCktAltDC,CisAjsCktAltLzIdx); + NCisAjsCktAlt,CisAjsCktAltIdx); + calculateQCACAQDC(QCisAjsCktAltQDC,LSLQ,w,NLSHam,NCisAjs, + NCisAjsCktAltDC,eleIdx,eleCfg,eleNum,eleProjCnt,e,ip,rbmCnt); } StopTimer(44); } @@ -517,12 +531,12 @@ for(i=0;i1) { /* QCisAjsQ, QCisAjsCktAltQ, LSLCisAjs */ //[TODO]: Check the value n - n = NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*NCisAjsCktAltDC + n = NLSHam*NLSHam*NCisAjs + NLSHam*NLSHam*(NCisAjsCktAltDC+NCisAjsCktAlt) + NLSHam*NCisAjs; vec = QCisAjsQ; #pragma omp parallel for default(shared) private(i) for(i=0;i 0) { - if(NLanczosMode <2) { - for (i = 0; i < NCisAjs; i++) { - fprintf(FileCisAjs, "%d %d %d %d % .18e % .18e \n", CisAjsIdx[i][0], CisAjsIdx[i][1], CisAjsIdx[i][2], - CisAjsIdx[i][3], creal(PhysCisAjs[i]), cimag(PhysCisAjs[i])); - } - } - else{ - int idx=0; - for (i = 0; i < NCisAjsLz; i++) { - idx = iOneBodyGIdx[CisAjsLzIdx[i][0] + CisAjsLzIdx[i][1] * Nsite][CisAjsLzIdx[i][2] + - CisAjsLzIdx[i][3] * Nsite]; - //fprintf(stdout, "Debug: idx= %d value= % .18e % .18e\n", idx, creal(PhysCisAjs[idx]), cimag(PhysCisAjs[idx])); - fprintf(FileCisAjs, "%d %d %d %d % .18e % .18e \n", CisAjsLzIdx[idx][0], CisAjsLzIdx[idx][1], - CisAjsLzIdx[idx][2], CisAjsLzIdx[idx][3], creal(PhysCisAjs[idx]), cimag(PhysCisAjs[idx])); - } - } + for (i = 0; i < NCisAjs; i++) + fprintf(FileCisAjs, "%d %d %d %d % .18e % .18e \n", CisAjsIdx[i][0], CisAjsIdx[i][1], CisAjsIdx[i][2], + CisAjsIdx[i][3], creal(PhysCisAjs[i]), cimag(PhysCisAjs[i])); fprintf(FileCisAjs, "\n"); } - /* zvo_cisajscktalt.dat */ + /* zvo_cisajscktaltex.dat */ if (NCisAjsCktAlt > 0) { for (i = 0; i < NCisAjsCktAlt; i++) fprintf(FileCisAjsCktAlt, "% .18e % .18e ", creal(PhysCisAjsCktAlt[i]), cimag(PhysCisAjsCktAlt[i])); fprintf(FileCisAjsCktAlt, "\n"); } - /* zvo_cisajscktaltdc.dat */ + /* zvo_cisajscktalt.dat */ if (NCisAjsCktAltDC > 0) { for (i = 0; i < NCisAjsCktAltDC; i++) { fprintf(FileCisAjsCktAltDC, "%d %d %d %d %d %d %d %d % .18e % .18e\n", @@ -702,16 +689,16 @@ void outputData() { if (NLanczosMode > 0) { if (AllComplexFlag == 0) { //real PhysCalLanczos_real( - QQQQ_real, QCisAjsQ_real, QCisAjsCktAltQ_real, - NLSHam, Nsite, NCisAjs, NCisAjsLz, iOneBodyGIdx, CisAjsLzIdx, NCisAjsCktAltDC, CisAjsCktAltDCIdx, NLanczosMode, - FileLS, FileLSQQQQ, FileLSQCisAjsQ, FileLSQCisAjsCktAltQ, - FileLSCisAjs, FileLSCisAjsCktAlt); + QQQQ_real, QCisAjsQ_real, QCisAjsCktAltQ_real, QCisAjsCktAltQDC_real, + NLSHam, Nsite, NCisAjs, NCisAjs, iOneBodyGIdx, CisAjsIdx, NCisAjsCktAlt, NCisAjsCktAltDC, CisAjsCktAltDCIdx, + NLanczosMode, FileLS, FileLSQQQQ, FileLSQCisAjsQ, FileLSQCisAjsCktAltQ, + FileLSCisAjs, FileLSCisAjsCktAlt, FileLSCisAjsCktAltDC); }else { //complex PhysCalLanczos_fcmp( - QQQQ, QCisAjsQ, QCisAjsCktAltQ, - NLSHam, Nsite, NCisAjs, NCisAjsLz, iOneBodyGIdx, CisAjsLzIdx, NCisAjsCktAltDC, CisAjsCktAltDCIdx, NLanczosMode, - FileLS, FileLSQQQQ, FileLSQCisAjsQ, FileLSQCisAjsCktAltQ, - FileLSCisAjs, FileLSCisAjsCktAlt); + QQQQ, QCisAjsQ, QCisAjsCktAltQ, QCisAjsCktAltQDC, + NLSHam, Nsite, NCisAjs, NCisAjs, iOneBodyGIdx, CisAjsIdx, NCisAjsCktAlt, NCisAjsCktAltDC, CisAjsCktAltDCIdx, + NLanczosMode, FileLS, FileLSQQQQ, FileLSQCisAjsQ, FileLSQCisAjsCktAltQ, + FileLSCisAjs, FileLSCisAjsCktAlt, FileLSCisAjsCktAltDC); } } }