From 2b598294c99dc1197d4f51556d925f181fec7ade Mon Sep 17 00:00:00 2001 From: Felix Erben Date: Thu, 18 Apr 2019 17:47:09 +0100 Subject: [PATCH] added distil source module --- Hadrons/Modules.hpp | 1 + Hadrons/Modules/MDistil/DistilSource.cc | 7 + Hadrons/Modules/MDistil/DistilSource.hpp | 221 +++++++++++++++++++++++ Hadrons/Modules/MDistil/PerambLight.hpp | 43 +---- Hadrons/modules.inc | 2 + 5 files changed, 236 insertions(+), 38 deletions(-) create mode 100644 Hadrons/Modules/MDistil/DistilSource.cc create mode 100644 Hadrons/Modules/MDistil/DistilSource.hpp diff --git a/Hadrons/Modules.hpp b/Hadrons/Modules.hpp index 94cfd316..dda60c3b 100644 --- a/Hadrons/Modules.hpp +++ b/Hadrons/Modules.hpp @@ -44,6 +44,7 @@ #include #include #include +#include #include #include #include diff --git a/Hadrons/Modules/MDistil/DistilSource.cc b/Hadrons/Modules/MDistil/DistilSource.cc new file mode 100644 index 00000000..9073185a --- /dev/null +++ b/Hadrons/Modules/MDistil/DistilSource.cc @@ -0,0 +1,7 @@ +#include + +using namespace Grid; +using namespace Hadrons; +using namespace MDistil; + +template class Grid::Hadrons::MDistil::TDistilSource; diff --git a/Hadrons/Modules/MDistil/DistilSource.hpp b/Hadrons/Modules/MDistil/DistilSource.hpp new file mode 100644 index 00000000..c9ce63df --- /dev/null +++ b/Hadrons/Modules/MDistil/DistilSource.hpp @@ -0,0 +1,221 @@ +#ifndef Hadrons_MDistil_DistilSource_hpp_ +#define Hadrons_MDistil_DistilSource_hpp_ + +#include +#include +#include +#include +#include +#include +#include + +// These are members of Distillation +#include + +BEGIN_HADRONS_NAMESPACE + + +/****************************************************************************** + * DistilSource * + ******************************************************************************/ +BEGIN_MODULE_NAMESPACE(MDistil) + +class DistilSourcePar: Serializable +{ +public: + GRID_SERIALIZABLE_CLASS_MEMBERS(DistilSourcePar, + std::string, noise, + std::string, eigenPack, + bool, multiFile, + int, tsrc, + int, nnoise, + int, LI, + int, SI, + int, TI, + int, nvec, + int, Ns, + int, Nt, + int, Nt_inv); +}; + +template +class TDistilSource: public Module +{ +public: + FERM_TYPE_ALIASES(FImpl,); +public: + // constructor + TDistilSource(const std::string name); + // destructor + virtual ~TDistilSource(void); + // dependency relation + virtual std::vector getInput(void); + virtual std::vector getOutput(void); + // setup + virtual void setup(void); + // execution + virtual void execute(void); +protected: + // These variables are created in setup() and freed in Cleanup() + GridCartesian * grid3d; // Owned by me, so I must delete it + GridCartesian * grid4d; // Owned by environment (so I won't delete it) +protected: + virtual void Cleanup(void); +}; + +MODULE_REGISTER_TMP(DistilSource, TDistilSource, MDistil); + +/****************************************************************************** + * TDistilSource implementation * + ******************************************************************************/ +// constructor ///////////////////////////////////////////////////////////////// +template +TDistilSource::TDistilSource(const std::string name) +: grid3d{nullptr}, grid4d{nullptr}, Module(name) +{} +// destructor +template +TDistilSource::~TDistilSource(void) +{ + Cleanup(); +}; + +// dependencies/products /////////////////////////////////////////////////////// +template +std::vector TDistilSource::getInput(void) +{ + std::vector in; + + in.push_back(par().noise); + in.push_back(par().eigenPack); + + return in; +} + +template +std::vector TDistilSource::getOutput(void) +{ + std::vector out = {getName()}; + + return out; +} + +// setup /////////////////////////////////////////////////////////////////////// +template +void TDistilSource::setup(void) +{ + Cleanup(); + auto &noise = envGet(std::vector, par().noise); + + int nnoise=par().nnoise; + int LI=par().LI; + int Ns=par().Ns; + int SI=par().SI; + int Nt_inv=par().Nt_inv; + + envCreate(std::vector, getName(), 1, + nnoise*LI*SI*Nt_inv, envGetGrid(FermionField)); + + grid4d = env().getGrid(); + std::vector latt_size = GridDefaultLatt(); + std::vector simd_layout = GridDefaultSimd(Nd, vComplex::Nsimd()); + std::vector mpi_layout = GridDefaultMpi(); + std::vector simd_layout_3 = GridDefaultSimd(Nd-1, vComplex::Nsimd()); + latt_size[Nd-1] = 1; + simd_layout_3.push_back( 1 ); + mpi_layout[Nd-1] = 1; + grid3d = MakeLowerDimGrid(grid4d); + + + envTmp(LatticeSpinColourVector, "tmp2",1,LatticeSpinColourVector(grid4d)); + envTmp(LatticeColourVector, "tmp_nospin",1,LatticeColourVector(grid4d)); + envTmp(LatticeSpinColourVector, "tmp3d",1,LatticeSpinColourVector(grid3d)); + envTmp(LatticeColourVector, "tmp3d_nospin",1,LatticeColourVector(grid3d)); + envTmp(LatticeSpinColourVector, "sink_tslice",1,LatticeSpinColourVector(grid3d)); + envTmp(LatticeColourVector, "evec3d",1,LatticeColourVector(grid3d)); +} + +// clean up any temporaries created by setup (that aren't stored in the environment) +template +void TDistilSource::Cleanup(void) +{ + if( grid3d != nullptr ) { + delete grid3d; + grid3d = nullptr; + } + grid4d = nullptr; +} + +// execution /////////////////////////////////////////////////////////////////// +template +void TDistilSource::execute(void) +{ + + auto &noise = envGet(std::vector, par().noise); + auto &epack = envGet(Grid::Hadrons::EigenPack, par().eigenPack); + auto &rho = envGet(std::vector, getName()); + + envGetTmp(LatticeSpinColourVector, tmp2); + envGetTmp(LatticeColourVector, tmp_nospin); + envGetTmp(LatticeSpinColourVector, tmp3d); + envGetTmp(LatticeColourVector, tmp3d_nospin); + envGetTmp(LatticeSpinColourVector, sink_tslice); + envGetTmp(LatticeColourVector, evec3d); + + + int Ntlocal = grid4d->LocalDimensions()[3]; + int Ntfirst = grid4d->LocalStarts()[3]; + + int tsrc=par().tsrc; + int nnoise=par().nnoise; + int LI=par().LI; + int Ns=par().Ns; + int Nt_inv=par().Nt_inv; // TODO: No input, but define through Nt, TI + int Nt=par().Nt; + int TI=par().TI; + int nvec=par().nvec; + int SI=par().SI; + + bool full_tdil=(TI==Nt); + + int vecindex; + int t_inv; + for (int inoise = 0; inoise < nnoise; inoise++) { + for (int dk = 0; dk < LI; dk++) { + for (int dt = 0; dt < Nt_inv; dt++) { + for (int ds = 0; ds < SI; ds++) { + vecindex = inoise + nnoise * dk + nnoise * LI * ds + nnoise *LI * SI*dt; + rho[vecindex] = zero; + tmp3d_nospin = zero; + for (int it = dt; it < Nt; it += TI){ + if (full_tdil) t_inv = tsrc; else t_inv = it; + if( t_inv >= Ntfirst && t_inv < Ntfirst + Ntlocal ) { + for (int ik = dk; ik < nvec; ik += LI){ + for (int is = ds; is < Ns; is += SI){ + ExtractSliceLocal(evec3d,epack.evec[ik],0,t_inv,3); + tmp3d_nospin = evec3d * noise[inoise + nnoise*(t_inv + Nt*(ik+nvec*is))]; + tmp3d=zero; + pokeSpin(tmp3d,tmp3d_nospin,is); + tmp2=zero; + InsertSliceLocal(tmp3d,tmp2,0,t_inv-Ntfirst,Grid::QCD::Tdir); + rho[vecindex] += tmp2; + } + } + } + } + } + } + } + } + + + // TEST TO SEE WHETHER THIS MIGHT BE THE MEMORY LEAK + Cleanup(); + +} + +END_MODULE_NAMESPACE + +END_HADRONS_NAMESPACE + +#endif // Hadrons_MDistil_DistilSource_hpp_ diff --git a/Hadrons/Modules/MDistil/PerambLight.hpp b/Hadrons/Modules/MDistil/PerambLight.hpp index 8dfff057..2b18e00c 100644 --- a/Hadrons/Modules/MDistil/PerambLight.hpp +++ b/Hadrons/Modules/MDistil/PerambLight.hpp @@ -25,6 +25,7 @@ class PerambLightPar: Serializable public: GRID_SERIALIZABLE_CLASS_MEMBERS(PerambLightPar, std::string, eigenPack, + std::string, noise, std::string, PerambFileName, //stem!!! std::string, UniqueIdentifier, bool, multiFile, @@ -86,6 +87,7 @@ std::vector TPerambLight::getInput(void) in.push_back(par().eigenPack); in.push_back(par().solver); + in.push_back(par().noise); return in; } @@ -93,7 +95,7 @@ std::vector TPerambLight::getInput(void) template std::vector TPerambLight::getOutput(void) { - std::vector out = {getName() + "_perambulator_light",getName() + "_noise",getName() + "_unsmeared_sink"}; + std::vector out = {getName(),getName() + "_unsmeared_sink"}; return out; } @@ -112,10 +114,8 @@ void TPerambLight::setup(void) const int Ns{Distil.Ns}; std::array sIndexNames{"Nt", "nvec", "LI", "nnoise", "Nt_inv", "SI"}; - envCreate(Perambulator, getName() + "_perambulator_light", 1, + envCreate(Perambulator, getName(), 1, sIndexNames,Distil.Nt,nvec,Distil.LI,Distil.nnoise,Distil.Nt_inv,Distil.SI); - envCreate(std::vector, getName() + "_noise", 1, - nvec*Distil.Ns*Distil.Nt*Distil.nnoise); envCreate(std::vector, getName() + "_unsmeared_sink", 1, nnoise*LI*Ns*Nt_inv, envGetGrid(FermionField)); @@ -176,50 +176,17 @@ void TPerambLight::execute(void) const std::string &UniqueIdentifier{par().UniqueIdentifier}; - auto &noise = envGet(std::vector, getName() + "_noise"); + auto &noise = envGet(std::vector, par().noise); auto &perambulator = envGet(Perambulator, getName() + "_perambulator_light"); auto &epack = envGet(Grid::Hadrons::EigenPack, par().eigenPack); auto &unsmeared_sink = envGet(std::vector, getName() + "_unsmeared_sink"); - //Create Noises - GridSerialRNG sRNG; - sRNG.SeedUniqueString(UniqueIdentifier + std::to_string(vm().getTrajectory())); //maybe add more?? - Real rn; - - for (int inoise=0;inoise 0.5) ? -1 : 1; - } - } - } - } - } // Load perambulator if it exists on disk instead of creating it // Not sure this is how we want it - rather specify an input flag 'read' // and assert that the file is there. const std::string &PerambFileName{par().PerambFileName}; -/* if( PerambFileName.length() ){ - bool bExists = false; - { - std::ifstream f(PerambFileName, std::ios::binary); - if( f.is_open() ) - bExists = true; - } - if( bExists ) { - perambulator.ReadBinary(PerambFileName); - return; - } - }*/ envGetTmp(LatticeSpinColourVector, dist_source); envGetTmp(LatticeSpinColourVector, tmp2); diff --git a/Hadrons/modules.inc b/Hadrons/modules.inc index 65e876c9..4cfc5e7e 100644 --- a/Hadrons/modules.inc +++ b/Hadrons/modules.inc @@ -39,6 +39,7 @@ modules_cc =\ Modules/MDistil/g5_multiply.cc \ Modules/MDistil/LapEvec.cc \ Modules/MDistil/Noises.cc \ + Modules/MDistil/DistilSource.cc \ Modules/MDistil/DistilVectors.cc \ Modules/MDistil/BContraction.cc \ Modules/MDistil/Baryon2pt.cc \ @@ -122,6 +123,7 @@ modules_hpp =\ Modules/MDistil/PerambFromSolve.hpp \ Modules/MDistil/Baryon2pt.hpp \ Modules/MDistil/LapEvec.hpp \ + Modules/MDistil/DistilSource.hpp \ Modules/MDistil/BContraction.hpp \ Modules/MDistil/DistilVectors.hpp \ Modules/MDistil/DistilSink.hpp \