From 842754bea9f8c4ba42c9295854342c6857f061bf Mon Sep 17 00:00:00 2001 From: Antonin Portelli Date: Wed, 13 Dec 2017 19:41:41 +0000 Subject: [PATCH] Hadrons: most modules ported to the new interface, compiles but untested --- extras/Hadrons/Module.hpp | 4 +- extras/Hadrons/Modules.hpp | 38 +++++------ extras/Hadrons/Modules/MAction/DWF.hpp | 11 ++-- extras/Hadrons/Modules/MAction/Wilson.hpp | 1 + .../Hadrons/Modules/MContraction/Baryon.hpp | 32 +++++++--- .../Hadrons/Modules/MContraction/DiscLoop.hpp | 17 +++-- .../Hadrons/Modules/MContraction/Gamma3pt.hpp | 21 +++++-- extras/Hadrons/Modules/MContraction/Meson.hpp | 11 ++-- .../Modules/MContraction/WardIdentity.hpp | 36 ++++++++--- .../Modules/MContraction/WeakHamiltonian.hpp | 1 + .../MContraction/WeakHamiltonianEye.cc | 41 ++++++++---- .../MContraction/WeakHamiltonianNonEye.cc | 44 ++++++++----- .../MContraction/WeakNeutral4ptDisc.cc | 39 ++++++++---- extras/Hadrons/Modules/MFermion/GaugeProp.hpp | 16 ++--- extras/Hadrons/Modules/MGauge/Load.cc | 19 ++++-- extras/Hadrons/Modules/MGauge/Load.hpp | 1 + extras/Hadrons/Modules/MGauge/Random.cc | 16 ++++- extras/Hadrons/Modules/MGauge/Random.hpp | 1 + extras/Hadrons/Modules/MGauge/StochEm.cc | 29 +++++---- extras/Hadrons/Modules/MGauge/StochEm.hpp | 1 + extras/Hadrons/Modules/MGauge/Unit.cc | 1 + extras/Hadrons/Modules/MLoop/NoiseLoop.hpp | 18 ++++-- extras/Hadrons/Modules/MSink/Point.hpp | 13 ++-- extras/Hadrons/Modules/MSink/Smear.hpp | 23 ++++--- extras/Hadrons/Modules/MSource/Point.hpp | 11 ++-- .../Hadrons/Modules/MSource/SeqConserved.hpp | 17 +++-- extras/Hadrons/Modules/MSource/SeqGamma.hpp | 53 +++++++++++----- extras/Hadrons/Modules/MSource/Wall.hpp | 49 ++++++++++----- extras/Hadrons/Modules/MSource/Z2.hpp | 38 ++++++++--- .../Modules/MUtilities/TestSeqConserved.hpp | 44 ++++++++----- .../Modules/MUtilities/TestSeqGamma.hpp | 28 ++++++--- extras/Hadrons/modules.inc | 63 +++++++++---------- tests/hadrons/Test_hadrons.hpp | 40 ++++++------ 33 files changed, 504 insertions(+), 273 deletions(-) diff --git a/extras/Hadrons/Module.hpp b/extras/Hadrons/Module.hpp index 25c0ac05..390573d8 100644 --- a/extras/Hadrons/Module.hpp +++ b/extras/Hadrons/Module.hpp @@ -92,8 +92,8 @@ static ns##mod##ModuleRegistrar ns##mod##ModuleRegistrarInstance; #define envGet(type, name)\ *env().template getObject(name) -#define envGetTmp(type, name)\ -*env().template getObject(getName() + "_tmp_" + name) +#define envGetTmp(type, var)\ +type &var = *env().template getObject(getName() + "_tmp_" + #var) #define envHasType(type, name)\ env().template isObjectOfType(name) diff --git a/extras/Hadrons/Modules.hpp b/extras/Hadrons/Modules.hpp index bb574a14..61a20058 100644 --- a/extras/Hadrons/Modules.hpp +++ b/extras/Hadrons/Modules.hpp @@ -30,31 +30,31 @@ See the full license in the file "LICENSE" in the top level distribution directo #include #include -// #include -// #include -// #include +#include +#include +#include #include -// #include -// #include -// #include -// #include -// #include +#include +#include +#include +#include +#include #include -// #include -// #include -// #include +#include +#include +#include #include -// #include +#include // #include // #include // #include #include -// #include +#include #include #include -// #include -// #include -// #include -// #include -// #include -// #include +#include +#include +#include +#include +#include +#include diff --git a/extras/Hadrons/Modules/MAction/DWF.hpp b/extras/Hadrons/Modules/MAction/DWF.hpp index 91e4ec94..0cb9a4cb 100644 --- a/extras/Hadrons/Modules/MAction/DWF.hpp +++ b/extras/Hadrons/Modules/MAction/DWF.hpp @@ -119,12 +119,13 @@ void TDWF::setup(void) << std::endl; LOG(Message) << "Fermion boundary conditions: " << par().boundary << std::endl; + env().createGrid(par().Ls); - auto &U = envGet(LatticeGaugeField, par().gauge); - auto &g4 = *env().getGrid(); - auto &grb4 = *env().getRbGrid(); - auto &g5 = *env().getGrid(par().Ls); - auto &grb5 = *env().getRbGrid(par().Ls); + auto &U = envGet(LatticeGaugeField, par().gauge); + auto &g4 = *env().getGrid(); + auto &grb4 = *env().getRbGrid(); + auto &g5 = *env().getGrid(par().Ls); + auto &grb5 = *env().getRbGrid(par().Ls); std::vector boundary = strToVec(par().boundary); typename DomainWallFermion::ImplParams implParams(boundary); envCreateDerived(FMat, DomainWallFermion, getName(), par().Ls, U, g5, diff --git a/extras/Hadrons/Modules/MAction/Wilson.hpp b/extras/Hadrons/Modules/MAction/Wilson.hpp index 1ca3bf59..a6b3f0d6 100644 --- a/extras/Hadrons/Modules/MAction/Wilson.hpp +++ b/extras/Hadrons/Modules/MAction/Wilson.hpp @@ -115,6 +115,7 @@ void TWilson::setup(void) << " using gauge field '" << par().gauge << "'" << std::endl; LOG(Message) << "Fermion boundary conditions: " << par().boundary << std::endl; + auto &U = envGet(LatticeGaugeField, par().gauge); auto &grid = *env().getGrid(); auto &gridRb = *env().getRbGrid(); diff --git a/extras/Hadrons/Modules/MContraction/Baryon.hpp b/extras/Hadrons/Modules/MContraction/Baryon.hpp index da927391..28f6aa51 100644 --- a/extras/Hadrons/Modules/MContraction/Baryon.hpp +++ b/extras/Hadrons/Modules/MContraction/Baryon.hpp @@ -71,8 +71,11 @@ public: virtual ~TBaryon(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: + // setup + virtual void setup(void); // execution virtual void execute(void); }; @@ -97,14 +100,29 @@ std::vector TBaryon::getInput(void) return input; } +template +std::vector TBaryon::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TBaryon::getOutput(void) { - std::vector out = {getName()}; + std::vector out = {}; return out; } +// setup /////////////////////////////////////////////////////////////////////// +template +void TBaryon::setup(void) +{ + envTmpLat(LatticeComplex, "c"); +} + // execution /////////////////////////////////////////////////////////////////// template void TBaryon::execute(void) @@ -113,12 +131,12 @@ void TBaryon::execute(void) << " quarks '" << par().q1 << "', '" << par().q2 << "', and '" << par().q3 << "'" << std::endl; - CorrWriter writer(par().output); - PropagatorField1 &q1 = *env().template getObject(par().q1); - PropagatorField2 &q2 = *env().template getObject(par().q2); - PropagatorField3 &q3 = *env().template getObject(par().q2); - LatticeComplex c(env().getGrid()); - Result result; + CorrWriter writer(par().output); + auto &q1 = envGet(PropagatorField1, par().q1); + auto &q2 = envGet(PropagatorField2, par().q2); + auto &q3 = envGet(PropagatorField3, par().q2); + envGetTmp(LatticeComplex, c); + Result result; // FIXME: do contractions diff --git a/extras/Hadrons/Modules/MContraction/DiscLoop.hpp b/extras/Hadrons/Modules/MContraction/DiscLoop.hpp index f8da3943..c0fbe296 100644 --- a/extras/Hadrons/Modules/MContraction/DiscLoop.hpp +++ b/extras/Hadrons/Modules/MContraction/DiscLoop.hpp @@ -67,6 +67,7 @@ public: virtual ~TDiscLoop(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -95,10 +96,18 @@ std::vector TDiscLoop::getInput(void) return in; } +template +std::vector TDiscLoop::getReference(void) +{ + std::vector out = {}; + + return out; +} + template std::vector TDiscLoop::getOutput(void) { - std::vector out = {getName()}; + std::vector out = {}; return out; } @@ -107,7 +116,7 @@ std::vector TDiscLoop::getOutput(void) template void TDiscLoop::setup(void) { - + envTmpLat(LatticeComplex, "c"); } // execution /////////////////////////////////////////////////////////////////// @@ -119,12 +128,12 @@ void TDiscLoop::execute(void) << " insertion." << std::endl; CorrWriter writer(par().output); - PropagatorField &q_loop = *env().template getObject(par().q_loop); - LatticeComplex c(env().getGrid()); + auto &q_loop = envGet(PropagatorField, par().q_loop); Gamma gamma(par().gamma); std::vector buf; Result result; + envGetTmp(LatticeComplex, c); c = trace(gamma*q_loop); sliceSum(c, buf, Tp); diff --git a/extras/Hadrons/Modules/MContraction/Gamma3pt.hpp b/extras/Hadrons/Modules/MContraction/Gamma3pt.hpp index a8653186..4a6baf3e 100644 --- a/extras/Hadrons/Modules/MContraction/Gamma3pt.hpp +++ b/extras/Hadrons/Modules/MContraction/Gamma3pt.hpp @@ -98,6 +98,7 @@ public: virtual ~TGamma3pt(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -126,10 +127,18 @@ std::vector TGamma3pt::getInput(void) return in; } +template +std::vector TGamma3pt::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TGamma3pt::getOutput(void) { - std::vector out = {getName()}; + std::vector out = {}; return out; } @@ -138,7 +147,7 @@ std::vector TGamma3pt::getOutput(void) template void TGamma3pt::setup(void) { - + envTmpLat(LatticeComplex, "c"); } // execution /////////////////////////////////////////////////////////////////// @@ -153,10 +162,9 @@ void TGamma3pt::execute(void) // Initialise variables. q2 and q3 are normal propagators, q1 may be // sink smeared. CorrWriter writer(par().output); - SlicedPropagator1 &q1 = *env().template getObject(par().q1); - PropagatorField2 &q2 = *env().template getObject(par().q2); - PropagatorField3 &q3 = *env().template getObject(par().q3); - LatticeComplex c(env().getGrid()); + auto &q1 = envGet(SlicedPropagator1, par().q1); + auto &q2 = envGet(PropagatorField2, par().q2); + auto &q3 = envGet(PropagatorField2, par().q3); Gamma g5(Gamma::Algebra::Gamma5); Gamma gamma(par().gamma); std::vector buf; @@ -165,6 +173,7 @@ void TGamma3pt::execute(void) // Extract relevant timeslice of sinked propagator q1, then contract & // sum over all spacial positions of gamma insertion. SitePropagator1 q1Snk = q1[par().tSnk]; + envGetTmp(LatticeComplex, c); c = trace(g5*q1Snk*adj(q2)*(g5*gamma)*q3); sliceSum(c, buf, Tp); diff --git a/extras/Hadrons/Modules/MContraction/Meson.hpp b/extras/Hadrons/Modules/MContraction/Meson.hpp index 3c179d44..1fd86d3a 100644 --- a/extras/Hadrons/Modules/MContraction/Meson.hpp +++ b/extras/Hadrons/Modules/MContraction/Meson.hpp @@ -161,6 +161,7 @@ void TMeson::parseGammaString(std::vector &gammaList) // Parse individual contractions from input string. gammaList = strToVec(par().gammas); } + envTmpLat(LatticeComplex, "c"); } // execution /////////////////////////////////////////////////////////////////// @@ -192,8 +193,8 @@ void TMeson::execute(void) if (envHasType(SlicedPropagator1, par().q1) and envHasType(SlicedPropagator2, par().q2)) { - SlicedPropagator1 &q1 = envGet(SlicedPropagator1, par().q1); - SlicedPropagator2 &q2 = envGet(SlicedPropagator2, par().q2); + auto &q1 = envGet(SlicedPropagator1, par().q1); + auto &q2 = envGet(SlicedPropagator2, par().q2); LOG(Message) << "(propagator already sinked)" << std::endl; for (unsigned int i = 0; i < result.size(); ++i) @@ -209,10 +210,10 @@ void TMeson::execute(void) } else { - PropagatorField1 &q1 = envGet(PropagatorField1, par().q1); - PropagatorField2 &q2 = envGet(PropagatorField2, par().q2); - LatticeComplex c(env().getGrid()); + auto &q1 = envGet(PropagatorField1, par().q1); + auto &q2 = envGet(PropagatorField2, par().q2); + envGetTmp(LatticeComplex, c); LOG(Message) << "(using sink '" << par().sink << "')" << std::endl; for (unsigned int i = 0; i < result.size(); ++i) { diff --git a/extras/Hadrons/Modules/MContraction/WardIdentity.hpp b/extras/Hadrons/Modules/MContraction/WardIdentity.hpp index 90922c27..c92c7243 100644 --- a/extras/Hadrons/Modules/MContraction/WardIdentity.hpp +++ b/extras/Hadrons/Modules/MContraction/WardIdentity.hpp @@ -73,6 +73,7 @@ public: virtual ~TWardIdentity(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -103,10 +104,18 @@ std::vector TWardIdentity::getInput(void) return in; } +template +std::vector TWardIdentity::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TWardIdentity::getOutput(void) { - std::vector out = {getName()}; + std::vector out = {}; return out; } @@ -120,6 +129,15 @@ void TWardIdentity::setup(void) { HADRON_ERROR(Size, "Ls mismatch between quark action and propagator"); } + envTmpLat(PropagatorField, "tmp"); + envTmpLat(PropagatorField, "vector_WI"); + if (par().test_axial) + { + envTmpLat(PropagatorField, "psi"); + envTmpLat(LatticeComplex, "PP"); + envTmpLat(LatticeComplex, "axial_defect"); + envTmpLat(LatticeComplex, "PJ5q"); + } } // execution /////////////////////////////////////////////////////////////////// @@ -129,12 +147,13 @@ void TWardIdentity::execute(void) LOG(Message) << "Performing Ward Identity checks for quark '" << par().q << "'." << std::endl; - PropagatorField tmp(env().getGrid()), vector_WI(env().getGrid()); - PropagatorField &q = *env().template getObject(par().q); - FMat &act = *(env().template getObject(par().action)); - Gamma g5(Gamma::Algebra::Gamma5); + auto &q = envGet(PropagatorField, par().q); + auto &act = envGet(FMat, par().action); + Gamma g5(Gamma::Algebra::Gamma5); // Compute D_mu V_mu, D here is backward derivative. + envGetTmp(PropagatorField, tmp); + envGetTmp(PropagatorField, vector_WI); vector_WI = zero; for (unsigned int mu = 0; mu < Nd; ++mu) { @@ -149,9 +168,10 @@ void TWardIdentity::execute(void) if (par().test_axial) { - PropagatorField psi(env().getGrid()); - LatticeComplex PP(env().getGrid()), axial_defect(env().getGrid()), - PJ5q(env().getGrid()); + envGetTmp(PropagatorField, psi); + envGetTmp(LatticeComplex, PP); + envGetTmp(LatticeComplex, axial_defect); + envGetTmp(LatticeComplex, PJ5q); std::vector axial_buf; // Compute , D is backwards derivative. diff --git a/extras/Hadrons/Modules/MContraction/WeakHamiltonian.hpp b/extras/Hadrons/Modules/MContraction/WeakHamiltonian.hpp index 7df40370..2b53c87a 100644 --- a/extras/Hadrons/Modules/MContraction/WeakHamiltonian.hpp +++ b/extras/Hadrons/Modules/MContraction/WeakHamiltonian.hpp @@ -99,6 +99,7 @@ public:\ virtual ~T##modname(void) = default;\ /* dependency relation */ \ virtual std::vector getInput(void);\ + virtual std::vector getReference(void);\ virtual std::vector getOutput(void);\ public:\ std::vector VA_label = {"V", "A"};\ diff --git a/extras/Hadrons/Modules/MContraction/WeakHamiltonianEye.cc b/extras/Hadrons/Modules/MContraction/WeakHamiltonianEye.cc index 314b080a..7a73a7e3 100644 --- a/extras/Hadrons/Modules/MContraction/WeakHamiltonianEye.cc +++ b/extras/Hadrons/Modules/MContraction/WeakHamiltonianEye.cc @@ -74,9 +74,16 @@ std::vector TWeakHamiltonianEye::getInput(void) return in; } +std::vector TWeakHamiltonianEye::getReference(void) +{ + std::vector out = {}; + + return out; +} + std::vector TWeakHamiltonianEye::getOutput(void) { - std::vector out = {getName()}; + std::vector out = {}; return out; } @@ -84,7 +91,15 @@ std::vector TWeakHamiltonianEye::getOutput(void) // setup /////////////////////////////////////////////////////////////////////// void TWeakHamiltonianEye::setup(void) { + unsigned int ndim = env().getNd(); + envTmpLat(LatticeComplex, "expbuf"); + envTmpLat(PropagatorField, "tmp1"); + envTmpLat(LatticeComplex, "tmp2"); + envTmp(std::vector, "S_body", 1, ndim, PropagatorField(env().getGrid())); + envTmp(std::vector, "S_loop", 1, ndim, PropagatorField(env().getGrid())); + envTmp(std::vector, "E_body", 1, ndim, LatticeComplex(env().getGrid())); + envTmp(std::vector, "E_loop", 1, ndim, LatticeComplex(env().getGrid())); } // execution /////////////////////////////////////////////////////////////////// @@ -96,22 +111,22 @@ void TWeakHamiltonianEye::execute(void) << "'." << std::endl; CorrWriter writer(par().output); - SlicedPropagator &q1 = *env().template getObject(par().q1); - PropagatorField &q2 = *env().template getObject(par().q2); - PropagatorField &q3 = *env().template getObject(par().q3); - PropagatorField &q4 = *env().template getObject(par().q4); - Gamma g5 = Gamma(Gamma::Algebra::Gamma5); - LatticeComplex expbuf(env().getGrid()); + auto &q1 = envGet(SlicedPropagator, par().q1); + auto &q2 = envGet(PropagatorField, par().q2); + auto &q3 = envGet(PropagatorField, par().q3); + auto &q4 = envGet(PropagatorField, par().q4); + Gamma g5 = Gamma(Gamma::Algebra::Gamma5); std::vector corrbuf; std::vector result(n_eye_diag); unsigned int ndim = env().getNd(); - PropagatorField tmp1(env().getGrid()); - LatticeComplex tmp2(env().getGrid()); - std::vector S_body(ndim, tmp1); - std::vector S_loop(ndim, tmp1); - std::vector E_body(ndim, tmp2); - std::vector E_loop(ndim, tmp2); + envGetTmp(LatticeComplex, expbuf); + envGetTmp(PropagatorField, tmp1); + envGetTmp(LatticeComplex, tmp2); + envGetTmp(std::vector, S_body); + envGetTmp(std::vector, S_loop); + envGetTmp(std::vector, E_body); + envGetTmp(std::vector, E_loop); // Get sink timeslice of q1. SitePropagator q1Snk = q1[par().tSnk]; diff --git a/extras/Hadrons/Modules/MContraction/WeakHamiltonianNonEye.cc b/extras/Hadrons/Modules/MContraction/WeakHamiltonianNonEye.cc index 2c4df68a..c333713d 100644 --- a/extras/Hadrons/Modules/MContraction/WeakHamiltonianNonEye.cc +++ b/extras/Hadrons/Modules/MContraction/WeakHamiltonianNonEye.cc @@ -74,9 +74,15 @@ std::vector TWeakHamiltonianNonEye::getInput(void) return in; } +std::vector TWeakHamiltonianNonEye::getReference(void) +{ + std::vector out = {}; + + return out; +} std::vector TWeakHamiltonianNonEye::getOutput(void) { - std::vector out = {getName()}; + std::vector out = {}; return out; } @@ -84,7 +90,15 @@ std::vector TWeakHamiltonianNonEye::getOutput(void) // setup /////////////////////////////////////////////////////////////////////// void TWeakHamiltonianNonEye::setup(void) { + unsigned int ndim = env().getNd(); + envTmpLat(LatticeComplex, "expbuf"); + envTmpLat(PropagatorField, "tmp1"); + envTmpLat(LatticeComplex, "tmp2"); + envTmp(std::vector, "C_i_side_loop", 1, ndim, PropagatorField(env().getGrid())); + envTmp(std::vector, "C_f_side_loop", 1, ndim, PropagatorField(env().getGrid())); + envTmp(std::vector, "W_i_side_loop", 1, ndim, LatticeComplex(env().getGrid())); + envTmp(std::vector, "W_f_side_loop", 1, ndim, LatticeComplex(env().getGrid())); } // execution /////////////////////////////////////////////////////////////////// @@ -95,23 +109,23 @@ void TWeakHamiltonianNonEye::execute(void) << par().q2 << ", '" << par().q3 << "' and '" << par().q4 << "'." << std::endl; - CorrWriter writer(par().output); - PropagatorField &q1 = *env().template getObject(par().q1); - PropagatorField &q2 = *env().template getObject(par().q2); - PropagatorField &q3 = *env().template getObject(par().q3); - PropagatorField &q4 = *env().template getObject(par().q4); - Gamma g5 = Gamma(Gamma::Algebra::Gamma5); - LatticeComplex expbuf(env().getGrid()); + CorrWriter writer(par().output); + auto &q1 = envGet(PropagatorField, par().q1); + auto &q2 = envGet(PropagatorField, par().q2); + auto &q3 = envGet(PropagatorField, par().q3); + auto &q4 = envGet(PropagatorField, par().q4); + Gamma g5 = Gamma(Gamma::Algebra::Gamma5); std::vector corrbuf; std::vector result(n_noneye_diag); - unsigned int ndim = env().getNd(); + unsigned int ndim = env().getNd(); - PropagatorField tmp1(env().getGrid()); - LatticeComplex tmp2(env().getGrid()); - std::vector C_i_side_loop(ndim, tmp1); - std::vector C_f_side_loop(ndim, tmp1); - std::vector W_i_side_loop(ndim, tmp2); - std::vector W_f_side_loop(ndim, tmp2); + envGetTmp(LatticeComplex, expbuf); + envGetTmp(PropagatorField, tmp1); + envGetTmp(LatticeComplex, tmp2); + envGetTmp(std::vector, C_i_side_loop); + envGetTmp(std::vector, C_f_side_loop); + envGetTmp(std::vector, W_i_side_loop); + envGetTmp(std::vector, W_f_side_loop); // Setup for C-type contractions. for (int mu = 0; mu < ndim; ++mu) diff --git a/extras/Hadrons/Modules/MContraction/WeakNeutral4ptDisc.cc b/extras/Hadrons/Modules/MContraction/WeakNeutral4ptDisc.cc index 6685f292..e0f07f6c 100644 --- a/extras/Hadrons/Modules/MContraction/WeakNeutral4ptDisc.cc +++ b/extras/Hadrons/Modules/MContraction/WeakNeutral4ptDisc.cc @@ -76,9 +76,16 @@ std::vector TWeakNeutral4ptDisc::getInput(void) return in; } +std::vector TWeakNeutral4ptDisc::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + std::vector TWeakNeutral4ptDisc::getOutput(void) { - std::vector out = {getName()}; + std::vector out = {}; return out; } @@ -86,7 +93,13 @@ std::vector TWeakNeutral4ptDisc::getOutput(void) // setup /////////////////////////////////////////////////////////////////////// void TWeakNeutral4ptDisc::setup(void) { + unsigned int ndim = env().getNd(); + envTmpLat(LatticeComplex, "expbuf"); + envTmpLat(PropagatorField, "tmp"); + envTmpLat(LatticeComplex, "curr"); + envTmp(std::vector, "meson", 1, ndim, PropagatorField(env().getGrid())); + envTmp(std::vector, "loop", 1, ndim, PropagatorField(env().getGrid())); } // execution /////////////////////////////////////////////////////////////////// @@ -97,21 +110,21 @@ void TWeakNeutral4ptDisc::execute(void) << par().q2 << ", '" << par().q3 << "' and '" << par().q4 << "'." << std::endl; - CorrWriter writer(par().output); - PropagatorField &q1 = *env().template getObject(par().q1); - PropagatorField &q2 = *env().template getObject(par().q2); - PropagatorField &q3 = *env().template getObject(par().q3); - PropagatorField &q4 = *env().template getObject(par().q4); - Gamma g5 = Gamma(Gamma::Algebra::Gamma5); - LatticeComplex expbuf(env().getGrid()); + CorrWriter writer(par().output); + auto &q1 = envGet(PropagatorField, par().q1); + auto &q2 = envGet(PropagatorField, par().q2); + auto &q3 = envGet(PropagatorField, par().q3); + auto &q4 = envGet(PropagatorField, par().q4); + Gamma g5 = Gamma(Gamma::Algebra::Gamma5); std::vector corrbuf; std::vector result(n_neut_disc_diag); - unsigned int ndim = env().getNd(); + unsigned int ndim = env().getNd(); - PropagatorField tmp(env().getGrid()); - std::vector meson(ndim, tmp); - std::vector loop(ndim, tmp); - LatticeComplex curr(env().getGrid()); + envGetTmp(LatticeComplex, expbuf); + envGetTmp(PropagatorField, tmp); + envGetTmp(LatticeComplex, curr); + envGetTmp(std::vector, meson); + envGetTmp(std::vector, loop); // Setup for type 1 contractions. for (int mu = 0; mu < ndim; ++mu) diff --git a/extras/Hadrons/Modules/MFermion/GaugeProp.hpp b/extras/Hadrons/Modules/MFermion/GaugeProp.hpp index 4d08841d..e77df287 100644 --- a/extras/Hadrons/Modules/MFermion/GaugeProp.hpp +++ b/extras/Hadrons/Modules/MFermion/GaugeProp.hpp @@ -154,21 +154,21 @@ void TGaugeProp::execute(void) LOG(Message) << "Computing quark propagator '" << getName() << "'" << std::endl; - FermionField &source = envGetTmp(FermionField, "source"); - FermionField &sol = envGetTmp(FermionField, "sol"); - FermionField &tmp = envGetTmp(FermionField, "tmp"); - std::string propName = (Ls_ == 1) ? getName() : (getName() + "_5d"); - PropagatorField &prop = envGet(PropagatorField, propName); - PropagatorField &fullSrc = envGet(PropagatorField, par().source); - SolverFn &solver = envGet(SolverFn, par().solver); + std::string propName = (Ls_ == 1) ? getName() : (getName() + "_5d"); + auto &prop = envGet(PropagatorField, propName); + auto &fullSrc = envGet(PropagatorField, par().source); + auto &solver = envGet(SolverFn, par().solver); + envGetTmp(FermionField, source); + envGetTmp(FermionField, sol); + envGetTmp(FermionField, tmp); LOG(Message) << "Inverting using solver '" << par().solver << "' on source '" << par().source << "'" << std::endl; for (unsigned int s = 0; s < Ns; ++s) for (unsigned int c = 0; c < Nc; ++c) { LOG(Message) << "Inversion for spin= " << s << ", color= " << c - << std::endl; + << std::endl; // source conversion for 4D sources if (!env().isObject5d(par().source)) { diff --git a/extras/Hadrons/Modules/MGauge/Load.cc b/extras/Hadrons/Modules/MGauge/Load.cc index 062e7e98..c2fd49de 100644 --- a/extras/Hadrons/Modules/MGauge/Load.cc +++ b/extras/Hadrons/Modules/MGauge/Load.cc @@ -49,6 +49,13 @@ std::vector TLoad::getInput(void) return in; } +std::vector TLoad::getReference(void) +{ + std::vector ref; + + return ref; +} + std::vector TLoad::getOutput(void) { std::vector out = {getName()}; @@ -59,19 +66,19 @@ std::vector TLoad::getOutput(void) // setup /////////////////////////////////////////////////////////////////////// void TLoad::setup(void) { - env().registerLattice(getName()); + envCreateLat(LatticeGaugeField, getName()); } // execution /////////////////////////////////////////////////////////////////// void TLoad::execute(void) { - FieldMetaData header; - std::string fileName = par().file + "." - + std::to_string(env().getTrajectory()); - + FieldMetaData header; + std::string fileName = par().file + "." + + std::to_string(vm().getTrajectory()); LOG(Message) << "Loading NERSC configuration from file '" << fileName << "'" << std::endl; - LatticeGaugeField &U = *env().createLattice(getName()); + + auto &U = envGet(LatticeGaugeField, getName()); NerscIO::readConfiguration(U, header, fileName); LOG(Message) << "NERSC header:" << std::endl; dump_meta_data(header, LOG(Message)); diff --git a/extras/Hadrons/Modules/MGauge/Load.hpp b/extras/Hadrons/Modules/MGauge/Load.hpp index a338af79..a967d714 100644 --- a/extras/Hadrons/Modules/MGauge/Load.hpp +++ b/extras/Hadrons/Modules/MGauge/Load.hpp @@ -57,6 +57,7 @@ public: virtual ~TLoad(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup diff --git a/extras/Hadrons/Modules/MGauge/Random.cc b/extras/Hadrons/Modules/MGauge/Random.cc index c10fdfc3..fdb0d145 100644 --- a/extras/Hadrons/Modules/MGauge/Random.cc +++ b/extras/Hadrons/Modules/MGauge/Random.cc @@ -44,7 +44,16 @@ TRandom::TRandom(const std::string name) // dependencies/products /////////////////////////////////////////////////////// std::vector TRandom::getInput(void) { - return std::vector(); + std::vector in; + + return in; +} + +std::vector TRandom::getReference(void) +{ + std::vector ref; + + return ref; } std::vector TRandom::getOutput(void) @@ -57,13 +66,14 @@ std::vector TRandom::getOutput(void) // setup /////////////////////////////////////////////////////////////////////// void TRandom::setup(void) { - env().registerLattice(getName()); + envCreateLat(LatticeGaugeField, getName()); } // execution /////////////////////////////////////////////////////////////////// void TRandom::execute(void) { LOG(Message) << "Generating random gauge configuration" << std::endl; - LatticeGaugeField &U = *env().createLattice(getName()); + + auto &U = envGet(LatticeGaugeField, getName()); SU3::HotConfiguration(*env().get4dRng(), U); } diff --git a/extras/Hadrons/Modules/MGauge/Random.hpp b/extras/Hadrons/Modules/MGauge/Random.hpp index a07130e4..30525113 100644 --- a/extras/Hadrons/Modules/MGauge/Random.hpp +++ b/extras/Hadrons/Modules/MGauge/Random.hpp @@ -50,6 +50,7 @@ public: virtual ~TRandom(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup diff --git a/extras/Hadrons/Modules/MGauge/StochEm.cc b/extras/Hadrons/Modules/MGauge/StochEm.cc index c7a9fc4f..a878ae2f 100644 --- a/extras/Hadrons/Modules/MGauge/StochEm.cc +++ b/extras/Hadrons/Modules/MGauge/StochEm.cc @@ -47,6 +47,13 @@ std::vector TStochEm::getInput(void) return in; } +std::vector TStochEm::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + std::vector TStochEm::getOutput(void) { std::vector out = {getName()}; @@ -57,32 +64,28 @@ std::vector TStochEm::getOutput(void) // setup /////////////////////////////////////////////////////////////////////// void TStochEm::setup(void) { - if (!env().hasRegisteredObject("_" + getName() + "_weight")) + if (!env().hasCreatedObject("_" + getName() + "_weight")) { - env().registerLattice("_" + getName() + "_weight"); + envCacheLat(EmComp, "_" + getName() + "_weight"); } - env().registerLattice(getName()); + envCreateLat(EmField, getName()); } // execution /////////////////////////////////////////////////////////////////// void TStochEm::execute(void) { + LOG(Message) << "Generating stochatic EM potential..." << std::endl; + PhotonR photon(par().gauge, par().zmScheme); - EmField &a = *env().createLattice(getName()); - EmComp *w; + auto &a = envGet(EmField, getName()); + auto &w = envGet(EmComp, "_" + getName() + "_weight"); if (!env().hasCreatedObject("_" + getName() + "_weight")) { LOG(Message) << "Caching stochatic EM potential weight (gauge: " << par().gauge << ", zero-mode scheme: " << par().zmScheme << ")..." << std::endl; - w = env().createLattice("_" + getName() + "_weight"); - photon.StochasticWeight(*w); + photon.StochasticWeight(w); } - else - { - w = env().getObject("_" + getName() + "_weight"); - } - LOG(Message) << "Generating stochatic EM potential..." << std::endl; - photon.StochasticField(a, *env().get4dRng(), *w); + photon.StochasticField(a, *env().get4dRng(), w); } diff --git a/extras/Hadrons/Modules/MGauge/StochEm.hpp b/extras/Hadrons/Modules/MGauge/StochEm.hpp index bacb5172..efc2e39b 100644 --- a/extras/Hadrons/Modules/MGauge/StochEm.hpp +++ b/extras/Hadrons/Modules/MGauge/StochEm.hpp @@ -59,6 +59,7 @@ public: virtual ~TStochEm(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup diff --git a/extras/Hadrons/Modules/MGauge/Unit.cc b/extras/Hadrons/Modules/MGauge/Unit.cc index bc05a785..af31f124 100644 --- a/extras/Hadrons/Modules/MGauge/Unit.cc +++ b/extras/Hadrons/Modules/MGauge/Unit.cc @@ -71,6 +71,7 @@ void TUnit::setup(void) void TUnit::execute(void) { LOG(Message) << "Creating unit gauge configuration" << std::endl; + auto &U = envGet(LatticeGaugeField, getName()); SU3::ColdConfiguration(*env().get4dRng(), U); } diff --git a/extras/Hadrons/Modules/MLoop/NoiseLoop.hpp b/extras/Hadrons/Modules/MLoop/NoiseLoop.hpp index 1f40dd48..0feb5efb 100644 --- a/extras/Hadrons/Modules/MLoop/NoiseLoop.hpp +++ b/extras/Hadrons/Modules/MLoop/NoiseLoop.hpp @@ -73,6 +73,7 @@ public: virtual ~TNoiseLoop(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -101,6 +102,15 @@ std::vector TNoiseLoop::getInput(void) return in; } + +template +std::vector TNoiseLoop::getReference(void) +{ + std::vector out = {}; + + return out; +} + template std::vector TNoiseLoop::getOutput(void) { @@ -113,16 +123,16 @@ std::vector TNoiseLoop::getOutput(void) template void TNoiseLoop::setup(void) { - env().template registerLattice(getName()); + envCreateLat(PropagatorField, getName()); } // execution /////////////////////////////////////////////////////////////////// template void TNoiseLoop::execute(void) { - PropagatorField &loop = *env().template createLattice(getName()); - PropagatorField &q = *env().template getObject(par().q); - PropagatorField &eta = *env().template getObject(par().eta); + auto &loop = envGet(PropagatorField, getName()); + auto &q = envGet(PropagatorField, par().q); + auto &eta = envGet(PropagatorField, par().eta); loop = q*adj(eta); } diff --git a/extras/Hadrons/Modules/MSink/Point.hpp b/extras/Hadrons/Modules/MSink/Point.hpp index 16b89434..42cae4f6 100644 --- a/extras/Hadrons/Modules/MSink/Point.hpp +++ b/extras/Hadrons/Modules/MSink/Point.hpp @@ -122,18 +122,19 @@ void TPoint::setup(void) // execution /////////////////////////////////////////////////////////////////// template void TPoint::execute(void) -{ - std::vector p = strToVec(par().mom); - LatticeComplex &ph = envGet(LatticeComplex, momphName_); - Complex i(0.0,1.0); - +{ LOG(Message) << "Setting up point sink function for momentum [" << par().mom << "]" << std::endl; + auto &ph = envGet(LatticeComplex, momphName_); + if (!hasPhase_) { - LatticeComplex &coor = envGetTmp(LatticeComplex, "coor"); + Complex i(0.0,1.0); + std::vector p; + envGetTmp(LatticeComplex, coor); + p = strToVec(par().mom); ph = zero; for(unsigned int mu = 0; mu < env().getNd(); mu++) { diff --git a/extras/Hadrons/Modules/MSink/Smear.hpp b/extras/Hadrons/Modules/MSink/Smear.hpp index b51d2f49..03cc861a 100644 --- a/extras/Hadrons/Modules/MSink/Smear.hpp +++ b/extras/Hadrons/Modules/MSink/Smear.hpp @@ -61,6 +61,7 @@ public: virtual ~TSmear(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -89,6 +90,14 @@ std::vector TSmear::getInput(void) return in; } +template +std::vector TSmear::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TSmear::getOutput(void) { @@ -101,9 +110,7 @@ std::vector TSmear::getOutput(void) template void TSmear::setup(void) { - unsigned int nt = env().getDim(Tp); - unsigned int size = nt * sizeof(SitePropagator); - env().registerObject(getName(), size); + envCreate(SlicedPropagator, getName(), 1, env().getDim(Tp)); } // execution /////////////////////////////////////////////////////////////////// @@ -114,11 +121,11 @@ void TSmear::execute(void) << "' using sink function '" << par().sink << "'." << std::endl; - SinkFn &sink = *env().template getObject(par().sink); - PropagatorField &q = *env().template getObject(par().q); - SlicedPropagator *out = new SlicedPropagator(env().getDim(Tp)); - *out = sink(q); - env().setObject(getName(), out); + auto &sink = envGet(SinkFn, par().sink); + auto &q = envGet(PropagatorField, par().q); + auto &out = envGet(SlicedPropagator, getName()); + + out = sink(q); } END_MODULE_NAMESPACE diff --git a/extras/Hadrons/Modules/MSource/Point.hpp b/extras/Hadrons/Modules/MSource/Point.hpp index 3fab41c0..6470c77f 100644 --- a/extras/Hadrons/Modules/MSource/Point.hpp +++ b/extras/Hadrons/Modules/MSource/Point.hpp @@ -128,12 +128,13 @@ void TPoint::setup(void) template void TPoint::execute(void) { - std::vector position = strToVec(par().position); - SitePropagator id; - LOG(Message) << "Creating point source at position [" << par().position - << "]" << std::endl; - PropagatorField &src = envGet(PropagatorField, getName()); + << "]" << std::endl; + + std::vector position = strToVec(par().position); + auto &src = envGet(PropagatorField, getName()); + SitePropagator id; + id = 1.; src = zero; pokeSite(id, src, position); diff --git a/extras/Hadrons/Modules/MSource/SeqConserved.hpp b/extras/Hadrons/Modules/MSource/SeqConserved.hpp index e8f91be1..9ccbee1b 100644 --- a/extras/Hadrons/Modules/MSource/SeqConserved.hpp +++ b/extras/Hadrons/Modules/MSource/SeqConserved.hpp @@ -82,6 +82,7 @@ public: virtual ~TSeqConserved(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -110,6 +111,14 @@ std::vector TSeqConserved::getInput(void) return in; } +template +std::vector TSeqConserved::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TSeqConserved::getOutput(void) { @@ -123,7 +132,7 @@ template void TSeqConserved::setup(void) { auto Ls_ = env().getObjectLs(par().action); - env().template registerLattice(getName(), Ls_); + envCreateLat(PropagatorField, getName(), Ls_); } // execution /////////////////////////////////////////////////////////////////// @@ -143,9 +152,9 @@ void TSeqConserved::execute(void) << par().mu << ") for " << par().tA << " <= t <= " << par().tB << std::endl; } - PropagatorField &src = *env().template createLattice(getName()); - PropagatorField &q = *env().template getObject(par().q); - FMat &mat = *(env().template getObject(par().action)); + auto &src = envGet(PropagatorField, getName()); + auto &q = envGet(PropagatorField, par().q); + auto &mat = envGet(FMat, par().action); std::vector mom = strToVec(par().mom); mat.SeqConservedCurrent(q, src, par().curr_type, par().mu, diff --git a/extras/Hadrons/Modules/MSource/SeqGamma.hpp b/extras/Hadrons/Modules/MSource/SeqGamma.hpp index 8f67f8fa..d2b3c958 100644 --- a/extras/Hadrons/Modules/MSource/SeqGamma.hpp +++ b/extras/Hadrons/Modules/MSource/SeqGamma.hpp @@ -80,12 +80,16 @@ public: virtual ~TSeqGamma(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup virtual void setup(void); // execution virtual void execute(void); +private: + bool hasPhase_{false}; + std::string momphName_, tName_; }; MODULE_REGISTER_NS(SeqGamma, TSeqGamma, MSource); @@ -97,6 +101,8 @@ MODULE_REGISTER_NS(SeqGamma, TSeqGamma, MSource); template TSeqGamma::TSeqGamma(const std::string name) : Module(name) +, momphName_ (name + "_momph") +, tName_ (name + "_t") {} // dependencies/products /////////////////////////////////////////////////////// @@ -108,6 +114,14 @@ std::vector TSeqGamma::getInput(void) return in; } +template +std::vector TSeqGamma::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TSeqGamma::getOutput(void) { @@ -120,7 +134,10 @@ std::vector TSeqGamma::getOutput(void) template void TSeqGamma::setup(void) { - env().template registerLattice(getName()); + envCreateLat(PropagatorField, getName()); + envCacheLat(Lattice>, tName_); + envCacheLat(LatticeComplex, momphName_); + envTmpLat(LatticeComplex, "coor"); } // execution /////////////////////////////////////////////////////////////////// @@ -138,23 +155,29 @@ void TSeqGamma::execute(void) << " sequential source for " << par().tA << " <= t <= " << par().tB << std::endl; } - PropagatorField &src = *env().template createLattice(getName()); - PropagatorField &q = *env().template getObject(par().q); - Lattice> t(env().getGrid()); - LatticeComplex ph(env().getGrid()), coor(env().getGrid()); - Gamma g(par().gamma); - std::vector p; - Complex i(0.0,1.0); + auto &src = envGet(PropagatorField, getName()); + auto &q = envGet(PropagatorField, par().q); + auto &ph = envGet(LatticeComplex, momphName_); + auto &t = envGet(Lattice>, tName_); + Gamma g(par().gamma); - p = strToVec(par().mom); - ph = zero; - for(unsigned int mu = 0; mu < env().getNd(); mu++) + if (!hasPhase_) { - LatticeCoordinate(coor, mu); - ph = ph + p[mu]*coor*((1./(env().getGrid()->_fdimensions[mu]))); + Complex i(0.0,1.0); + std::vector p; + + envGetTmp(LatticeComplex, coor); + p = strToVec(par().mom); + ph = zero; + for(unsigned int mu = 0; mu < env().getNd(); mu++) + { + LatticeCoordinate(coor, mu); + ph = ph + (p[mu]/env().getGrid()->_fdimensions[mu])*coor; + } + ph = exp((Real)(2*M_PI)*i*ph); + LatticeCoordinate(t, Tp); + hasPhase_ = true; } - ph = exp((Real)(2*M_PI)*i*ph); - LatticeCoordinate(t, Tp); src = where((t >= par().tA) and (t <= par().tB), ph*(g*q), 0.*q); } diff --git a/extras/Hadrons/Modules/MSource/Wall.hpp b/extras/Hadrons/Modules/MSource/Wall.hpp index 57dee06d..d9814d9e 100644 --- a/extras/Hadrons/Modules/MSource/Wall.hpp +++ b/extras/Hadrons/Modules/MSource/Wall.hpp @@ -72,12 +72,16 @@ public: virtual ~TWall(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup virtual void setup(void); // execution virtual void execute(void); +private: + bool hasPhase_{false}; + std::string momphName_, tName_; }; MODULE_REGISTER_NS(Wall, TWall, MSource); @@ -89,17 +93,27 @@ MODULE_REGISTER_NS(Wall, TWall, MSource); template TWall::TWall(const std::string name) : Module(name) +, momphName_ (name + "_momph") +, tName_ (name + "_t") {} // dependencies/products /////////////////////////////////////////////////////// template std::vector TWall::getInput(void) { - std::vector in; + std::vector in = {}; return in; } +template +std::vector TWall::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TWall::getOutput(void) { @@ -112,7 +126,7 @@ std::vector TWall::getOutput(void) template void TWall::setup(void) { - env().template registerLattice(getName()); + envCreateLat(PropagatorField, getName()); } // execution /////////////////////////////////////////////////////////////////// @@ -122,21 +136,28 @@ void TWall::execute(void) LOG(Message) << "Generating wall source at t = " << par().tW << " with momentum " << par().mom << std::endl; - PropagatorField &src = *env().template createLattice(getName()); - Lattice> t(env().getGrid()); - LatticeComplex ph(env().getGrid()), coor(env().getGrid()); - std::vector p; - Complex i(0.0,1.0); + auto &src = envGet(PropagatorField, getName()); + auto &ph = envGet(LatticeComplex, momphName_); + auto &t = envGet(Lattice>, tName_); - p = strToVec(par().mom); - ph = zero; - for(unsigned int mu = 0; mu < Nd; mu++) + if (!hasPhase_) { - LatticeCoordinate(coor, mu); - ph = ph + p[mu]*coor*((1./(env().getGrid()->_fdimensions[mu]))); + Complex i(0.0,1.0); + std::vector p; + + envGetTmp(LatticeComplex, coor); + p = strToVec(par().mom); + ph = zero; + for(unsigned int mu = 0; mu < env().getNd(); mu++) + { + LatticeCoordinate(coor, mu); + ph = ph + (p[mu]/env().getGrid()->_fdimensions[mu])*coor; + } + ph = exp((Real)(2*M_PI)*i*ph); + LatticeCoordinate(t, Tp); + hasPhase_ = true; } - ph = exp((Real)(2*M_PI)*i*ph); - LatticeCoordinate(t, Tp); + src = 1.; src = where((t == par().tW), src*ph, 0.*src); } diff --git a/extras/Hadrons/Modules/MSource/Z2.hpp b/extras/Hadrons/Modules/MSource/Z2.hpp index e2cc4f34..2e864ff0 100644 --- a/extras/Hadrons/Modules/MSource/Z2.hpp +++ b/extras/Hadrons/Modules/MSource/Z2.hpp @@ -75,12 +75,16 @@ public: virtual ~TZ2(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup virtual void setup(void); // execution virtual void execute(void); +private: + bool hasT_{false}; + std::string tName_; }; MODULE_REGISTER_NS(Z2, TZ2, MSource); @@ -93,6 +97,7 @@ MODULE_REGISTER_NS(ScalarZ2, TZ2, MSource); template TZ2::TZ2(const std::string name) : Module(name) +, tName_ (name + "_t") {} // dependencies/products /////////////////////////////////////////////////////// @@ -104,6 +109,14 @@ std::vector TZ2::getInput(void) return in; } +template +std::vector TZ2::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TZ2::getOutput(void) { @@ -116,29 +129,36 @@ std::vector TZ2::getOutput(void) template void TZ2::setup(void) { - env().template registerLattice(getName()); + envCreateLat(PropagatorField, getName()); + envCacheLat(Lattice>, tName_); + envTmpLat(LatticeComplex, "eta"); } // execution /////////////////////////////////////////////////////////////////// template void TZ2::execute(void) { - Lattice> t(env().getGrid()); - LatticeComplex eta(env().getGrid()); - Complex shift(1., 1.); - if (par().tA == par().tB) { LOG(Message) << "Generating Z_2 wall source at t= " << par().tA - << std::endl; + << std::endl; } else { LOG(Message) << "Generating Z_2 band for " << par().tA << " <= t <= " - << par().tB << std::endl; + << par().tB << std::endl; } - PropagatorField &src = *env().template createLattice(getName()); - LatticeCoordinate(t, Tp); + + auto &src = envGet(PropagatorField, getName()); + auto &t = envGet(Lattice>, getName()); + Complex shift(1., 1.); + + if (!hasT_) + { + LatticeCoordinate(t, Tp); + hasT_ = true; + } + envGetTmp(LatticeComplex, eta); bernoulli(*env().get4dRng(), eta); eta = (2.*eta - shift)*(1./::sqrt(2.)); eta = where((t >= par().tA) and (t <= par().tB), eta, 0.*eta); diff --git a/extras/Hadrons/Modules/MUtilities/TestSeqConserved.hpp b/extras/Hadrons/Modules/MUtilities/TestSeqConserved.hpp index b0f2846f..081d2911 100644 --- a/extras/Hadrons/Modules/MUtilities/TestSeqConserved.hpp +++ b/extras/Hadrons/Modules/MUtilities/TestSeqConserved.hpp @@ -79,6 +79,7 @@ public: virtual ~TTestSeqConserved(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -107,6 +108,14 @@ std::vector TTestSeqConserved::getInput(void) return in; } +template +std::vector TTestSeqConserved::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TTestSeqConserved::getOutput(void) { @@ -124,36 +133,37 @@ void TTestSeqConserved::setup(void) { HADRON_ERROR(Size, "Ls mismatch between quark action and propagator"); } + envTmpLat(PropagatorField, "tmp"); + envTmpLat(LatticeComplex, "c"); } // execution /////////////////////////////////////////////////////////////////// template void TTestSeqConserved::execute(void) { - PropagatorField tmp(env().getGrid()); - PropagatorField &q = *env().template getObject(par().q); - PropagatorField &qSeq = *env().template getObject(par().qSeq); - FMat &act = *(env().template getObject(par().action)); - Gamma g5(Gamma::Algebra::Gamma5); - Gamma::Algebra gA = (par().curr == Current::Axial) ? - Gamma::Algebra::Gamma5 : - Gamma::Algebra::Identity; - Gamma g(gA); - SitePropagator qSite; - Complex test_S, test_V, check_S, check_V; - std::vector check_buf; - LatticeComplex c(env().getGrid()); - // Check sequential insertion of current gives same result as conserved // current sink upon contraction. Assume q uses a point source. - std::vector siteCoord; + + auto &q = envGet(PropagatorField, par().q); + auto &qSeq = envGet(PropagatorField, par().qSeq); + auto &act = envGet(FMat, par().action); + Gamma g5(Gamma::Algebra::Gamma5); + Gamma::Algebra gA = (par().curr == Current::Axial) ? + Gamma::Algebra::Gamma5 : + Gamma::Algebra::Identity; + Gamma g(gA); + SitePropagator qSite; + Complex test_S, test_V, check_S, check_V; + std::vector check_buf; + std::vector siteCoord; + + envGetTmp(PropagatorField, tmp); + envGetTmp(LatticeComplex, c); siteCoord = strToVec(par().origin); peekSite(qSite, qSeq, siteCoord); test_S = trace(qSite*g); test_V = trace(qSite*g*Gamma::gmu[par().mu]); - act.ContractConservedCurrent(q, q, tmp, par().curr, par().mu); - c = trace(tmp*g); sliceSum(c, check_buf, Tp); check_S = TensorRemove(check_buf[par().t_J]); diff --git a/extras/Hadrons/Modules/MUtilities/TestSeqGamma.hpp b/extras/Hadrons/Modules/MUtilities/TestSeqGamma.hpp index 9736ab54..30bd4b69 100644 --- a/extras/Hadrons/Modules/MUtilities/TestSeqGamma.hpp +++ b/extras/Hadrons/Modules/MUtilities/TestSeqGamma.hpp @@ -63,6 +63,7 @@ public: virtual ~TTestSeqGamma(void) = default; // dependency relation virtual std::vector getInput(void); + virtual std::vector getReference(void); virtual std::vector getOutput(void); protected: // setup @@ -91,6 +92,14 @@ std::vector TTestSeqGamma::getInput(void) return in; } +template +std::vector TTestSeqGamma::getReference(void) +{ + std::vector ref = {}; + + return ref; +} + template std::vector TTestSeqGamma::getOutput(void) { @@ -103,26 +112,27 @@ std::vector TTestSeqGamma::getOutput(void) template void TTestSeqGamma::setup(void) { - + envTmpLat(LatticeComplex, "c"); } // execution /////////////////////////////////////////////////////////////////// template void TTestSeqGamma::execute(void) { - PropagatorField &q = *env().template getObject(par().q); - PropagatorField &qSeq = *env().template getObject(par().qSeq); - LatticeComplex c(env().getGrid()); - Gamma g5(Gamma::Algebra::Gamma5); - Gamma g(par().gamma); - SitePropagator qSite; - Complex test, check; + auto &q = envGet(PropagatorField, par().q); + auto &qSeq = envGet(PropagatorField, par().qSeq); + Gamma g5(Gamma::Algebra::Gamma5); + Gamma g(par().gamma); + SitePropagator qSite; + Complex test, check; std::vector check_buf; + std::vector siteCoord; // Check sequential insertion of gamma matrix gives same result as // insertion of gamma at sink upon contraction. Assume q uses a point // source. - std::vector siteCoord; + + envGetTmp(LatticeComplex, c); siteCoord = strToVec(par().origin); peekSite(qSite, qSeq, siteCoord); test = trace(g*qSite); diff --git a/extras/Hadrons/modules.inc b/extras/Hadrons/modules.inc index 5ce2435f..2f4d183e 100644 --- a/extras/Hadrons/modules.inc +++ b/extras/Hadrons/modules.inc @@ -1,45 +1,38 @@ modules_cc =\ - Modules/MGauge/Unit.cc - # Modules/MContraction/WeakHamiltonianEye.cc \ - # Modules/MContraction/WeakHamiltonianNonEye.cc \ - # Modules/MContraction/WeakNeutral4ptDisc.cc \ - # Modules/MGauge/Load.cc \ - # Modules/MGauge/Random.cc \ - # Modules/MGauge/StochEm.cc \ - # Modules/MScalar/ChargedProp.cc \ - # Modules/MScalar/FreeProp.cc + Modules/MGauge/Unit.cc \ + Modules/MContraction/WeakHamiltonianEye.cc \ + Modules/MContraction/WeakHamiltonianNonEye.cc \ + Modules/MContraction/WeakNeutral4ptDisc.cc \ + Modules/MGauge/Load.cc \ + Modules/MGauge/Random.cc \ + Modules/MGauge/StochEm.cc modules_hpp =\ Modules/MAction/DWF.hpp \ Modules/MAction/Wilson.hpp \ Modules/MSink/Point.hpp \ Modules/MSource/Point.hpp \ + Modules/MGauge/Load.hpp \ + Modules/MGauge/Random.hpp \ + Modules/MGauge/StochEm.hpp \ Modules/MGauge/Unit.hpp \ Modules/MSolver/RBPrecCG.hpp \ Modules/MFermion/GaugeProp.hpp \ - Modules/MContraction/Meson.hpp - - # Modules/MContraction/Baryon.hpp \ - # Modules/MContraction/DiscLoop.hpp \ - # Modules/MContraction/Gamma3pt.hpp \ - # Modules/MContraction/WardIdentity.hpp \ - # Modules/MContraction/WeakHamiltonian.hpp \ - # Modules/MContraction/WeakHamiltonianEye.hpp \ - # Modules/MContraction/WeakHamiltonianNonEye.hpp \ - # Modules/MContraction/WeakNeutral4ptDisc.hpp \ - # Modules/MFermion/GaugeProp.hpp \ - # Modules/MGauge/Load.hpp \ - # Modules/MGauge/Random.hpp \ - # Modules/MGauge/StochEm.hpp \ - # Modules/MLoop/NoiseLoop.hpp \ - # Modules/MScalar/ChargedProp.hpp \ - # Modules/MScalar/FreeProp.hpp \ - # Modules/MScalar/Scalar.hpp \ - # Modules/MSink/Smear.hpp \ - # Modules/MSolver/RBPrecCG.hpp \ - # Modules/MSource/SeqConserved.hpp \ - # Modules/MSource/SeqGamma.hpp \ - # Modules/MSource/Wall.hpp \ - # Modules/MSource/Z2.hpp \ - # Modules/MUtilities/TestSeqConserved.hpp \ - # Modules/MUtilities/TestSeqGamma.hpp + Modules/MContraction/Baryon.hpp \ + Modules/MContraction/DiscLoop.hpp \ + Modules/MContraction/Gamma3pt.hpp \ + Modules/MContraction/Meson.hpp \ + Modules/MContraction/WardIdentity.hpp \ + Modules/MContraction/WeakHamiltonian.hpp \ + Modules/MContraction/WeakHamiltonianEye.hpp \ + Modules/MContraction/WeakHamiltonianNonEye.hpp \ + Modules/MContraction/WeakNeutral4ptDisc.hpp \ + Modules/MLoop/NoiseLoop.hpp \ + Modules/MSink/Smear.hpp \ + Modules/MSolver/RBPrecCG.hpp \ + Modules/MSource/SeqConserved.hpp \ + Modules/MSource/SeqGamma.hpp \ + Modules/MSource/Wall.hpp \ + Modules/MSource/Z2.hpp \ + Modules/MUtilities/TestSeqConserved.hpp \ + Modules/MUtilities/TestSeqGamma.hpp diff --git a/tests/hadrons/Test_hadrons.hpp b/tests/hadrons/Test_hadrons.hpp index 9bd3ee0a..0265f5a6 100644 --- a/tests/hadrons/Test_hadrons.hpp +++ b/tests/hadrons/Test_hadrons.hpp @@ -118,7 +118,7 @@ inline void makeWilsonAction(Application &application, std::string actionName, std::string &gaugeField, double mass, std::string boundary = "1 1 1 -1") { - if (!(Environment::getInstance().hasModule(actionName))) + if (!(VirtualMachine::getInstance().hasModule(actionName))) { MAction::Wilson::Par actionPar; actionPar.gauge = gaugeField; @@ -144,7 +144,7 @@ inline void makeDWFAction(Application &application, std::string actionName, std::string &gaugeField, double mass, double M5, unsigned int Ls, std::string boundary = "1 1 1 -1") { - if (!(Environment::getInstance().hasModule(actionName))) + if (!(VirtualMachine::getInstance().hasModule(actionName))) { MAction::DWF::Par actionPar; actionPar.gauge = gaugeField; @@ -173,7 +173,7 @@ inline void makeDWFAction(Application &application, std::string actionName, inline void makeRBPrecCGSolver(Application &application, std::string &solverName, std::string &actionName, double residual = 1e-8) { - if (!(Environment::getInstance().hasModule(solverName))) + if (!(VirtualMachine::getInstance().hasModule(solverName))) { MSolver::RBPrecCG::Par solverPar; solverPar.action = actionName; @@ -195,7 +195,7 @@ inline void makePointSource(Application &application, std::string srcName, std::string pos) { // If the source already exists, don't make the module again. - if (!(Environment::getInstance().hasModule(srcName))) + if (!(VirtualMachine::getInstance().hasModule(srcName))) { MSource::Point::Par pointPar; pointPar.position = pos; @@ -219,7 +219,7 @@ inline void makeSequentialSource(Application &application, std::string srcName, std::string mom = ZERO_MOM) { // If the source already exists, don't make the module again. - if (!(Environment::getInstance().hasModule(srcName))) + if (!(VirtualMachine::getInstance().hasModule(srcName))) { MSource::SeqGamma::Par seqPar; seqPar.q = qSrc; @@ -255,7 +255,7 @@ inline void makeConservedSequentialSource(Application &application, std::string mom = ZERO_MOM) { // If the source already exists, don't make the module again. - if (!(Environment::getInstance().hasModule(srcName))) + if (!(VirtualMachine::getInstance().hasModule(srcName))) { MSource::SeqConserved::Par seqPar; seqPar.q = qSrc; @@ -280,7 +280,7 @@ inline void makeConservedSequentialSource(Application &application, inline void makeNoiseSource(Application &application, std::string &srcName, unsigned int tA, unsigned int tB) { - if (!(Environment::getInstance().hasModule(srcName))) + if (!(VirtualMachine::getInstance().hasModule(srcName))) { MSource::Z2::Par noisePar; noisePar.tA = tA; @@ -302,7 +302,7 @@ inline void makeWallSource(Application &application, std::string &srcName, unsigned int tW, std::string mom = ZERO_MOM) { // If the source already exists, don't make the module again. - if (!(Environment::getInstance().hasModule(srcName))) + if (!(VirtualMachine::getInstance().hasModule(srcName))) { MSource::Wall::Par wallPar; wallPar.tW = tW; @@ -324,7 +324,7 @@ inline void makePointSink(Application &application, std::string &sinkFnct, std::string mom = ZERO_MOM) { // If the sink function already exists, don't make it again. - if (!(Environment::getInstance().hasModule(sinkFnct))) + if (!(VirtualMachine::getInstance().hasModule(sinkFnct))) { MSink::Point::Par pointPar; pointPar.mom = mom; @@ -345,7 +345,7 @@ inline void sinkSmear(Application &application, std::string &sinkFnct, std::string &propName, std::string &smearedProp) { // If the propagator has already been smeared, don't smear it again. - if (!(Environment::getInstance().hasModule(smearedProp))) + if (!(VirtualMachine::getInstance().hasModule(smearedProp))) { MSink::Smear::Par smearPar; smearPar.q = propName; @@ -367,7 +367,7 @@ inline void makePropagator(Application &application, std::string &propName, std::string &srcName, std::string &solver) { // If the propagator already exists, don't make the module again. - if (!(Environment::getInstance().hasModule(propName))) + if (!(VirtualMachine::getInstance().hasModule(propName))) { MFermion::GaugeProp::Par quarkPar; quarkPar.source = srcName; @@ -390,7 +390,7 @@ inline void makeLoop(Application &application, std::string &propName, std::string &srcName, std::string &resName) { // If the loop propagator already exists, don't make the module again. - if (!(Environment::getInstance().hasModule(propName))) + if (!(VirtualMachine::getInstance().hasModule(propName))) { MLoop::NoiseLoop::Par loopPar; loopPar.q = resName; @@ -421,7 +421,7 @@ inline void mesonContraction(Application &application, std::string &sink, std::string gammas = "") { - if (!(Environment::getInstance().hasModule(modName))) + if (!(VirtualMachine::getInstance().hasModule(modName))) { MContraction::Meson::Par mesPar; mesPar.output = output; @@ -453,7 +453,7 @@ inline void gamma3ptContraction(Application &application, unsigned int npt, Gamma::Algebra gamma = Gamma::Algebra::Identity) { std::string modName = std::to_string(npt) + "pt_" + label; - if (!(Environment::getInstance().hasModule(modName))) + if (!(VirtualMachine::getInstance().hasModule(modName))) { MContraction::Gamma3pt::Par gamma3ptPar; gamma3ptPar.output = std::to_string(npt) + "pt/" + label; @@ -487,7 +487,7 @@ inline void weakContraction##top(Application &application, unsigned int npt,\ std::string &label, unsigned int tSnk = 0)\ {\ std::string modName = std::to_string(npt) + "pt_" + label;\ - if (!(Environment::getInstance().hasModule(modName)))\ + if (!(VirtualMachine::getInstance().hasModule(modName)))\ {\ MContraction::WeakHamiltonian##top::Par weakPar;\ weakPar.output = std::to_string(npt) + "pt/" + label;\ @@ -521,7 +521,7 @@ inline void disc0Contraction(Application &application, std::string &label) { std::string modName = "4pt_" + label; - if (!(Environment::getInstance().hasModule(modName))) + if (!(VirtualMachine::getInstance().hasModule(modName))) { MContraction::WeakNeutral4ptDisc::Par disc0Par; disc0Par.output = "4pt/" + label; @@ -547,7 +547,7 @@ inline void discLoopContraction(Application &application, std::string &q_loop, std::string &modName, Gamma::Algebra gamma = Gamma::Algebra::Identity) { - if (!(Environment::getInstance().hasModule(modName))) + if (!(VirtualMachine::getInstance().hasModule(modName))) { MContraction::DiscLoop::Par discPar; discPar.output = "disc/" + modName; @@ -574,7 +574,7 @@ inline void makeWITest(Application &application, std::string &modName, std::string &propName, std::string &actionName, double mass, unsigned int Ls = 1, bool test_axial = false) { - if (!(Environment::getInstance().hasModule(modName))) + if (!(VirtualMachine::getInstance().hasModule(modName))) { MContraction::WardIdentity::Par wiPar; if (Ls > 1) @@ -613,7 +613,7 @@ inline void makeSeqCurrComparison(Application &application, std::string &modName std::string &actionName, std::string &origin, unsigned int t_J, unsigned int mu, Current curr) { - if (!(Environment::getInstance().hasModule(modName))) + if (!(VirtualMachine::getInstance().hasModule(modName))) { MUtilities::TestSeqConserved::Par seqPar; seqPar.q = propName; @@ -646,7 +646,7 @@ inline void makeSeqGamComparison(Application &application, std::string &modName, std::string &origin, Gamma::Algebra gamma, unsigned int t_g) { - if (!(Environment::getInstance().hasModule(modName))) + if (!(VirtualMachine::getInstance().hasModule(modName))) { MUtilities::TestSeqGamma::Par seqPar; seqPar.q = propName;