From c645d33db5a6bab5e3255c0dd2287e7f26bf849c Mon Sep 17 00:00:00 2001 From: James Harrison Date: Fri, 3 Nov 2017 10:59:26 +0000 Subject: [PATCH] QedFVol: Redo optimisation of charged propagator, and fix I/O bug --- extras/Hadrons/Modules/MScalar/ChargedProp.cc | 131 +++++++----------- 1 file changed, 49 insertions(+), 82 deletions(-) diff --git a/extras/Hadrons/Modules/MScalar/ChargedProp.cc b/extras/Hadrons/Modules/MScalar/ChargedProp.cc index 1b901bf1..cb7e6c79 100644 --- a/extras/Hadrons/Modules/MScalar/ChargedProp.cc +++ b/extras/Hadrons/Modules/MScalar/ChargedProp.cc @@ -151,23 +151,25 @@ void TChargedProp::execute(void) buf = GFSrc; momD1(buf, fft); buf = -G*buf; - fft.FFT_all_dim(propQ, buf, FFT::backward); + fft.FFT_dim(propQ, buf, env().getNd()-1, FFT::backward); // G*momD1*G*momD1*G*F*Src (here buf = G*momD1*G*F*Src) buf = -buf; momD1(buf, fft); propSun = G*buf; - fft.FFT_all_dim(propSun, propSun, FFT::backward); + fft.FFT_dim(propSun, propSun, env().getNd()-1, FFT::backward); // -G*momD2*G*F*Src (momD2 = F*D2*Finv) buf = GFSrc; momD2(buf, fft); buf = -G*buf; - fft.FFT_all_dim(propTad, buf, FFT::backward); - - // full charged scalar propagator - prop = (*prop0_) + q*propQ + q*q*propSun + q*q*propTad; + fft.FFT_dim(propTad, buf, env().getNd()-1, FFT::backward); + // full charged scalar propagator + buf = GFSrc; + fft.FFT_dim(buf, buf, env().getNd()-1, FFT::backward); + prop = buf + q*propQ + q*q*propSun + q*q*propTad; + // OUTPUT IF NECESSARY if (!par().output.empty()) { @@ -184,94 +186,59 @@ void TChargedProp::execute(void) << filename << "'..." << std::endl; CorrWriter writer(filename); - std::vector vecBuf; - std::vector result; + // std::vector vecBuf; + std::vector result, result0, resultQ, resultSun, resultTad; + result.resize(env().getGrid()->_fdimensions[env().getNd()-1]); + result0.resize(env().getGrid()->_fdimensions[env().getNd()-1]); + resultQ.resize(env().getGrid()->_fdimensions[env().getNd()-1]); + resultSun.resize(env().getGrid()->_fdimensions[env().getNd()-1]); + resultTad.resize(env().getGrid()->_fdimensions[env().getNd()-1]); write(writer, "charge", q); write(writer, "mass", par().mass); - // Write full propagator - buf = prop; + TComplex site; + std::vector whichmom; + whichmom.resize(env().getNd()); + for (unsigned int j = 0; j < env().getNd()-1; ++j) { - for (unsigned int momcount = 0; momcount < mom[j]; ++momcount) - { - buf = buf*adj(*phase_[j]); - } + whichmom[j] = mom[j]; } - sliceSum(buf, vecBuf, Tp); - result.resize(vecBuf.size()); - for (unsigned int t = 0; t < vecBuf.size(); ++t) + + for (unsigned int t = 0; t < env().getGrid()->_fdimensions[env().getNd()-1]; ++t) { - result[t] = TensorRemove(vecBuf[t]); + whichmom[env().getNd()-1] = t; + // Write full propagator + peekSite(site, prop, whichmom); + result[t]=TensorRemove(site); + // Write free propagator + peekSite(site, buf, whichmom); + result0[t]=TensorRemove(site); + // Write propagator O(q) term + peekSite(site, propQ, whichmom); + resultQ[t]=TensorRemove(site); + // Write propagator sunset term + peekSite(site, propSun, whichmom); + resultSun[t]=TensorRemove(site); + // Write propagator tadpole term + peekSite(site, propTad, whichmom); + resultTad[t]=TensorRemove(site); } 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); + write(writer, "prop_0", result0); + write(writer, "prop_Q", resultQ); + write(writer, "prop_Sun", resultSun); + write(writer, "prop_Tad", resultTad); } } + + std::vector mask(env().getNd(),1); + mask[env().getNd()-1] = 0; + fft.FFT_dim_mask(prop, prop, mask, FFT::backward); + fft.FFT_dim_mask(propQ, propQ, mask, FFT::backward); + fft.FFT_dim_mask(propSun, propSun, mask, FFT::backward); + fft.FFT_dim_mask(propTad, propTad, mask, FFT::backward); } void TChargedProp::momD1(ScalarField &s, FFT &fft)