Skip to content

Commit 4caec30

Browse files
三澤貴宏三澤貴宏
authored andcommitted
Revert "Revert "Merge pull request #76 from k-ido/dev-cg""
This reverts commit b833dbf.
1 parent b833dbf commit 4caec30

19 files changed

Lines changed: 562 additions & 23 deletions

File tree

config/oneapi.cmake

Lines changed: 26 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,26 @@
1+
# for Intel Compiler
2+
3+
set(CMAKE_C_COMPILER "icx" CACHE STRING "" FORCE)
4+
set(CMAKE_CXX_COMPILER "icpx" CACHE STRING "" FORCE)
5+
set(CMAKE_C_FLAGS_DEBUG "-O0 -g -Wall -Wformat -Werror=format-security")
6+
set(CMAKE_C_FLAGS_RELEASE "-Wno-unknown-pragmas -O2 -DNDEBUG -xHost" CACHE STRING "" FORCE)
7+
set(CMAKE_CXX_FLAGS_RELEASE " -Wno-unknown-pragmas -O2 -DNDEBUG -xHost" CACHE STRING "" FORCE)
8+
9+
set(CMAKE_Fortran_COMPILER "ifx" CACHE STRING "" FORCE)
10+
set(CMAKE_Fortran_FLAGS_RELEASE "-O3 -DNDEBUG -xHost" CACHE STRING "" FORCE)
11+
12+
# for Intel MKL
13+
set(BLA_VENDOR "Intel10_64lp" CACHE STRING "" FORCE)
14+
15+
# for BLIS & block pfaffian
16+
set(BLIS_ARTIFACT_CONFIG "intel64")
17+
18+
if(USE_SCALAPACK)
19+
if(SCALAPACK_LIBRARIES MATCHES "")
20+
set(SCALAPACK_LIBRARIES "-qmkl=cluster")
21+
#set(SCALAPACK_LIBRARIES "-lmkl_scalapack_lp64 -lmkl_blacs_intelmpi_lp64 -lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -liomp5 -lpthread -lm")
22+
endif(SCALAPACK_LIBRARIES MATCHES "")
23+
24+
message(STATUS "SCALAPACK_LIBRARIES is ${SCALAPACK_LIBRARIES}")
25+
endif(USE_SCALAPACK)
26+

doc/en/source/expert.rst

Lines changed: 12 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -634,6 +634,18 @@ Keywords and parameters
634634
:math:`O(N_\text{p}) + O(N_\text{p}N_\text{MCS})` when
635635
:math:`N_\text{p} > N_\text{MCS}`.
636636

637+
- ``useDiagScale``
638+
639+
**Type :** int-type (0 or 1, default value: 0)
640+
641+
**Description :** The option of using the point Jacobi method (scaling by diagonal elements of :math:`S` matrix) when solving the linear equation :math:`Sx=g` in the SR method by CG method (0: off, 1: on, ``NSRCG`` must be 1).
642+
643+
- ``NneuronGeneral``
644+
645+
**Type :** int-type (default value: 0)
646+
647+
**Description :** The number of neurons :math:`N_{\rm General RBM}` in the hidden layer of RBM.
648+
637649
LocSpin file (locspn.def)
638650
~~~~~~~~~~~~~~~~~~~~~~~~~
639651

doc/ja/source/expert.rst

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -615,6 +615,12 @@ ModParaファイル (modpara.def)
615615
を陽に構築せずに解くことでメモリを削減する [4]_ オプション[NeuscammanUmrigarChan_ ](1で機能On,
616616
``NStore`` は1に固定されます)。
617617

618+
- ``useDiagScale``
619+
620+
**形式 :** int型 (0もしくは1、デフォルト値=0)
621+
622+
**説明 :** SR法での連立一次方程式 :math:`Sx=g` をCG法により解く際に、 Point Jacobi法 (:math:`S` 行列の対角スケーリング)による前処理付きCG法を使用するオプション(1で機能ON, ``NSRCG=1`` である必要がある)。
623+
618624
- ``NneuronGeneral``
619625

620626
**形式 :** int型 (デフォルト値=0)

src/mVMC/include/global.h

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -359,4 +359,7 @@ int Counter_max = 6;
359359
/* 0: hopping, 1: hopping accept, 2: exchange try, 3: exchange accept */
360360
/* 4: local spin flip try, 5 local spin flip accept*/
361361

362+
int useDiagScale=0;
363+
int RescaleSmat=0;
364+
362365
#endif /* _INCLUDE_GLOBAL */

src/mVMC/readdef.c

Lines changed: 43 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -200,6 +200,7 @@ int ReadGreen(char *xNameListFile, int Nca, int **caIdx, int Ncacadc, int **caca
200200
int i, info = 0;
201201

202202
cFileNameListFile = malloc(sizeof(char) * D_CharTmpReadDef * KWIdxInt_end);
203+
// free at vmcmain.c
203204
fprintf(stdout, " Read File %s .\n", xNameListFile);
204205
if (GetFileName(xNameListFile, cFileNameListFile) != 0) {
205206
fprintf(stderr, " error: Definition files(*.def) are incomplete.\n");
@@ -664,6 +665,8 @@ int ReadDefFileNInt(char *xNameListFile, MPI_Comm comm) {
664665
MPI_Bcast(&FlagRBM, 1, MPI_INT, 0, comm);
665666
MPI_Bcast(&NStoreO, 1, MPI_INT, 0, comm); // for NStoreO
666667
MPI_Bcast(&NSRCG, 1, MPI_INT, 0, comm); // for NCG
668+
MPI_Bcast(&RescaleSmat, 1, MPI_INT, 0, comm); // for Rescale S matrix
669+
MPI_Bcast(&useDiagScale, 1, MPI_INT, 0, comm); // for Jacobi preconditioned CG
667670
MPI_Bcast(&AllComplexFlag, 1, MPI_INT, 0, comm); // for Real
668671
MPI_Bcast(&iFlgOrbitalGeneral, 1, MPI_INT, 0, comm); // for fsz
669672
MPI_Bcast(bufDouble, nBufDouble, MPI_DOUBLE, 0, comm);
@@ -746,6 +749,35 @@ int ReadDefFileNInt(char *xNameListFile, MPI_Comm comm) {
746749
APFlag = 0;
747750
}
748751

752+
if (NSRCG == 2){
753+
useDiagScale = 1;
754+
NSRCG = 1;
755+
// if (rank == 0) printf("remark: use preconditioned CG (Diag Scale)\n");
756+
} //else {
757+
// useDiagScale = 0;
758+
//}
759+
if (useDiagScale){
760+
if(NSRCG == 1){
761+
if (rank == 0) printf("remark: use preconditioned CG (Diag Scale)\n");
762+
}else{
763+
if (rank == 0) printf("remark: not use preconditioned CG (Diag Scale) because NSRCG=%d != 1. Use direct method instead.\n",NSRCG);
764+
}
765+
}
766+
767+
if (RescaleSmat){
768+
if (!(NSRCG == 1 || NStoreO != 0)) {
769+
if (rank == 0) {
770+
fprintf(stderr,
771+
"error: invalid RescaleSmat configuration. "
772+
"RescaleSmat=1 requires NSRCG==1 or NStoreO!=0, "
773+
"but got NSRCG=%d, NStoreO=%d.\n",
774+
NSRCG, NStoreO);
775+
}
776+
MPI_Abort(comm, EXIT_FAILURE);
777+
}
778+
if (rank == 0) printf("remark: rescale S matrix \n");
779+
}
780+
749781
if (DSROptStepDt < 0) {
750782
SRFlag = 1; /* diagonalization */
751783
if (rank == 0) fprintf(stderr, "remark: Diagonalization Mode\n");
@@ -1820,6 +1852,9 @@ void SetDefaultValuesModPara(int *bufInt, double *bufDouble) {
18201852
bufDouble[IdxSROptCGTol] = 1.0e-10;
18211853
NStoreO = 1;
18221854
NSRCG = 0;
1855+
RescaleSmat = 0;
1856+
useDiagScale = 0;
1857+
18231858
}
18241859

18251860
int GetInfoFromModPara(int *bufInt, double *bufDouble) {
@@ -1933,6 +1968,10 @@ int GetInfoFromModPara(int *bufInt, double *bufDouble) {
19331968
NStoreO = (int) dtmp;
19341969
} else if (CheckWords(ctmp, "NSRCG") == 0) {
19351970
NSRCG = (int) dtmp;
1971+
} else if (CheckWords(ctmp, "RescaleSmat") == 0) {
1972+
RescaleSmat = (int) dtmp;
1973+
} else if (CheckWords(ctmp, "useDiagScale") == 0) {
1974+
useDiagScale = (int) dtmp;
19361975
//RBM
19371976
} else if (CheckWords(ctmp, "Nneuron") == 0) {
19381977
bufInt[IdxNneuron] = (int) dtmp;
@@ -1959,9 +1998,10 @@ int GetInfoFromModPara(int *bufInt, double *bufDouble) {
19591998
default:
19601999
break;
19612000
}
2001+
fclose(fp);
19622002
}
19632003
}
1964-
fclose(fp);
2004+
//fclose(fp);
19652005
fprintf(stdout, "End: Read ModPara File .\n");
19662006
return iret;
19672007
}
@@ -2105,6 +2145,8 @@ int GetInfoOpt(FILE *fp, int *ArrayOpt, int iComplxFlag, int *iTotalOptCount, in
21052145
fscanf(fp, "%d\n", &(ArrayOpt[2 * fidx])); // TBC real
21062146
if(iComplxFlag>0){
21072147
ArrayOpt[2 * fidx + 1] = ArrayOpt[2 * fidx]; // TBC imaginary
2148+
}else{
2149+
ArrayOpt[2 * fidx + 1] = 0;
21082150
}
21092151
fidx++;
21102152
(iLocalOptCount)++;

src/mVMC/setmemory.c

Lines changed: 38 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -271,11 +271,15 @@ void FreeMemoryDef() {
271271
free(CisAjsIdx);
272272
free(QPTransSgn);
273273
free(QPTrans);
274+
free(QPTransInv);
274275
free(OrbitalIdx);
276+
free(OrbitalSgn);
275277
free(DoublonHolon4siteIdx);
276278
free(DoublonHolon2siteIdx);
277279
free(JastrowIdx);
278280
free(ParaTransfer);
281+
free(ParaCoulombIntra);
282+
free(ParaQPTrans);
279283
free(ExchangeCoupling);
280284
free(PairHopping);
281285
free(HundCoupling);
@@ -359,7 +363,11 @@ void SetMemory() {
359363
InvM = (double complex*)malloc( sizeof(double complex)*(NQPFull*(Nsize*Nsize+1)) );
360364
PfM = InvM + NQPFull*Nsize*Nsize;
361365
// for real TBC
362-
SlaterElm_real = (double*)malloc(sizeof(double)*(NQPFull*(2*Nsite)*(2*Nsite)) );
366+
if (AllComplexFlag == 0){
367+
SlaterElm_real = (double*)malloc(sizeof(double)*(NQPFull*(2*Nsite)*(2*Nsite)) );
368+
InvM_real = (double*)malloc(sizeof(double)*(NQPFull*(Nsize*Nsize+1)) );
369+
PfM_real = InvM_real + NQPFull*Nsize*Nsize;
370+
}
363371

364372
if (FlagRBM) {
365373
SlaterElmBF_real = (double*)malloc( sizeof(double)*(NQPFull*(2*Nsite)*(2*Nsite)) );
@@ -376,8 +384,6 @@ void SetMemory() {
376384
BFSubIdx[i] = (int*)malloc(sizeof(int)*NrangeIdx);
377385
}
378386
}
379-
InvM_real = (double*)malloc(sizeof(double)*(NQPFull*(Nsize*Nsize+1)) );
380-
PfM_real = InvM_real + NQPFull*Nsize*Nsize;
381387

382388
/***** Quantum Projection *****/
383389
QPFullWeight = (double complex*)malloc(sizeof(double complex)*(NQPFull+NQPFix+5*NSPGaussLeg));
@@ -402,15 +408,17 @@ void SetMemory() {
402408
SROptO = SROptHO + 2*SROptSize; //TBC
403409
}
404410
//for real
405-
if(NSRCG==0){
406-
SROptOO_real = (double*)malloc( sizeof(double )*SROptSize*(SROptSize+2)) ; //TBC
407-
SROptHO_real = SROptOO_real + (SROptSize)*(SROptSize); //TBC
408-
SROptO_real = SROptHO_real + (SROptSize); //TBC
409-
}else{
410-
// OO contains only <O_i> and <O_i O_i> in SR-CG
411-
SROptOO_real = (double*)malloc( sizeof(double )*SROptSize*4) ; //TBC
412-
SROptHO_real = SROptOO_real + SROptSize*2; //TBC
413-
SROptO_real = SROptHO_real + SROptSize; //TBC
411+
if (AllComplexFlag == 0){
412+
if(NSRCG==0){
413+
SROptOO_real = (double*)malloc( sizeof(double )*SROptSize*(SROptSize+2)) ; //TBC
414+
SROptHO_real = SROptOO_real + (SROptSize)*(SROptSize); //TBC
415+
SROptO_real = SROptHO_real + (SROptSize); //TBC
416+
}else{
417+
// OO contains only <O_i> and <O_i O_i> in SR-CG
418+
SROptOO_real = (double*)malloc( sizeof(double )*SROptSize*4) ; //TBC
419+
SROptHO_real = SROptOO_real + SROptSize*2; //TBC
420+
SROptO_real = SROptHO_real + SROptSize; //TBC
421+
}
414422
}
415423

416424
if(NSRCG==1 || NStoreO!=0){
@@ -481,18 +489,36 @@ void FreeMemory() {
481489
if(NVMCCalMode==0){
482490
free(SROptData);
483491
free(SROptOO);
492+
//for real
493+
if (AllComplexFlag == 0){
494+
free(SROptOO_real);
495+
}
496+
if(NSRCG==1 || NStoreO!=0){
497+
if(AllComplexFlag==0){ //real & sz=0
498+
free(SROptO_Store_real);
499+
}else{
500+
free(SROptO_Store);
501+
}
502+
}
484503
}
485504

486505
free(QPFullWeight);
487506

488507
free(InvM);
489508
free(SlaterElm);
490509

510+
if (AllComplexFlag == 0){
511+
free(InvM_real);
512+
free(SlaterElm_real);
513+
}
514+
491515
free(BurnEleIdx);
492516
free(TmpEleIdx);
493517
free(logSqPfFullSlater);
494518
free(EleProjCnt);
495519
free(EleIdx);
520+
free(EleNum);
521+
free(EleSpn);
496522
free(EleCfg);
497523

498524
free(Para);

src/mVMC/stcopt_cg.c

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -32,6 +32,19 @@ inline double xdot(const int n, const double * const p, const double * const q)
3232
}
3333
return z;
3434
}
35+
36+
inline double get_absmax(int idx_start, int idx_end, const double * const p) {
37+
int pi;
38+
if (idx_end <= idx_start) return 0.0;
39+
double absval_max=fabs(p[idx_start]);
40+
41+
#pragma loop noalias
42+
for ( pi=idx_start+1; pi < idx_end; pi ++ ) {
43+
if ( absval_max < fabs(p[pi]) ) absval_max = fabs(p[pi]) ;
44+
}
45+
return absval_max;
46+
}
47+
3548
extern inline double xdot(const int n, const double * const p, const double * const q);
3649

3750
#define MVMC_SRCG_REAL
@@ -43,8 +56,10 @@ int StochasticOptCG(MPI_Comm comm)
4356
{
4457
int ret=0;
4558
if(AllComplexFlag==0){
59+
if(RescaleSmat) Rescale4SRCG_real(comm);
4660
ret = StochasticOptCG_real(comm);
4761
}else{
62+
if(RescaleSmat) Rescale4SRCG_fcmp(comm);
4863
ret = StochasticOptCG_fcmp(comm);
4964
}
5065
return ret;

0 commit comments

Comments
 (0)