Switch to std::unique_ptr<GridCartesian> grid3d;

Remove hand-coded reference to pi - switch to <math.h> definition

Hadrons/Modules/MDistil/Distil.hpp

@@ -30,6 +30,9 @@
 #ifndef Hadrons_MDistil_Distil_hpp_
 #define Hadrons_MDistil_Distil_hpp_
 
+#define _USE_MATH_DEFINES
+#include <math.h>
+
 #include <Hadrons/NamedTensor.hpp>
 #include <Hadrons/Module.hpp>
 #include <Hadrons/ModuleFactory.hpp>
@@ -43,7 +46,7 @@ BEGIN_MODULE_NAMESPACE(MDistil)
 /******************************************************************************
  Distillation code that is common across modules
 
- Documentation on how t use this code available at
+ Documentation on how to use this code available at
 
  *  https://aportelli.github.io/Hadrons-doc/#/mdistil  *
  
@@ -119,7 +122,7 @@ inline void RotateEigen(std::vector<LatticeColourVector> & evec)
         //const Grid::Complex phase = std::conj(cplx0 / cplx0_mag);
         const Real argphase = std::arg(phase);
 #endif
-        std::cout << GridLogMessage << "RotateEigen() : Site 0 : |" << cplx0 << "|=" << cplx0_mag << " => phase=" << (argphase / 3.14159265) << " pi" << std::endl;
+        std::cout << GridLogMessage << "RotateEigen() : Site 0 : |" << cplx0 << "|=" << cplx0_mag << " => phase=" << (argphase / M_PI) << " pi" << std::endl;
         {
             // TODO: Only really needed on the master slice
             for( int k = 0 ; k < evec.size() ; k++ )
@@ -137,4 +140,4 @@ inline void RotateEigen(std::vector<LatticeColourVector> & evec)
 
 END_MODULE_NAMESPACE
 END_HADRONS_NAMESPACE
-#endif // Hadrons_MDistil_Distil_hpp_
+#endif

Hadrons/Modules/MDistil/DistilPar.hpp

@@ -33,11 +33,11 @@
 #include <Hadrons/Modules/MDistil/Distil.hpp>
 
 BEGIN_HADRONS_NAMESPACE
+BEGIN_MODULE_NAMESPACE(MDistil)
 
 /******************************************************************************
  *                         DistilPar                                 *
  ******************************************************************************/
-BEGIN_MODULE_NAMESPACE(MDistil)
 
 template <typename FImpl>
 class TDistilPar: public Module<DistilParameters>
@@ -63,9 +63,7 @@ MODULE_REGISTER_TMP(DistilPar, TDistilPar<FIMPL>, MDistil);
  ******************************************************************************/
 // constructor /////////////////////////////////////////////////////////////////
 template <typename FImpl>
-TDistilPar<FImpl>::TDistilPar(const std::string name)
+TDistilPar<FImpl>::TDistilPar(const std::string name) : Module<DistilParameters>(name) {}
-: Module<DistilParameters>(name)
-{}
 
 // dependencies/products ///////////////////////////////////////////////////////
 template <typename FImpl>
@@ -95,7 +93,5 @@ void TDistilPar<FImpl>::execute(void)
 }
 
 END_MODULE_NAMESPACE
-
 END_HADRONS_NAMESPACE
-
+#endif
-#endif // Hadrons_MDistil_DistilPar_hpp_

Hadrons/Modules/MDistil/DistilVectors.hpp

@@ -49,7 +49,7 @@ public:
                                     std::string, lapevec,
                                     std::string, rho,
                                     std::string, phi,
-                                    std::string, DistilPar);
+                                    std::string, DistilParams);
 };
 
 template <typename FImpl>
@@ -69,18 +69,9 @@ public:
     // execution
     virtual void execute(void);
 protected:
-    // These variables are created in setup() and freed in Cleanup()
+    std::unique_ptr<GridCartesian> grid3d; // Owned by me, so I must delete it
-    GridCartesian * grid3d; // Owned by me, so I must delete it
-    GridCartesian * grid4d; // Owned by environment (so I won't delete it)
-    virtual void Cleanup(void);
 public:
     // These variables contain parameters
-    std::string PerambulatorName;
-    std::string NoiseVectorName;
-    std::string LapEvecName;
-    std::string DParName;
-    bool bMakeRho;
-    bool bMakePhi;
     std::string RhoName;
     std::string PhiName;
 };
@@ -92,47 +83,31 @@ MODULE_REGISTER_TMP(DistilVectors, TDistilVectors<FIMPL>, MDistil);
  ******************************************************************************/
 // constructor /////////////////////////////////////////////////////////////////
 template <typename FImpl>
-TDistilVectors<FImpl>::TDistilVectors(const std::string name)
+TDistilVectors<FImpl>::TDistilVectors(const std::string name) : Module<DistilVectorsPar>(name) {}
-:  grid3d{nullptr}, grid4d{nullptr}, Module<DistilVectorsPar>(name)
-{}
-// destructor
-template <typename FImpl>
-TDistilVectors<FImpl>::~TDistilVectors(void)
-{
-    Cleanup();
-};
 
 // dependencies/products ///////////////////////////////////////////////////////
 template <typename FImpl>
 std::vector<std::string> TDistilVectors<FImpl>::getInput(void)
 {
-    PerambulatorName = par().perambulator;
+    return {par().noise,par().perambulator,par().lapevec,par().DistilParams};
-    NoiseVectorName = par().noise;
-    LapEvecName = par().lapevec;
-    DParName = par().DistilPar;
-    return { PerambulatorName, NoiseVectorName, LapEvecName, DParName };
 }
 
 template <typename FImpl>
 std::vector<std::string> TDistilVectors<FImpl>::getOutput(void)
 {
-    RhoName  = par().rho;
+    RhoName = par().rho;
-    PhiName  = par().phi;
+    PhiName = par().phi;
-    bMakeRho = ( RhoName.size() > 0 );
+    if (RhoName.empty() && PhiName.empty())
-    bMakePhi = ( PhiName.size() > 0 );
-    if (!bMakeRho && !bMakePhi)
     {
         RhoName = getName();
         PhiName = RhoName;
         RhoName.append("_rho");
         PhiName.append("_phi");
-        bMakeRho = true;
-        bMakePhi = true;
     }
     std::vector<std::string> out;
-    if (bMakeRho)
+    if (!RhoName.empty())
         out.push_back(RhoName);
-    if (bMakePhi)
+    if (!PhiName.empty())
         out.push_back(PhiName);
     return out;
 }
@@ -141,64 +116,45 @@ std::vector<std::string> TDistilVectors<FImpl>::getOutput(void)
 template <typename FImpl>
 void TDistilVectors<FImpl>::setup(void)
 {
-    Cleanup();
-    auto &noise        = envGet(NoiseTensor,  NoiseVectorName);
-    auto &perambulator = envGet(PerambTensor, PerambulatorName);
-    auto &DPar         = envGet(MDistil::DistilParameters,  DParName);
-    
     // We expect the perambulator to have been created with these indices
-    assert( perambulator.ValidateIndexNames() && "Perambulator index names bad" );
+    auto &perambulator = envGet(PerambTensor, par().perambulator);
+    if (!perambulator.ValidateIndexNames())
+        HADRONS_ERROR(Range,"Perambulator index names bad");
+
+    const DistilParameters &dp{envGet(DistilParameters, par().DistilParams)};
+    const int Nt{env().getDim(Tdir)};
+    const bool full_tdil{ dp.TI == Nt };
+    const int Nt_inv{ full_tdil ? 1 : dp.TI };
     
-    const int Nt{env().getDim(Tdir)}; 
+    if (!RhoName.empty())
-    const int nvec{DPar.nvec}; 
+        envCreate(std::vector<FermionField>, RhoName, 1, dp.nnoise*dp.LI*dp.SI*Nt_inv, envGetGrid(FermionField));
-    const int nnoise{DPar.nnoise}; 
+    if (!PhiName.empty())
-    const int TI{DPar.TI}; 
+        envCreate(std::vector<FermionField>, PhiName, 1, dp.nnoise*dp.LI*dp.SI*Nt_inv, envGetGrid(FermionField));
-    const int LI{DPar.LI}; 
-    const int SI{DPar.SI}; 
-    const bool full_tdil{ TI == Nt }; 
-    const int Nt_inv{ full_tdil ? 1 : TI };
     
-    if (bMakeRho)
-        envCreate(std::vector<FermionField>, RhoName, 1, nnoise*LI*SI*Nt_inv, envGetGrid(FermionField));
-    if (bMakePhi)
-        envCreate(std::vector<FermionField>,   PhiName, 1, nnoise*LI*SI*Nt_inv, envGetGrid(FermionField));
-    
-    grid4d = env().getGrid();
     Coordinate latt_size   = GridDefaultLatt();
     Coordinate mpi_layout  = GridDefaultMpi();
     Coordinate 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);
+    GridCartesian * const grid4d{env().getGrid()};
+    grid3d.reset(MakeLowerDimGrid(grid4d));
     
     envTmp(LatticeSpinColourVector, "source4d",1,LatticeSpinColourVector(grid4d));
-    envTmp(LatticeSpinColourVector, "source3d",1,LatticeSpinColourVector(grid3d));
+    envTmp(LatticeSpinColourVector, "source3d",1,LatticeSpinColourVector(grid3d.get()));
-    envTmp(LatticeColourVector, "source3d_nospin",1,LatticeColourVector(grid3d));
+    envTmp(LatticeColourVector, "source3d_nospin",1,LatticeColourVector(grid3d.get()));
-    envTmp(LatticeSpinColourVector, "sink3d",1,LatticeSpinColourVector(grid3d));
+    envTmp(LatticeSpinColourVector, "sink3d",1,LatticeSpinColourVector(grid3d.get()));
-    envTmp(LatticeColourVector, "evec3d",1,LatticeColourVector(grid3d));
+    envTmp(LatticeColourVector, "evec3d",1,LatticeColourVector(grid3d.get()));
-}
-
-// clean up any temporaries created by setup (that aren't stored in the environment)
-template <typename FImpl>
-void TDistilVectors<FImpl>::Cleanup(void)
-{
-    if ( grid3d != nullptr)
-    {
-        delete grid3d;
-        grid3d = nullptr;
-    }
-    grid4d = nullptr;
 }
 
 // execution ///////////////////////////////////////////////////////////////////
 template <typename FImpl>
 void TDistilVectors<FImpl>::execute(void)
 {
-    auto &noise        = envGet(NoiseTensor, NoiseVectorName);
+    auto &noise        = envGet(NoiseTensor,  par().noise);
-    auto &perambulator = envGet(PerambTensor, PerambulatorName);
+    auto &perambulator = envGet(PerambTensor, par().perambulator);
-    auto &epack        = envGet(Grid::Hadrons::EigenPack<LatticeColourVector>, LapEvecName);
+    auto &epack        = envGet(Grid::Hadrons::EigenPack<LatticeColourVector>, par().lapevec);
-    auto &DPar         = envGet(MDistil::DistilParameters,  DParName);
+    const DistilParameters &DPar{envGet(DistilParameters, par().DistilParams)};
     
     envGetTmp(LatticeSpinColourVector, source4d);
     envGetTmp(LatticeSpinColourVector, source3d);
@@ -206,6 +162,7 @@ void TDistilVectors<FImpl>::execute(void)
     envGetTmp(LatticeSpinColourVector, sink3d);
     envGetTmp(LatticeColourVector,     evec3d);
     
+    GridCartesian * const grid4d{env().getGrid()};
     const int Ntlocal{ grid4d->LocalDimensions()[3] };
     const int Ntfirst{ grid4d->LocalStarts()[3] };
     
@@ -221,7 +178,7 @@ void TDistilVectors<FImpl>::execute(void)
     
     int vecindex;
     int t_inv;
-    if (bMakeRho)
+    if (!RhoName.empty())
     {
         auto &rho = envGet(std::vector<FermionField>, RhoName);
         for (int inoise = 0; inoise < nnoise; inoise++) {
@@ -252,7 +209,7 @@ void TDistilVectors<FImpl>::execute(void)
             }
         }
     }
-    if (bMakePhi) {
+    if (!PhiName.empty()) {
         auto &phi = envGet(std::vector<FermionField>, PhiName);
         for (int inoise = 0; inoise < nnoise; inoise++) {
             for (int dk = 0; dk < LI; dk++) {

Hadrons/Modules/MDistil/LapEvec.hpp

@@ -111,12 +111,7 @@ public:
     // execution
     virtual void execute(void);
 protected:
-    // These variables are created in setup() and freed in Cleanup()
+    std::unique_ptr<GridCartesian> gridLD; // Owned by me, so I must delete it
-    GridCartesian * gridLD; // Owned by me, so I must delete it
-    GridCartesian * gridHD; // Owned by environment (so I won't delete it)
-    std::string sGaugeName;
-protected:
-    virtual void Cleanup(void);
 };
 
 MODULE_REGISTER_TMP(LapEvec, TLapEvec<GIMPL>, MDistil);
@@ -127,61 +122,36 @@ MODULE_REGISTER_TMP(LapEvec, TLapEvec<GIMPL>, MDistil);
 
 // constructor /////////////////////////////////////////////////////////////////
 template <typename GImpl>
-TLapEvec<GImpl>::TLapEvec(const std::string name) : gridLD{nullptr}, Module<LapEvecPar>(name)
+TLapEvec<GImpl>::TLapEvec(const std::string name) : Module<LapEvecPar>(name) {}
-{
-}
-
-// destructor /////////////////////////////////////////////////////////////////
-template <typename GImpl>
-TLapEvec<GImpl>::~TLapEvec()
-{
-    Cleanup();
-}
 
 // dependencies/products ///////////////////////////////////////////////////////
 template <typename GImpl>
 std::vector<std::string> TLapEvec<GImpl>::getInput(void)
 {
-    sGaugeName = par().gauge;
+    return std::vector<std::string>{par().gauge};
-    return std::vector<std::string>{ sGaugeName };
 }
 
 template <typename GImpl>
 std::vector<std::string> TLapEvec<GImpl>::getOutput(void)
 {
-    std::vector<std::string> out = {getName()}; // This is the higher dimensional eigenpack
+    return {getName()}; // This is the higher dimensional eigenpack
-    return out;
 }
 
 // setup ///////////////////////////////////////////////////////////////////////
 template <typename GImpl>
 void TLapEvec<GImpl>::setup(void)
 {
-    Cleanup();
+    GridCartesian * gridHD = env().getGrid();
-    Environment & e{env()};
+    gridLD.reset(MakeLowerDimGrid(gridHD));
-    gridHD = e.getGrid();
-    gridLD = MakeLowerDimGrid( gridHD );
     const int Ntlocal{gridHD->LocalDimensions()[Tdir]};
     // Temporaries
     envTmpLat(GaugeField, "Umu_stout");
     envTmpLat(GaugeField, "Umu_smear");
-    envTmp(LatticeGaugeField, "UmuNoTime",1,LatticeGaugeField(gridLD));
+    envTmp(LatticeGaugeField, "UmuNoTime",1,LatticeGaugeField(gridLD.get()));
-    envTmp(LatticeColourVector, "src",1,LatticeColourVector(gridLD));
+    envTmp(LatticeColourVector, "src",1,LatticeColourVector(gridLD.get()));
     envTmp(std::vector<LapEvecs>, "eig",1,std::vector<LapEvecs>(Ntlocal));
     // Output objects
-    envCreate(LapEvecs, getName(), 1, par().Lanczos.Nvec, gridHD );
+    envCreate(LapEvecs, getName(), 1, par().Lanczos.Nvec, gridHD);
-}
-
-// clean up any temporaries created by setup (that aren't stored in the environment)
-template <typename GImpl>
-void TLapEvec<GImpl>::Cleanup(void)
-{
-    if (gridLD != nullptr)
-    {
-        delete gridLD;
-        gridLD = nullptr;
-    }
-    gridHD = nullptr;
 }
 
 /*************************************************************************************
@@ -260,7 +230,7 @@ void TLapEvec<GImpl>::execute(void)
     
     // Stout smearing
     envGetTmp(GaugeField, Umu_smear);
-    Umu_smear = envGet(GaugeField, sGaugeName); // The smeared field starts off as the Gauge field
+    Umu_smear = envGet(GaugeField, par().gauge); // The smeared field starts off as the Gauge field
     LOG(Message) << "Initial plaquette: " << WilsonLoops<PeriodicGimplR>::avgPlaquette(Umu_smear) << std::endl;
     const StoutParameters &Stout{par().Stout};
     if( Stout.steps )
@@ -282,6 +252,7 @@ void TLapEvec<GImpl>::execute(void)
     envGetTmp(std::vector<LapEvecs>, eig);   // Eigenpack for each timeslice
     envGetTmp(LatticeGaugeField, UmuNoTime); // Gauge field without time dimension
     envGetTmp(LatticeColourVector, src);
+    GridCartesian * gridHD = env().getGrid();
     const int Ntlocal{gridHD->LocalDimensions()[Tdir]};
     const int Ntfirst{gridHD->LocalStarts()[Tdir]};
     uint32_t ConvergenceErrors{0};
@@ -290,7 +261,7 @@ void TLapEvec<GImpl>::execute(void)
         LOG(Message) << "------------------------------------------------------------" << std::endl;
         LOG(Message) << " Compute eigenpack, local timeslice = " << t << " / " << Ntlocal << std::endl;
         LOG(Message) << "------------------------------------------------------------" << std::endl;
-        eig[t].resize(LPar.Nk+LPar.Np,gridLD);
+        eig[t].resize(LPar.Nk+LPar.Np,gridLD.get());
         
         // Construct smearing operator
         ExtractSliceLocal(UmuNoTime,Umu_smear,0,t,Tdir); // switch to 3d/4d objects
@@ -317,7 +288,7 @@ void TLapEvec<GImpl>::execute(void)
             LOG(Error) << "MDistil::LapEvec : Not enough eigenvectors converged. If this occurs in practice, we should modify the eigensolver to iterate once more to ensure the second convergence test does not take us below the requested number of eigenvectors" << std::endl;
         }
         if( Nconv != LPar.Nvec )
-            eig[t].resize( LPar.Nvec, gridLD );
+            eig[t].resize(LPar.Nvec, gridLD.get());
         RotateEigen( eig[t].evec ); // Rotate the eigenvectors into our phase convention
         
         for (int i=0;i<LPar.Nvec;i++){

Hadrons/Modules/MDistil/Noises.hpp

@@ -43,7 +43,7 @@ class NoisesPar: Serializable
 {
 public:
     GRID_SERIALIZABLE_CLASS_MEMBERS(NoisesPar,
-                                    std::string, DistilPar,
+                                    std::string, DistilParams,
                                     std::string, NoiseFileName)
 };
 
@@ -62,8 +62,6 @@ public:
     virtual void setup(void);
     // execution
     virtual void execute(void);
-protected:
-    std::string DParName;
 };
 
 MODULE_REGISTER_TMP(Noises, TNoises<FIMPL>, MDistil);
@@ -73,16 +71,13 @@ MODULE_REGISTER_TMP(Noises, TNoises<FIMPL>, MDistil);
  ******************************************************************************/
 // constructor /////////////////////////////////////////////////////////////////
 template <typename FImpl>
-TNoises<FImpl>::TNoises(const std::string name)
+TNoises<FImpl>::TNoises(const std::string name) : Module<NoisesPar>(name) {}
-: Module<NoisesPar>(name)
-{}
 
 // dependencies/products ///////////////////////////////////////////////////////
 template <typename FImpl>
 std::vector<std::string> TNoises<FImpl>::getInput(void)
 {
-    DParName = par().DistilPar;
+    return {par().DistilParams};
-    return { DParName };
 }
 
 template <typename FImpl>
@@ -96,32 +91,26 @@ std::vector<std::string> TNoises<FImpl>::getOutput(void)
 template <typename FImpl>
 void TNoises<FImpl>::setup(void)
 {
-    auto &DPar         = envGet(MDistil::DistilParameters,  DParName);
+    const DistilParameters &dp{envGet(DistilParameters, par().DistilParams)};
     const int Nt{env().getDim(Tdir)}; 
-    const int nvec{DPar.nvec}; 
+    envCreate(NoiseTensor, getName(), 1, dp.nnoise, Nt, dp.nvec, Ns);
-    const int nnoise{DPar.nnoise}; 
-    envCreate(NoiseTensor, getName(), 1, nnoise, Nt, nvec, Ns);
 }
 
 // execution ///////////////////////////////////////////////////////////////////
 template <typename FImpl>
 void TNoises<FImpl>::execute(void)
 {
-    auto &DPar         = envGet(MDistil::DistilParameters,  DParName);
+    const DistilParameters &dp{envGet(DistilParameters, par().DistilParams)};
     const int Nt{env().getDim(Tdir)};
-    const int nnoise{DPar.nnoise};
+    const bool full_tdil{ dp.TI == Nt };
-    const int nvec{DPar.nvec};
+    const bool exact_distillation{ full_tdil && dp.LI == dp.nvec };
-    const int TI{DPar.TI};
-    const int LI{DPar.LI};
-    const bool full_tdil{ TI == Nt }; 
-    const bool exact_distillation{ full_tdil && LI == nvec }; 
     
     // We use our own seeds so we can specify different noises per quark
     Real rn;
     auto &noise = envGet(NoiseTensor, getName());
-    for (int inoise = 0; inoise < nnoise; inoise++) {
+    for (int inoise = 0; inoise < dp.nnoise; inoise++) {
         for (int t = 0; t < Nt; t++) {
-            for (int ivec = 0; ivec < nvec; ivec++) {
+            for (int ivec = 0; ivec < dp.nvec; ivec++) {
                 for (int is = 0; is < Ns; is++) {
                     if (exact_distillation)
                         noise.tensor(inoise, t, ivec, is) = 1.;
@@ -146,4 +135,4 @@ void TNoises<FImpl>::execute(void)
 
 END_MODULE_NAMESPACE
 END_HADRONS_NAMESPACE
-#endif // Hadrons_MDistil_Noises_hpp_
+#endif

Hadrons/Modules/MDistil/PerambFromSolve.hpp

@@ -51,9 +51,9 @@ public:
                                     std::string, eigenPack,
                                     std::string, PerambFileName,
                                     std::string, solve,
-                                    std::string, nvec_reduced,
+                                    int, nvec_reduced,
-                                    std::string, LI_reduced,
+                                    int, LI_reduced,
-                                    std::string, DistilPar);
+                                    std::string, DistilParams);
 };
 
 template <typename FImpl>
@@ -73,13 +73,7 @@ public:
     // execution
     virtual void execute(void);
 protected:
-    GridCartesian * grid3d; // Owned by me, so I must delete it
+    std::unique_ptr<GridCartesian> grid3d; // Owned by me, so I must delete it
-    GridCartesian * grid4d;
-    //other members
-    std::string DParName;
-protected:
-    virtual void Cleanup(void);
-    
 };
 
 MODULE_REGISTER_TMP(PerambFromSolve, TPerambFromSolve<FIMPL>, MDistil);
@@ -89,67 +83,37 @@ MODULE_REGISTER_TMP(PerambFromSolve, TPerambFromSolve<FIMPL>, MDistil);
  ******************************************************************************/
 // constructor /////////////////////////////////////////////////////////////////
 template <typename FImpl>
-TPerambFromSolve<FImpl>::TPerambFromSolve(const std::string name)
+TPerambFromSolve<FImpl>::TPerambFromSolve(const std::string name) : Module<PerambFromSolvePar>(name){}
-:grid3d{nullptr}, grid4d{nullptr}, Module<PerambFromSolvePar>(name)
-{}
-//destructor
-template <typename FImpl>
-TPerambFromSolve<FImpl>::~TPerambFromSolve(void)
-{
-    Cleanup();
-};
 
 // dependencies/products ///////////////////////////////////////////////////////
 template <typename FImpl>
 std::vector<std::string> TPerambFromSolve<FImpl>::getInput(void)
 {
-    //return std::vector<std::string>{ par().solve, par().eigenPack };
+    return {par().solve, par().eigenPack, par().DistilParams};
-    DParName = par().DistilPar;
-    return {par().solve, par().eigenPack, DParName };
 }
 
 template <typename FImpl>
 std::vector<std::string> TPerambFromSolve<FImpl>::getOutput(void)
 {
-    return std::vector<std::string>{ getName() };
+    return std::vector<std::string>{getName()};
 }
 
 // setup ///////////////////////////////////////////////////////////////////////
 template <typename FImpl>
 void TPerambFromSolve<FImpl>::setup(void)
 {
-    Cleanup();
+    const DistilParameters & dp{envGet(MDistil::DistilParameters, par().DistilParams)};
-    auto &DPar         = envGet(MDistil::DistilParameters,  DParName);
+    const int Nt{env().getDim(Tdir)};
-    const int Nt{env().getDim(Tdir)}; 
+    const bool full_tdil{ dp.TI == Nt };
-    const int nvec{DPar.nvec}; 
+    const int Nt_inv{ full_tdil ? 1 : dp.TI };
-    const int nnoise{DPar.nnoise}; 
+    grid3d.reset( MakeLowerDimGrid( env().getGrid() ) );
-    const int LI{DPar.LI}; 
+    envCreate(PerambTensor, getName(), 1, Nt,par().nvec_reduced,par().LI_reduced,dp.nnoise,Nt_inv,dp.SI);
-    const int TI{DPar.TI}; 
+    envCreate(NoiseTensor, getName() + "_noise", 1, dp.nnoise, Nt, dp.nvec, Ns );
-    const int SI{DPar.SI}; 
+    envTmp(LatticeColourVector, "result3d_nospin",1,LatticeColourVector(grid3d.get()));
-    const bool full_tdil{ TI == Nt }; 
+    envTmp(LatticeColourVector, "evec3d",1,LatticeColourVector(grid3d.get()));
-    const int Nt_inv{ full_tdil ? 1 : TI };
-    const int nvec_reduced{ par().nvec_reduced.empty() ? nvec:std::stoi(par().nvec_reduced)};
-    const int LI_reduced{ par().LI_reduced.empty() ? LI:std::stoi(par().LI_reduced)};
-    grid4d = env().getGrid();
-    grid3d = MakeLowerDimGrid(grid4d);
-    envCreate(PerambTensor, getName(), 1, Nt,nvec_reduced,LI_reduced,nnoise,Nt_inv,SI);
-    envCreate(NoiseTensor, getName() + "_noise", 1, nnoise, Nt, nvec, Ns );
-    envTmp(LatticeColourVector, "result3d_nospin",1,LatticeColourVector(grid3d));
-    envTmp(LatticeColourVector, "evec3d",1,LatticeColourVector(grid3d));
     envTmpLat(LatticeColourVector, "result4d_nospin");
 }
 
-template <typename FImpl>
-void TPerambFromSolve<FImpl>::Cleanup(void)
-{
-    if (grid3d != nullptr)
-    {
-        delete grid3d;
-        grid3d = nullptr;
-    }
-    grid4d = nullptr;
-}
-
 // execution ///////////////////////////////////////////////////////////////////
 template <typename FImpl>
 void TPerambFromSolve<FImpl>::execute(void)
@@ -157,17 +121,12 @@ void TPerambFromSolve<FImpl>::execute(void)
     GridCartesian * grid4d = env().getGrid();
     const int Ntlocal{grid4d->LocalDimensions()[3]};
     const int Ntfirst{grid4d->LocalStarts()[3]};
-    auto &DPar         = envGet(MDistil::DistilParameters,  DParName);
+    const DistilParameters &dp{envGet(DistilParameters, par().DistilParams)};
-    const int Nt{env().getDim(Tdir)}; 
+    const int Nt{env().getDim(Tdir)};
-    const int nvec{DPar.nvec}; 
+    const bool full_tdil{ dp.TI == Nt };
-    const int nnoise{DPar.nnoise}; 
+    const int Nt_inv{ full_tdil ? 1 : dp.TI };
-    const int TI{DPar.TI}; 
+    const int nvec_reduced{par().nvec_reduced};
-    const int LI{DPar.LI}; 
+    const int LI_reduced{par().LI_reduced};
-    const int SI{DPar.SI}; 
-    const bool full_tdil{ TI == Nt }; 
-    const int Nt_inv{ full_tdil ? 1 : TI };
-    const int nvec_reduced{ par().nvec_reduced.empty() ? nvec:std::stoi(par().nvec_reduced)};
-    const int LI_reduced{ par().LI_reduced.empty() ? LI:std::stoi(par().LI_reduced)};
     auto &perambulator  = envGet(PerambTensor, getName());
     auto &solve         = envGet(std::vector<FermionField>, par().solve);
     auto &epack         = envGet(Grid::Hadrons::EigenPack<LatticeColourVector>, par().eigenPack);
@@ -176,17 +135,17 @@ void TPerambFromSolve<FImpl>::execute(void)
     envGetTmp(LatticeColourVector, result3d_nospin);
     envGetTmp(LatticeColourVector, evec3d);
     
-    for (int inoise = 0; inoise < nnoise; inoise++)
+    for (int inoise = 0; inoise < dp.nnoise; inoise++)
     {
         for (int dk = 0; dk < LI_reduced; dk++)
        	{
             for (int dt = 0; dt < Nt_inv; dt++)
 	    {
-                for (int ds = 0; ds < SI; ds++)
+                for (int ds = 0; ds < dp.SI; ds++)
 	       	{
                     for (int is = 0; is < Ns; is++)
 		    {
-                        result4d_nospin = peekSpin(solve[inoise+nnoise*(dk+LI*(dt+Nt_inv*ds))],is);
+                        result4d_nospin = peekSpin(solve[inoise+dp.nnoise*(dk+dp.LI*(dt+Nt_inv*ds))],is);
                         for (int t = Ntfirst; t < Ntfirst + Ntlocal; t++)
 		       	{
                             ExtractSliceLocal(result3d_nospin,result4d_nospin,0,t-Ntfirst,Tdir);

Hadrons/Modules/MDistil/Perambulator.hpp

@@ -48,7 +48,7 @@ public:
                                     std::string, noise,
                                     std::string, PerambFileName,
                                     std::string, UnsmearedSinkFileName,
-                                    std::string, DistilPar);
+                                    std::string, DistilParams);
 };
 
 template <typename FImpl>
@@ -69,12 +69,7 @@ public:
     // execution
     virtual void execute(void);
 protected:
-    virtual void Cleanup(void);
+    std::unique_ptr<GridCartesian> grid3d; // Owned by me, so I must delete it
-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)
-    // Other members
     unsigned int Ls_;
 };
 
@@ -85,22 +80,13 @@ MODULE_REGISTER_TMP(Perambulator, TPerambulator<FIMPL>, MDistil);
  ******************************************************************************/
 // constructor /////////////////////////////////////////////////////////////////
 template <typename FImpl>
-TPerambulator<FImpl>::TPerambulator(const std::string name)
+TPerambulator<FImpl>::TPerambulator(const std::string name) : Module<PerambulatorPar>(name) {}
-: grid3d{nullptr}, grid4d{nullptr}, Module<PerambulatorPar>(name)
-{}
-
-// destructor
-template <typename FImpl>
-TPerambulator<FImpl>::~TPerambulator(void)
-{
-    Cleanup();
-};
 
 // dependencies/products ///////////////////////////////////////////////////////
 template <typename FImpl>
 std::vector<std::string> TPerambulator<FImpl>::getInput(void)
 {
-    return {par().lapevec, par().solver, par().noise, par().DistilPar};
+    return {par().lapevec, par().solver, par().noise, par().DistilParams};
 }
 
 template <typename FImpl>
@@ -113,10 +99,8 @@ std::vector<std::string> TPerambulator<FImpl>::getOutput(void)
 template <typename FImpl>
 void TPerambulator<FImpl>::setup(void)
 {
-    Cleanup();
+    grid3d.reset( MakeLowerDimGrid( env().getGrid() ) );
-    grid4d = env().getGrid();
+    const DistilParameters &dp = envGet(DistilParameters, par().DistilParams);
-    grid3d = MakeLowerDimGrid(grid4d);
-    const DistilParameters &dp = envGet(DistilParameters, par().DistilPar);
     const int  Nt{env().getDim(Tdir)};
     const bool full_tdil{ dp.TI == Nt };
     const int  Nt_inv{ full_tdil ? 1 : dp.TI };
@@ -127,12 +111,12 @@ void TPerambulator<FImpl>::setup(void)
     
     envTmpLat(LatticeSpinColourVector,   "dist_source");
     envTmpLat(LatticeSpinColourVector,   "source4d");
-    envTmp(LatticeSpinColourVector,      "source3d",1,LatticeSpinColourVector(grid3d));
+    envTmp(LatticeSpinColourVector,      "source3d",1,LatticeSpinColourVector(grid3d.get()));
-    envTmp(LatticeColourVector,          "source3d_nospin",1,LatticeColourVector(grid3d));
+    envTmp(LatticeColourVector,          "source3d_nospin",1,LatticeColourVector(grid3d.get()));
     envTmpLat(LatticeSpinColourVector,   "result4d");
     envTmpLat(LatticeColourVector,       "result4d_nospin");
-    envTmp(LatticeColourVector,          "result3d_nospin",1,LatticeColourVector(grid3d));
+    envTmp(LatticeColourVector,          "result3d_nospin",1,LatticeColourVector(grid3d.get()));
-    envTmp(LatticeColourVector,          "evec3d",1,LatticeColourVector(grid3d));
+    envTmp(LatticeColourVector,          "evec3d",1,LatticeColourVector(grid3d.get()));
     
     Ls_ = env().getObjectLs(par().solver);
     envTmpLat(FermionField, "v4dtmp");
@@ -140,23 +124,11 @@ void TPerambulator<FImpl>::setup(void)
     envTmpLat(FermionField, "v5dtmp_sol", Ls_);
 }
 
-// clean up any temporaries created by setup (that aren't stored in the environment)
-template <typename FImpl>
-void TPerambulator<FImpl>::Cleanup(void)
-{
-    if (grid3d != nullptr)
-    {
-        delete grid3d;
-        grid3d = nullptr;
-    }
-    grid4d = nullptr;
-}
-
 // execution ///////////////////////////////////////////////////////////////////
 template <typename FImpl>
 void TPerambulator<FImpl>::execute(void)
 {
-    const DistilParameters &DPar{ envGet(DistilParameters, par().DistilPar) };
+    const DistilParameters &DPar{ envGet(DistilParameters, par().DistilParams) };
     const int Nt{env().getDim(Tdir)};
     const int nvec{DPar.nvec}; 
     const int nnoise{DPar.nnoise}; 
@@ -184,6 +156,7 @@ void TPerambulator<FImpl>::execute(void)
     envGetTmp(LatticeColourVector, result4d_nospin);
     envGetTmp(LatticeColourVector, result3d_nospin);
     envGetTmp(LatticeColourVector, evec3d);
+    GridCartesian * const grid4d{ env().getGrid() }; // Owned by environment (so I won't delete it)
     const int Ntlocal{grid4d->LocalDimensions()[3]};
     const int Ntfirst{grid4d->LocalStarts()[3]};
     const std::string UnsmearedSinkFileName{ par().UnsmearedSinkFileName };
@@ -286,10 +259,8 @@ void TPerambulator<FImpl>::execute(void)
     //Save the unsmeared sinks if filename specified
     if (!UnsmearedSinkFileName.empty())
     {
-       // bool bMulti = ( Hadrons::MDistil::DistilParameters::ParameterDefault( par().UnsmearedSinkMultiFile, 1, false ) != 0 );
-        bool bMulti =0;
         LOG(Message) << "Writing unsmeared sink to " << UnsmearedSinkFileName << std::endl;
-        A2AVectorsIo::write(UnsmearedSinkFileName, unsmeared_sink, bMulti, vm().getTrajectory());
+        A2AVectorsIo::write(UnsmearedSinkFileName, unsmeared_sink, false, vm().getTrajectory());
     }
 }