Merge branch 'develop' of https://github.com/paboyle/Grid into develop

extras/Hadrons/Modules/MAction/DWF.hpp

@@ -48,7 +48,8 @@ public:
                                     std::string, gauge,
                                     unsigned int, Ls,
                                     double      , mass,
-                                    double      , M5);
+                                    double      , M5,
+                                    std::string , boundary);
 };
 
 template <typename FImpl>
@@ -116,14 +117,19 @@ void TDWF<FImpl>::execute(void)
                  << par().mass << ", M5= " << par().M5 << " and Ls= "
                  << par().Ls << " using gauge field '" << par().gauge << "'"
                  << std::endl;
+    LOG(Message) << "Fermion boundary conditions: " << par().boundary 
+                 << std::endl;
     env().createGrid(par().Ls);
     auto &U      = *env().template getObject<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<Complex> boundary = strToVec<Complex>(par().boundary);
+    typename DomainWallFermion<FImpl>::ImplParams implParams(boundary);
     FMat *fMatPt = new DomainWallFermion<FImpl>(U, g5, grb5, g4, grb4,
-                                                par().mass, par().M5);
+                                                par().mass, par().M5,
+                                                implParams);
     env().setObject(getName(), fMatPt);
 }
 

extras/Hadrons/Modules/MAction/Wilson.hpp

@@ -46,7 +46,8 @@ class WilsonPar: Serializable
 public:
     GRID_SERIALIZABLE_CLASS_MEMBERS(WilsonPar,
                                     std::string, gauge,
-                                    double     , mass);
+                                    double     , mass,
+                                    std::string, boundary);
 };
 
 template <typename FImpl>
@@ -112,10 +113,15 @@ void TWilson<FImpl>::execute()
 {
     LOG(Message) << "Setting up TWilson fermion matrix with m= " << par().mass
                  << " using gauge field '" << par().gauge << "'" << std::endl;
+    LOG(Message) << "Fermion boundary conditions: " << par().boundary 
+                 << std::endl;
     auto &U      = *env().template getObject<LatticeGaugeField>(par().gauge);
     auto &grid   = *env().getGrid();
     auto &gridRb = *env().getRbGrid();
-    FMat *fMatPt = new WilsonFermion<FImpl>(U, grid, gridRb, par().mass);
+    std::vector<Complex> boundary = strToVec<Complex>(par().boundary);
+    typename WilsonFermion<FImpl>::ImplParams implParams(boundary);
+    FMat *fMatPt = new WilsonFermion<FImpl>(U, grid, gridRb, par().mass,
+                                            implParams);
     env().setObject(getName(), fMatPt);
 }
 

extras/Hadrons/Modules/MContraction/Meson.hpp

@@ -131,12 +131,11 @@ std::vector<std::string> TMeson<FImpl1, FImpl2>::getOutput(void)
 template <typename FImpl1, typename FImpl2>
 void TMeson<FImpl1, FImpl2>::parseGammaString(std::vector<GammaPair> &gammaList)
 {
+    gammaList.clear();
     // Determine gamma matrices to insert at source/sink.
     if (par().gammas.compare("all") == 0)
     {
         // Do all contractions.
-        unsigned int n_gam = Ns * Ns;
-        gammaList.resize(n_gam*n_gam);
         for (unsigned int i = 1; i < Gamma::nGamma; i += 2)
         {
             for (unsigned int j = 1; j < Gamma::nGamma; j += 2)

tests/hadrons/Test_hadrons_meson_3pt.cc

@@ -61,6 +61,10 @@ int main(int argc, char *argv[])
     
     // gauge field
     application.createModule<MGauge::Unit>("gauge");
+    
+    // set fermion boundary conditions to be periodic space, antiperiodic time.
+    std::string boundary = "1 1 1 -1";
+
     for (unsigned int i = 0; i < flavour.size(); ++i)
     {
         // actions
@@ -69,6 +73,7 @@ int main(int argc, char *argv[])
         actionPar.Ls    = 12;
         actionPar.M5    = 1.8;
         actionPar.mass  = mass[i];
+        actionPar.boundary = boundary;
         application.createModule<MAction::DWF>("DWF_" + flavour[i], actionPar);
         
         // solvers

tests/hadrons/Test_hadrons_rarekaon.cc

@@ -98,6 +98,10 @@ int main(int argc, char *argv[])
         gaugePar.file = configStem;
         application.createModule<MGauge::Load>("gauge", gaugePar);
     }
+    
+    // set fermion boundary conditions to be periodic space, antiperiodic time.
+    std::string boundary = "1 1 1 -1";
+
     for (unsigned int i = 0; i < flavour.size(); ++i)
     {
         // actions
@@ -106,6 +110,7 @@ int main(int argc, char *argv[])
         actionPar.Ls    = 16;
         actionPar.M5    = 1.8;
         actionPar.mass  = mass[i];
+        actionPar.boundary = boundary;
         application.createModule<MAction::DWF>("DWF_" + flavour[i], actionPar);
 
         // solvers

tests/hadrons/Test_hadrons_spectrum.cc

@@ -63,6 +63,10 @@ int main(int argc, char *argv[])
     MSource::Point::Par ptPar;
     ptPar.position = "0 0 0 0";
     application.createModule<MSource::Point>("pt", ptPar);
+    
+    // set fermion boundary conditions to be periodic space, antiperiodic time.
+    std::string boundary = "1 1 1 -1";
+
     for (unsigned int i = 0; i < flavour.size(); ++i)
     {
         // actions
@@ -71,6 +75,7 @@ int main(int argc, char *argv[])
         actionPar.Ls    = 12;
         actionPar.M5    = 1.8;
         actionPar.mass  = mass[i];
+        actionPar.boundary = boundary;
         application.createModule<MAction::DWF>("DWF_" + flavour[i], actionPar);
         
         // solvers