From 026e8f7c0696f8048ffbc3abcf4663925311462c Mon Sep 17 00:00:00 2001 From: Ye Luo Date: Tue, 9 Jan 2018 16:43:06 -0600 Subject: [PATCH] Replace 1 nel det with 2 nel/2 dets. --- src/QMCWaveFunctions/Determinant.h | 13 +++++--- src/QMCWaveFunctions/DeterminantRef.h | 13 +++++--- src/QMCWaveFunctions/WaveFunction.cpp | 42 +++++++++++++++++------- src/QMCWaveFunctions/WaveFunction.h | 8 +++-- src/QMCWaveFunctions/WaveFunctionRef.cpp | 42 +++++++++++++++++------- 5 files changed, 84 insertions(+), 34 deletions(-) diff --git a/src/QMCWaveFunctions/Determinant.h b/src/QMCWaveFunctions/Determinant.h index 9f8b547a4..862057fa5 100644 --- a/src/QMCWaveFunctions/Determinant.h +++ b/src/QMCWaveFunctions/Determinant.h @@ -178,7 +178,8 @@ void checkDiff(const MT1 &a, const MT2 &b, const std::string &tag) struct DiracDeterminant : public WaveFunctionComponentBase { - DiracDeterminant(int nels, RandomGenerator RNG) + DiracDeterminant(int nels, RandomGenerator RNG, int First=0) + : FirstIndex(First) { psiMinv.resize(nels, nels); psiV.resize(nels); @@ -222,7 +223,7 @@ struct DiracDeterminant : public WaveFunctionComponentBase GradType evalGrad(ParticleSet &P, int iat) {} - ValueType ratioGrad(ParticleSet &P, int iat, GradType &grad) {} + ValueType ratioGrad(ParticleSet &P, int iat, GradType &grad) { return ratio(P, iat); } void evaluateGL(ParticleSet &P, ParticleSet::ParticleGradient_t &G, ParticleSet::ParticleLaplacian_t &L, bool fromscratch = false) {} @@ -246,7 +247,7 @@ struct DiracDeterminant : public WaveFunctionComponentBase constexpr double czero(0); for (int j = 0; j < nels; ++j) psiV[j] = myRandom() - shift; - curRatio = inner_product_n(psiV.data(), psiMinv[iel], nels, czero); + curRatio = inner_product_n(psiV.data(), psiMinv[iel-FirstIndex], nels, czero); return curRatio; } @@ -254,8 +255,8 @@ struct DiracDeterminant : public WaveFunctionComponentBase inline void acceptMove(ParticleSet &P, int iel) { const int nels = psiV.size(); - updateRow(psiMinv.data(), psiV.data(), nels, nels, iel, curRatio); - std::copy_n(psiV.data(), nels, psiMsave[iel]); + updateRow(psiMinv.data(), psiV.data(), nels, nels, iel-FirstIndex, curRatio); + std::copy_n(psiV.data(), nels, psiMsave[iel-FirstIndex]); } /** accessor functions for checking */ @@ -269,6 +270,8 @@ struct DiracDeterminant : public WaveFunctionComponentBase double curRatio; /// workspace size int LWork; + /// initial particle index + const int FirstIndex; /// inverse matrix to be update Matrix psiMinv; /// a SPO set for the row update diff --git a/src/QMCWaveFunctions/DeterminantRef.h b/src/QMCWaveFunctions/DeterminantRef.h index e1531eb54..a766e5485 100644 --- a/src/QMCWaveFunctions/DeterminantRef.h +++ b/src/QMCWaveFunctions/DeterminantRef.h @@ -182,7 +182,8 @@ struct DiracDeterminantRef : public qmcplusplus::WaveFunctionComponentBase { using ParticleSet = qmcplusplus::ParticleSet; - DiracDeterminantRef(int nels, qmcplusplus::RandomGenerator RNG) + DiracDeterminantRef(int nels, qmcplusplus::RandomGenerator RNG, int First=0) + : FirstIndex(First) { psiMinv.resize(nels, nels); psiV.resize(nels); @@ -225,7 +226,7 @@ struct DiracDeterminantRef : public qmcplusplus::WaveFunctionComponentBase // FIXME do we want remainder of evaluateLog? } GradType evalGrad(ParticleSet &P, int iat) {} - ValueType ratioGrad(ParticleSet &P, int iat, GradType &grad) {} + ValueType ratioGrad(ParticleSet &P, int iat, GradType &grad) { return ratio(P, iat); } void evaluateGL(ParticleSet &P, ParticleSet::ParticleGradient_t &G, ParticleSet::ParticleLaplacian_t &L, bool fromscratch = false) {} @@ -248,7 +249,7 @@ struct DiracDeterminantRef : public qmcplusplus::WaveFunctionComponentBase constexpr double czero(0); for (int j = 0; j < nels; ++j) psiV[j] = myRandom() - shift; - curRatio = inner_product_n(psiV.data(), psiMinv[iel], nels, czero); + curRatio = inner_product_n(psiV.data(), psiMinv[iel-FirstIndex], nels, czero); return curRatio; } @@ -256,8 +257,8 @@ struct DiracDeterminantRef : public qmcplusplus::WaveFunctionComponentBase inline void acceptMove(ParticleSet &P, int iel) { const int nels = psiV.size(); - updateRow(psiMinv.data(), psiV.data(), nels, nels, iel, curRatio); - std::copy_n(psiV.data(), nels, psiMsave[iel]); + updateRow(psiMinv.data(), psiV.data(), nels, nels, iel-FirstIndex, curRatio); + std::copy_n(psiV.data(), nels, psiMsave[iel-FirstIndex]); } /** accessor functions for checking */ @@ -271,6 +272,8 @@ struct DiracDeterminantRef : public qmcplusplus::WaveFunctionComponentBase double curRatio; /// workspace size int LWork; + /// initial particle index + const int FirstIndex; /// inverse matrix to be update qmcplusplus::Matrix psiMinv; /// a SPO set for the row update diff --git a/src/QMCWaveFunctions/WaveFunction.cpp b/src/QMCWaveFunctions/WaveFunction.cpp index d4f602e27..cadf628c7 100644 --- a/src/QMCWaveFunctions/WaveFunction.cpp +++ b/src/QMCWaveFunctions/WaveFunction.cpp @@ -33,7 +33,9 @@ WaveFunction::WaveFunction(ParticleSet &ions, ParticleSet &els, d_ie = DistanceTable::add(ions, els, DT_SOA); // determinant component - Det = new DetType(els.getTotalNum(), RNG); + nelup = els.getTotalNum()/2; + Det_up = new DetType(nelup, RNG, 0); + Det_dn = new DetType(els.getTotalNum()-nelup, RNG, nelup); // J1 component J1 = new J1OrbType(ions, els); @@ -48,7 +50,14 @@ WaveFunction::WaveFunction(ParticleSet &ions, ParticleSet &els, buildJeeI(*J3, els.Lattice.WignerSeitzRadius); } -WaveFunction::~WaveFunction() { delete J2; } +WaveFunction::~WaveFunction() +{ + delete Det_up; + delete Det_dn; + delete J1; + delete J2; + delete J3; +} void WaveFunction::evaluateLog(ParticleSet &P) { @@ -57,7 +66,8 @@ void WaveFunction::evaluateLog(ParticleSet &P) { P.G = czero; P.L = czero; - LogValue = Det->evaluateLog(P, P.G, P.L); + LogValue = Det_up->evaluateLog(P, P.G, P.L); + LogValue += Det_dn->evaluateLog(P, P.G, P.L); LogValue += J1->evaluateLog(P, P.G, P.L); LogValue += J2->evaluateLog(P, P.G, P.L); LogValue += J3->evaluateLog(P, P.G, P.L); @@ -67,26 +77,35 @@ void WaveFunction::evaluateLog(ParticleSet &P) WaveFunctionBase::posT WaveFunction::evalGrad(ParticleSet &P, int iat) { - return Det->evalGrad(P, iat) + J1->evalGrad(P, iat) + J2->evalGrad(P, iat) + - J3->evalGrad(P, iat); + return ( iatevalGrad(P, iat) : Det_dn->evalGrad(P, iat) ) + + J1->evalGrad(P, iat) + + J2->evalGrad(P, iat) + + J3->evalGrad(P, iat); } WaveFunctionBase::valT WaveFunction::ratioGrad(ParticleSet &P, int iat, posT &grad) { - return Det->ratioGrad(P, iat, grad) + J1->ratioGrad(P, iat, grad) + - J2->ratioGrad(P, iat, grad) + J3->ratioGrad(P, iat, grad); + return ( iatratioGrad(P, iat, grad) : Det_dn->ratioGrad(P, iat, grad) ) + * J1->ratioGrad(P, iat, grad) + * J2->ratioGrad(P, iat, grad) + * J3->ratioGrad(P, iat, grad); } WaveFunctionBase::valT WaveFunction::ratio(ParticleSet &P, int iat) { - return Det->ratio(P, iat) * J1->ratio(P, iat) * J2->ratio(P, iat) * - J3->ratio(P, iat); + return ( iatratio(P, iat) : Det_dn->ratio(P, iat) ) + * J1->ratio(P, iat) + * J2->ratio(P, iat) + * J3->ratio(P, iat); } void WaveFunction::acceptMove(ParticleSet &P, int iat) { - Det->acceptMove(P, iat); + if(iatacceptMove(P, iat); + else + Det_dn->acceptMove(P, iat); J1->acceptMove(P, iat); J2->acceptMove(P, iat); J3->acceptMove(P, iat); @@ -99,7 +118,8 @@ void WaveFunction::evaluateGL(ParticleSet &P) constexpr valT czero(0); P.G = czero; P.L = czero; - Det->evaluateGL(P, P.G, P.L); + Det_up->evaluateGL(P, P.G, P.L); + Det_dn->evaluateGL(P, P.G, P.L); J1->evaluateGL(P, P.G, P.L); J2->evaluateGL(P, P.G, P.L); J3->evaluateGL(P, P.G, P.L); diff --git a/src/QMCWaveFunctions/WaveFunction.h b/src/QMCWaveFunctions/WaveFunction.h index f8209c902..1c9a25fba 100644 --- a/src/QMCWaveFunctions/WaveFunction.h +++ b/src/QMCWaveFunctions/WaveFunction.h @@ -73,7 +73,9 @@ struct WaveFunction : public WaveFunctionBase J1OrbType *J1; J2OrbType *J2; J3OrbType *J3; - DetType *Det; + int nelup; + DetType *Det_up; + DetType *Det_dn; WaveFunction(ParticleSet &ions, ParticleSet &els, RandomGenerator RNG); @@ -106,7 +108,9 @@ struct WaveFunctionRef : public qmcplusplus::WaveFunctionBase J1OrbType *J1; J2OrbType *J2; J3OrbType *J3; - DetType *Det; + int nelup; + DetType *Det_up; + DetType *Det_dn; PooledData Buffer; WaveFunctionRef(ParticleSet &ions, ParticleSet &els, diff --git a/src/QMCWaveFunctions/WaveFunctionRef.cpp b/src/QMCWaveFunctions/WaveFunctionRef.cpp index a0e9e9f41..4e2a878f3 100644 --- a/src/QMCWaveFunctions/WaveFunctionRef.cpp +++ b/src/QMCWaveFunctions/WaveFunctionRef.cpp @@ -39,7 +39,9 @@ WaveFunctionRef::WaveFunctionRef(ParticleSet &ions, ParticleSet &els, d_ie = DistanceTable::add(ions, els, DT_SOA); // determinant component - Det = new DetType(els.getTotalNum(), RNG); + nelup = els.getTotalNum()/2; + Det_up = new DetType(nelup, RNG, 0); + Det_dn = new DetType(els.getTotalNum()-nelup, RNG, nelup); // J1 component J1 = new J1OrbType(ions, els); @@ -54,7 +56,14 @@ WaveFunctionRef::WaveFunctionRef(ParticleSet &ions, ParticleSet &els, buildJeeI(*J3, els.Lattice.WignerSeitzRadius); } -WaveFunctionRef::~WaveFunctionRef() { delete J2; } +WaveFunctionRef::~WaveFunctionRef() +{ + delete Det_up; + delete Det_dn; + delete J1; + delete J2; + delete J3; +} void WaveFunctionRef::evaluateLog(ParticleSet &P) { @@ -63,7 +72,8 @@ void WaveFunctionRef::evaluateLog(ParticleSet &P) { P.G = czero; P.L = czero; - LogValue = Det->evaluateLog(P, P.G, P.L); + LogValue = Det_up->evaluateLog(P, P.G, P.L); + LogValue += Det_dn->evaluateLog(P, P.G, P.L); LogValue += J1->evaluateLog(P, P.G, P.L); LogValue += J2->evaluateLog(P, P.G, P.L); LogValue += J3->evaluateLog(P, P.G, P.L); @@ -73,26 +83,35 @@ void WaveFunctionRef::evaluateLog(ParticleSet &P) WaveFunctionBase::posT WaveFunctionRef::evalGrad(ParticleSet &P, int iat) { - return Det->evalGrad(P, iat) + J1->evalGrad(P, iat) + J2->evalGrad(P, iat) + - J3->evalGrad(P, iat); + return ( iatevalGrad(P, iat) : Det_dn->evalGrad(P, iat) ) + + J1->evalGrad(P, iat) + + J2->evalGrad(P, iat) + + J3->evalGrad(P, iat); } WaveFunctionBase::valT WaveFunctionRef::ratioGrad(ParticleSet &P, int iat, posT &grad) { - return Det->ratioGrad(P, iat, grad) + J1->ratioGrad(P, iat, grad) + - J2->ratioGrad(P, iat, grad) + J3->ratioGrad(P, iat, grad); + return ( iatratioGrad(P, iat, grad) : Det_dn->ratioGrad(P, iat, grad) ) + * J1->ratioGrad(P, iat, grad) + * J2->ratioGrad(P, iat, grad) + * J3->ratioGrad(P, iat, grad); } WaveFunctionBase::valT WaveFunctionRef::ratio(ParticleSet &P, int iat) { - return Det->ratio(P, iat) * J1->ratio(P, iat) * J2->ratio(P, iat) * - J3->ratio(P, iat); + return ( iatratio(P, iat) : Det_dn->ratio(P, iat) ) + * J1->ratio(P, iat) + * J2->ratio(P, iat) + * J3->ratio(P, iat); } void WaveFunctionRef::acceptMove(ParticleSet &P, int iat) { - Det->acceptMove(P, iat); + if(iatacceptMove(P, iat); + else + Det_dn->acceptMove(P, iat); J1->acceptMove(P, iat); J2->acceptMove(P, iat); J3->acceptMove(P, iat); @@ -105,7 +124,8 @@ void WaveFunctionRef::evaluateGL(ParticleSet &P) constexpr valT czero(0); P.G = czero; P.L = czero; - Det->evaluateGL(P, P.G, P.L); + Det_up->evaluateGL(P, P.G, P.L); + Det_dn->evaluateGL(P, P.G, P.L); J1->evaluateGL(P, P.G, P.L); J2->evaluateGL(P, P.G, P.L); J3->evaluateGL(P, P.G, P.L);