From 0d4e31ca58ff63903b8832838ab7eba3b0f6a4a6 Mon Sep 17 00:00:00 2001 From: James Harrison Date: Mon, 30 Oct 2017 15:46:50 +0000 Subject: [PATCH] QedFVol: Calculate phase factors for momentum projections once per configuration only. --- extras/Hadrons/Modules/MScalar/ChargedProp.cc | 81 ----------- extras/Hadrons/Modules/MScalar/ScalarVP.cc | 131 ++++-------------- 2 files changed, 28 insertions(+), 184 deletions(-) diff --git a/extras/Hadrons/Modules/MScalar/ChargedProp.cc b/extras/Hadrons/Modules/MScalar/ChargedProp.cc index 6ffda844..3aaf3d5e 100644 --- a/extras/Hadrons/Modules/MScalar/ChargedProp.cc +++ b/extras/Hadrons/Modules/MScalar/ChargedProp.cc @@ -230,87 +230,6 @@ void TChargedProp::execute(void) write(writer, "prop_Q", resultQ); write(writer, "prop_Sun", resultSun); write(writer, "prop_Tad", resultTad); - - // // Write full propagator - // buf = prop; - // for (unsigned int j = 0; j < env().getNd()-1; ++j) - // { - // for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - // { - // buf = buf*adj(*phase_[j]); - // } - // } - // sliceSum(buf, vecBuf, Tp); - // result.resize(vecBuf.size()); - // for (unsigned int t = 0; t < vecBuf.size(); ++t) - // { - // result[t] = TensorRemove(vecBuf[t]); - // } - // write(writer, "prop", result); - - // // Write free propagator - // buf = *prop0_; - // for (unsigned int j = 0; j < env().getNd()-1; ++j) - // { - // for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - // { - // buf = buf*adj(*phase_[j]); - // } - // } - // sliceSum(buf, vecBuf, Tp); - // for (unsigned int t = 0; t < vecBuf.size(); ++t) - // { - // result[t] = TensorRemove(vecBuf[t]); - // } - // write(writer, "prop_0", result); - - // // Write propagator O(q) term - // buf = propQ; - // for (unsigned int j = 0; j < env().getNd()-1; ++j) - // { - // for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - // { - // buf = buf*adj(*phase_[j]); - // } - // } - // sliceSum(buf, vecBuf, Tp); - // for (unsigned int t = 0; t < vecBuf.size(); ++t) - // { - // result[t] = TensorRemove(vecBuf[t]); - // } - // write(writer, "prop_Q", result); - - // // Write propagator sunset term - // buf = propSun; - // for (unsigned int j = 0; j < env().getNd()-1; ++j) - // { - // for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - // { - // buf = buf*adj(*phase_[j]); - // } - // } - // sliceSum(buf, vecBuf, Tp); - // for (unsigned int t = 0; t < vecBuf.size(); ++t) - // { - // result[t] = TensorRemove(vecBuf[t]); - // } - // write(writer, "prop_Sun", result); - - // // Write propagator tadpole term - // buf = propTad; - // for (unsigned int j = 0; j < env().getNd()-1; ++j) - // { - // for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - // { - // buf = buf*adj(*phase_[j]); - // } - // } - // sliceSum(buf, vecBuf, Tp); - // for (unsigned int t = 0; t < vecBuf.size(); ++t) - // { - // result[t] = TensorRemove(vecBuf[t]); - // } - // write(writer, "prop_Tad", result); } } diff --git a/extras/Hadrons/Modules/MScalar/ScalarVP.cc b/extras/Hadrons/Modules/MScalar/ScalarVP.cc index 4d923802..297a823d 100644 --- a/extras/Hadrons/Modules/MScalar/ScalarVP.cc +++ b/extras/Hadrons/Modules/MScalar/ScalarVP.cc @@ -144,13 +144,19 @@ void TScalarVP::execute(void) } // Open output files if necessary + std::vector vecBuf; + std::vector result; + ScalarField vpPhase(env().getGrid()); std::vector writer, writer0, writerD; + std::vector momphases; if (!par().output.empty()) { + LOG(Message) << "Preparing output files..." << std::endl; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { std::vector mom = strToVec(par().outputMom[i_p]); + // Open output files std::string filename = par().output + "_" + std::to_string(mom[0]) + std::to_string(mom[1]) + std::to_string(mom[2]) @@ -180,11 +186,20 @@ void TScalarVP::execute(void) write(*writer0[i_p], "mass", static_cast(env().getModule(par().scalarProp))->par().mass); write(*writerD[i_p], "charge", q); write(*writerD[i_p], "mass", static_cast(env().getModule(par().scalarProp))->par().mass); + + // Calculate phase factors + vpPhase = Complex(1.0,0.0); + for (unsigned int j = 0; j < env().getNd()-1; ++j) + { + for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) + { + vpPhase = vpPhase*(*phase_[j]); + } + } + vpPhase = adj(vpPhase); + momphases.push_back(vpPhase); } } - std::vector vecBuf; - std::vector result; - ScalarField vpPhase(env().getGrid()); // Do contractions for (unsigned int nu = 0; nu < env().getNd(); ++nu) @@ -207,18 +222,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = freeVpTensor[mu][nu]; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = freeVpTensor[mu][nu]*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -244,18 +250,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -281,18 +278,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -316,18 +304,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -350,18 +329,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -385,18 +355,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -419,18 +380,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -457,18 +409,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -494,18 +437,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = tmp_vp; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = tmp_vp*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t) @@ -521,18 +455,9 @@ void TScalarVP::execute(void) // Output if necessary if (!par().output.empty()) { - std::vector mom; for (unsigned int i_p = 0; i_p < par().outputMom.size(); ++i_p) { - mom = strToVec(par().outputMom[i_p]); - vpPhase = vpTensor[mu][nu]; - for (unsigned int j = 0; j < env().getNd()-1; ++j) - { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - vpPhase = vpPhase*adj(*phase_[j]); - } - } + vpPhase = vpTensor[mu][nu]*momphases[i_p]; sliceSum(vpPhase, vecBuf, Tp); result.resize(vecBuf.size()); for (unsigned int t = 0; t < vecBuf.size(); ++t)