diff --git a/src/QMCDrivers/WFOpt/QMCCostFunction.cpp b/src/QMCDrivers/WFOpt/QMCCostFunction.cpp index cc1f4ba3c9..cfbfec7375 100644 --- a/src/QMCDrivers/WFOpt/QMCCostFunction.cpp +++ b/src/QMCDrivers/WFOpt/QMCCostFunction.cpp @@ -110,7 +110,7 @@ void QMCCostFunction::GradCost(std::vector& PGradient, ltz = false; Return_rt delE = std::pow(std::abs(eloc_new - EtargetEff), PowerE); Return_rt ddelE = PowerE * std::pow(std::abs(eloc_new - EtargetEff), PowerE - 1); - const Return_t* Dsaved = (*DerivRecords[ip])[iw]; + const Return_t* Dsaved = (*DerivRecords[ip])[iw]; const Return_rt* HDsaved = (*HDerivRecords[ip])[iw]; for (int pm = 0; pm < NumOptimizables; pm++) { @@ -137,7 +137,7 @@ void QMCCostFunction::GradCost(std::vector& PGradient, Return_rt eloc_new = saved[ENERGY_NEW]; Return_rt delta_l = (eloc_new - curAvg_w); Return_rt sigma_l = delta_l * delta_l; - const Return_t* Dsaved = (*DerivRecords[ip])[iw]; + const Return_t* Dsaved = (*DerivRecords[ip])[iw]; const Return_rt* HDsaved = (*HDerivRecords[ip])[iw]; for (int pm = 0; pm < NumOptimizables; pm++) { @@ -643,17 +643,7 @@ QMCCostFunction::Return_rt QMCCostFunction::fillOverlapHamiltonianMatrices(Matri { ScopedTimer tmp_timer(fill_timer_); - RealType b1, b2; - if (GEVType == "H2") - { - b1 = w_beta; - b2 = 0; - } - else - { - b2 = w_beta; - b1 = 0; - } + RealType b2(w_beta); Right = 0.0; Left = 0.0; @@ -661,10 +651,7 @@ QMCCostFunction::Return_rt QMCCostFunction::fillOverlapHamiltonianMatrices(Matri // resetPsi(); curAvg_w = SumValue[SUM_E_WGT] / SumValue[SUM_WGT]; Return_rt curAvg2_w = SumValue[SUM_ESQ_WGT] / SumValue[SUM_WGT]; - // RealType H2_avg = 1.0/curAvg2_w; - RealType H2_avg = 1.0 / (curAvg_w * curAvg_w); - // RealType H2_avg = 1.0/std::sqrt(curAvg_w*curAvg_w*curAvg2_w); - RealType V_avg = curAvg2_w - curAvg_w * curAvg_w; + RealType V_avg = curAvg2_w - curAvg_w * curAvg_w; std::vector D_avg(getNumParams(), 0.0); Return_rt wgtinv = 1.0 / SumValue[SUM_WGT]; for (int ip = 0; ip < NumThreads; ip++) @@ -674,7 +661,7 @@ QMCCostFunction::Return_rt QMCCostFunction::fillOverlapHamiltonianMatrices(Matri { const Return_rt* restrict saved = (*RecordsOnNode[ip])[iw]; Return_rt weight = saved[REWEIGHT] * wgtinv; - const Return_t* Dsaved = (*DerivRecords[ip])[iw]; + const Return_t* Dsaved = (*DerivRecords[ip])[iw]; for (int pm = 0; pm < getNumParams(); pm++) { D_avg[pm] += Dsaved[pm] * weight; @@ -692,19 +679,16 @@ QMCCostFunction::Return_rt QMCCostFunction::fillOverlapHamiltonianMatrices(Matri const Return_rt* restrict saved = (*RecordsOnNode[ip])[iw]; Return_rt weight = saved[REWEIGHT] * wgtinv; Return_rt eloc_new = saved[ENERGY_NEW]; - const Return_t* Dsaved = (*DerivRecords[ip])[iw]; + const Return_t* Dsaved = (*DerivRecords[ip])[iw]; const Return_rt* HDsaved = (*HDerivRecords[ip])[iw]; #pragma omp parallel for for (int pm = 0; pm < getNumParams(); pm++) { - Return_t wfe = (HDsaved[pm] + (Dsaved[pm] - D_avg[pm]) * eloc_new) * weight; - Return_t wfd = (Dsaved[pm] - D_avg[pm]) * weight; - Return_t vterm = - HDsaved[pm] * (eloc_new - curAvg_w) + (Dsaved[pm] - D_avg[pm]) * eloc_new * (eloc_new - RealType(2.0) * curAvg_w); + Return_t wfe = (HDsaved[pm] + (Dsaved[pm] - D_avg[pm]) * eloc_new) * weight; + Return_t wfd = (Dsaved[pm] - D_avg[pm]) * weight; + Return_t vterm = HDsaved[pm] * (eloc_new - curAvg_w) + + (Dsaved[pm] - D_avg[pm]) * eloc_new * (eloc_new - RealType(2.0) * curAvg_w); // Return_t vterm = (HDsaved[pm]+(Dsaved[pm]-D_avg[pm])*eloc_new -curAvg_w)*(eloc_new-curAvg_w); - // H2 - Right(0, pm + 1) += b1 * H2_avg * std::real(vterm) * weight; - Right(pm + 1, 0) += b1 * H2_avg * std::real(vterm) * weight; // Variance Left(0, pm + 1) += b2 * std::real(vterm) * weight; Left(pm + 1, 0) += b2 * std::real(vterm) * weight; @@ -714,18 +698,18 @@ QMCCostFunction::Return_rt QMCCostFunction::fillOverlapHamiltonianMatrices(Matri for (int pm2 = 0; pm2 < getNumParams(); pm2++) { // Hamiltonian - Left(pm + 1, pm2 + 1) += std::real((1 - b2) * std::conj(wfd) * (HDsaved[pm2] + (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); + Left(pm + 1, pm2 + 1) += + std::real((1 - b2) * std::conj(wfd) * (HDsaved[pm2] + (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); // Overlap RealType ovlij = std::real(std::conj(wfd) * (Dsaved[pm2] - D_avg[pm2])); Right(pm + 1, pm2 + 1) += ovlij; // Variance - RealType varij = weight * std::real((HDsaved[pm] - RealType(2.0) * std::conj(Dsaved[pm] - D_avg[pm]) * eloc_new) * - (HDsaved[pm2] - RealType(2.0) * (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); + RealType varij = weight * + std::real((HDsaved[pm] - RealType(2.0) * std::conj(Dsaved[pm] - D_avg[pm]) * eloc_new) * + (HDsaved[pm2] - RealType(2.0) * (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); // RealType varij=weight*(HDsaved[pm] +(Dsaved[pm]-D_avg[pm])*eloc_new-curAvg_w)* // (HDsaved[pm2] + (Dsaved[pm2]-D_avg[pm2])*eloc_new-curAvg_w); Left(pm + 1, pm2 + 1) += b2 * (varij + V_avg * ovlij); - // H2 - Right(pm + 1, pm2 + 1) += b1 * H2_avg * varij; } } } @@ -733,9 +717,7 @@ QMCCostFunction::Return_rt QMCCostFunction::fillOverlapHamiltonianMatrices(Matri myComm->allreduce(Right); myComm->allreduce(Left); Left(0, 0) = (1 - b2) * curAvg_w + b2 * V_avg; - Right(0, 0) = 1.0 + b1 * H2_avg * V_avg; - if (GEVType == "H2") - return H2_avg; + Right(0, 0) = 1.0; return 1.0; } diff --git a/src/QMCDrivers/WFOpt/QMCCostFunctionBase.cpp b/src/QMCDrivers/WFOpt/QMCCostFunctionBase.cpp index fc6e7a6483..e22a10d102 100644 --- a/src/QMCDrivers/WFOpt/QMCCostFunctionBase.cpp +++ b/src/QMCDrivers/WFOpt/QMCCostFunctionBase.cpp @@ -47,7 +47,6 @@ QMCCostFunctionBase::QMCCostFunctionBase(ParticleSet& w, TrialWaveFunction& psi, w_w(0.0), MaxWeight(1e6), w_beta(0.0), - GEVType("mixed"), vmc_or_dmc(2.0), needGrads(true), targetExcitedStr("no"), @@ -318,6 +317,7 @@ bool QMCCostFunctionBase::put(xmlNodePtr q) std::string includeNonlocalH; std::string writeXmlPerStep("no"); std::string computeNLPPderiv; + std::string GEVType; astring variational_subset_str; ParameterSet m_param; m_param.add(writeXmlPerStep, "dumpXML"); @@ -326,7 +326,7 @@ bool QMCCostFunctionBase::put(xmlNodePtr q) m_param.add(includeNonlocalH, "nonlocalpp", {}, TagStatus::DEPRECATED); m_param.add(computeNLPPderiv, "use_nonlocalpp_deriv", {}, TagStatus::DEPRECATED); m_param.add(w_beta, "beta"); - m_param.add(GEVType, "GEVMethod"); + m_param.add(GEVType, "GEVMethod", {}, TagStatus::DEPRECATED); m_param.add(targetExcitedStr, "targetExcited"); m_param.add(omega_shift, "omega"); m_param.add(do_override_output, "output_vp_override", {true}); diff --git a/src/QMCDrivers/WFOpt/QMCCostFunctionBase.h b/src/QMCDrivers/WFOpt/QMCCostFunctionBase.h index 7394803cb2..c8ec41ee9a 100644 --- a/src/QMCDrivers/WFOpt/QMCCostFunctionBase.h +++ b/src/QMCDrivers/WFOpt/QMCCostFunctionBase.h @@ -221,7 +221,6 @@ class QMCCostFunctionBase : public CostFunctionBase, public Return_rt curVar_abs; Return_rt w_beta; - std::string GEVType; Return_rt vmc_or_dmc; bool needGrads; ///whether we are targeting an excited state diff --git a/src/QMCDrivers/WFOpt/QMCCostFunctionBatched.cpp b/src/QMCDrivers/WFOpt/QMCCostFunctionBatched.cpp index 1e260bf26c..bb71af64c4 100644 --- a/src/QMCDrivers/WFOpt/QMCCostFunctionBatched.cpp +++ b/src/QMCDrivers/WFOpt/QMCCostFunctionBatched.cpp @@ -112,12 +112,12 @@ void QMCCostFunctionBatched::GradCost(std::vector& PGradient, ltz = false; Return_rt delE = std::pow(std::abs(eloc_new - EtargetEff), PowerE); Return_rt ddelE = PowerE * std::pow(std::abs(eloc_new - EtargetEff), PowerE - 1); - const Return_t* Dsaved = DerivRecords_[iw]; + const Return_t* Dsaved = DerivRecords_[iw]; const Return_rt* HDsaved = HDerivRecords_[iw]; for (int pm = 0; pm < NumOptimizables; pm++) { //From Toulouse J. Chem. Phys. 126, 084102 (2007), this is H_0j+H_j0, which are independent - //estimates of 1/2 the energy gradient g. So g1+g2 is an estimate of g. + //estimates of 1/2 the energy gradient g. So g1+g2 is an estimate of g. EDtotals_w[pm] += weight * (HDsaved[pm] + 2.0 * std::real(Dsaved[pm]) * delta_l); URV[pm] += 2.0 * (eloc_new * HDsaved[pm] - curAvg * HD_avg[pm]); if (ltz) @@ -139,7 +139,7 @@ void QMCCostFunctionBatched::GradCost(std::vector& PGradient, Return_rt eloc_new = saved[ENERGY_NEW]; Return_rt delta_l = (eloc_new - curAvg_w); Return_rt sigma_l = delta_l * delta_l; - const Return_t* Dsaved = DerivRecords_[iw]; + const Return_t* Dsaved = DerivRecords_[iw]; const Return_rt* HDsaved = HDerivRecords_[iw]; for (int pm = 0; pm < NumOptimizables; pm++) { @@ -193,7 +193,7 @@ void QMCCostFunctionBatched::getConfigurations(const std::string& aroot) } } - /** Compute number of batches and final batch size given the number of samples +/** Compute number of batches and final batch size given the number of samples * and a batch size. * \param[in] sample_size number of samples to process. * \param[in] batch_size process samples in batch_size at a time (typically the number of walkers in a crowd). @@ -354,7 +354,7 @@ void QMCCostFunctionBatched::checkConfigurations(EngineHandle& handle) for (int j = 0; j < nparams; j++) { //dlogpsi is in general complex if psi is complex. - DerivRecords[is][j] = dlogpsi_array[ib][j]; + DerivRecords[is][j] = dlogpsi_array[ib][j]; //but E_L and d E_L/dc are real if c is real. HDerivRecords[is][j] = std::real(dhpsioverpsi_array[ib][j]); } @@ -594,7 +594,7 @@ QMCCostFunctionBatched::EffectiveWeight QMCCostFunctionBatched::correlatedSampli if (optVars.recompute(j)) { //In general, dlogpsi is complex. - DerivRecords[is][j] = dlogpsi_array[ib][j]; + DerivRecords[is][j] = dlogpsi_array[ib][j]; //However, E_L is always real, and so d E_L/dc is real, provided c is real. HDerivRecords[is][j] = std::real(dhpsioverpsi_array[ib][j]); } @@ -700,28 +700,12 @@ QMCCostFunctionBatched::Return_rt QMCCostFunctionBatched::fillOverlapHamiltonian { ScopedTimer tmp_timer(fill_timer_); - RealType b1, b2; - if (GEVType == "H2") - { - b1 = w_beta; - b2 = 0; - } - else - { - b2 = w_beta; - b1 = 0; - } - Right = 0.0; Left = 0.0; - // resetPsi(); curAvg_w = SumValue[SUM_E_WGT] / SumValue[SUM_WGT]; Return_rt curAvg2_w = SumValue[SUM_ESQ_WGT] / SumValue[SUM_WGT]; - // RealType H2_avg = 1.0/curAvg2_w; - RealType H2_avg = 1.0 / (curAvg_w * curAvg_w); - // RealType H2_avg = 1.0/std::sqrt(curAvg_w*curAvg_w*curAvg2_w); - RealType V_avg = curAvg2_w - curAvg_w * curAvg_w; + RealType V_avg = curAvg2_w - curAvg_w * curAvg_w; std::vector D_avg(getNumParams(), 0.0); Return_rt wgtinv = 1.0 / SumValue[SUM_WGT]; @@ -729,7 +713,7 @@ QMCCostFunctionBatched::Return_rt QMCCostFunctionBatched::fillOverlapHamiltonian { const Return_rt* restrict saved = RecordsOnNode_[iw]; Return_rt weight = saved[REWEIGHT] * wgtinv; - const Return_t* Dsaved = DerivRecords_[iw]; + const Return_t* Dsaved = DerivRecords_[iw]; for (int pm = 0; pm < getNumParams(); pm++) { D_avg[pm] += Dsaved[pm] * weight; @@ -743,7 +727,7 @@ QMCCostFunctionBatched::Return_rt QMCCostFunctionBatched::fillOverlapHamiltonian const Return_rt* restrict saved = RecordsOnNode_[iw]; Return_rt weight = saved[REWEIGHT] * wgtinv; Return_rt eloc_new = saved[ENERGY_NEW]; - const Return_t* Dsaved = DerivRecords_[iw]; + const Return_t* Dsaved = DerivRecords_[iw]; const Return_rt* HDsaved = HDerivRecords_[iw]; size_t opt_num_crowds = walkers_per_crowd_.size(); @@ -752,21 +736,18 @@ QMCCostFunctionBatched::Return_rt QMCCostFunctionBatched::fillOverlapHamiltonian auto constructMatrices = [](int crowd_id, std::vector& crowd_ranges, int numParams, const Return_t* Dsaved, - const Return_rt* HDsaved, Return_rt weight, Return_rt eloc_new, RealType H2_avg, - RealType V_avg, std::vector& D_avg, RealType b1, RealType b2, - RealType curAvg_w, Matrix& Left, Matrix& Right) { + const Return_rt* HDsaved, Return_rt weight, Return_rt eloc_new, RealType V_avg, + std::vector& D_avg, RealType b2, RealType curAvg_w, Matrix& Left, + Matrix& Right) { int local_pm_start = crowd_ranges[crowd_id]; int local_pm_end = crowd_ranges[crowd_id + 1]; for (int pm = local_pm_start; pm < local_pm_end; pm++) { - Return_t wfe = (HDsaved[pm] + (Dsaved[pm] - D_avg[pm]) * eloc_new) * weight; - Return_t wfd = (Dsaved[pm] - D_avg[pm]) * weight; - Return_t vterm = - HDsaved[pm] * (eloc_new - curAvg_w) + (Dsaved[pm] - D_avg[pm]) * eloc_new * (eloc_new - RealType(2.0) * curAvg_w); - // H2 - Right(0, pm + 1) += b1 * H2_avg * std::real(vterm) * weight; - Right(pm + 1, 0) += b1 * H2_avg * std::real(vterm) * weight; + Return_t wfe = (HDsaved[pm] + (Dsaved[pm] - D_avg[pm]) * eloc_new) * weight; + Return_t wfd = (Dsaved[pm] - D_avg[pm]) * weight; + Return_t vterm = HDsaved[pm] * (eloc_new - curAvg_w) + + (Dsaved[pm] - D_avg[pm]) * eloc_new * (eloc_new - RealType(2.0) * curAvg_w); // Variance Left(0, pm + 1) += b2 * std::real(vterm) * weight; Left(pm + 1, 0) += b2 * std::real(vterm) * weight; @@ -776,30 +757,28 @@ QMCCostFunctionBatched::Return_rt QMCCostFunctionBatched::fillOverlapHamiltonian for (int pm2 = 0; pm2 < numParams; pm2++) { // Hamiltonian - Left(pm + 1, pm2 + 1) += std::real((1 - b2) * std::conj(wfd) * (HDsaved[pm2] + (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); + Left(pm + 1, pm2 + 1) += + std::real((1 - b2) * std::conj(wfd) * (HDsaved[pm2] + (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); // Overlap RealType ovlij = std::real(std::conj(wfd) * (Dsaved[pm2] - D_avg[pm2])); Right(pm + 1, pm2 + 1) += ovlij; // Variance - RealType varij = weight * std::real((HDsaved[pm] - RealType(2.0) * std::conj(Dsaved[pm] - D_avg[pm]) * eloc_new) * - (HDsaved[pm2] - RealType(2.0) * (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); + RealType varij = weight * + std::real((HDsaved[pm] - RealType(2.0) * std::conj(Dsaved[pm] - D_avg[pm]) * eloc_new) * + (HDsaved[pm2] - RealType(2.0) * (Dsaved[pm2] - D_avg[pm2]) * eloc_new)); Left(pm + 1, pm2 + 1) += b2 * (varij + V_avg * ovlij); - // H2 - Right(pm + 1, pm2 + 1) += b1 * H2_avg * varij; } } }; ParallelExecutor<> crowd_tasks; crowd_tasks(opt_num_crowds, constructMatrices, params_per_crowd, getNumParams(), Dsaved, HDsaved, weight, eloc_new, - H2_avg, V_avg, D_avg, b1, b2, curAvg_w, Left, Right); + V_avg, D_avg, w_beta, curAvg_w, Left, Right); } myComm->allreduce(Right); myComm->allreduce(Left); - Left(0, 0) = (1 - b2) * curAvg_w + b2 * V_avg; - Right(0, 0) = 1.0 + b1 * H2_avg * V_avg; - if (GEVType == "H2") - return H2_avg; + Left(0, 0) = (1 - w_beta) * curAvg_w + w_beta * V_avg; + Right(0, 0) = 1.0; return 1.0; } diff --git a/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.cpp b/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.cpp index 5c65f09cd7..6a759dd780 100644 --- a/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.cpp +++ b/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.cpp @@ -59,9 +59,7 @@ QMCFixedSampleLinearOptimize::QMCFixedSampleLinearOptimize(const ProjectData& pr bigChange(50), exp0(-16), stepsize(0.25), - GEVtype("mixed"), StabilizerMethod("best"), - GEVSplit("no"), bestShift_i(-1.0), bestShift_s(-1.0), shift_i_input(0.01), @@ -176,8 +174,6 @@ QMCFixedSampleLinearOptimize::QMCFixedSampleLinearOptimize(const ProjectData& pr // m_param.add(quadstep,"quadstep"); // m_param.add(stepsize,"stepsize"); // m_param.add(exp1,"exp1"); - // m_param.add(GEVtype,"GEVMethod"); - // m_param.add(GEVSplit,"GEVSplit"); // m_param.add(StabilizerMethod,"StabilizerMethod"); // m_param.add(LambdaMax,"LambdaMax"); //Set parameters for line minimization: diff --git a/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.h b/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.h index 54749d0684..2a4001f208 100644 --- a/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.h +++ b/src/QMCDrivers/WFOpt/QMCFixedSampleLinearOptimize.h @@ -136,7 +136,7 @@ class QMCFixedSampleLinearOptimize : public QMCDriver, public LinearMethod, priv int nstabilizers; RealType stabilizerScale, bigChange, exp0, exp1, stepsize, savedQuadstep; - std::string GEVtype, StabilizerMethod, GEVSplit; + std::string StabilizerMethod; RealType w_beta; /// number of previous steps to orthogonalize to. int eigCG;