From 2bd42339198d96e3df8fa938d5e62db0bf125e56 Mon Sep 17 00:00:00 2001 From: Guido Cossu Date: Wed, 7 Dec 2016 04:56:37 +0000 Subject: [PATCH] Completed testing of the HMC for Ls vectorised version (on AVX2) --- lib/cartesian/Cartesian_base.h | 1 + lib/lattice/Lattice_rng.h | 28 +-- lib/qcd/action/fermion/FermionOperatorImpl.h | 5 +- .../Test_hmc_EODWFRatioLsVectorised_Binary.cc | 170 ++++++++++++++++++ tests/hmc/Test_hmc_EODWFRatio_Binary.cc | 27 +-- 5 files changed, 193 insertions(+), 38 deletions(-) create mode 100644 tests/hmc/Test_hmc_EODWFRatioLsVectorised_Binary.cc diff --git a/lib/cartesian/Cartesian_base.h b/lib/cartesian/Cartesian_base.h index 24ff5d4a..71563b2b 100644 --- a/lib/cartesian/Cartesian_base.h +++ b/lib/cartesian/Cartesian_base.h @@ -8,6 +8,7 @@ Author: Peter Boyle Author: paboyle +Author: Guido Cossu This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by diff --git a/lib/lattice/Lattice_rng.h b/lib/lattice/Lattice_rng.h index c5607ed8..2d0e87c7 100644 --- a/lib/lattice/Lattice_rng.h +++ b/lib/lattice/Lattice_rng.h @@ -148,7 +148,7 @@ namespace Grid { ss<<_generators[gen]; ss.seekg(0,ss.beg); for(int i=0;i>saved[i]; + ss>>saved[i]; } } void SetState(std::vector & saved,int gen){ @@ -315,26 +315,26 @@ namespace Grid { for(int gidx=0;gidxGlobalIndexToGlobalCoor(gidx,gcoor); - _grid->GlobalCoorToRankIndex(rank,o_idx,i_idx,gcoor); + int rank,o_idx,i_idx; + _grid->GlobalIndexToGlobalCoor(gidx,gcoor); + _grid->GlobalCoorToRankIndex(rank,o_idx,i_idx,gcoor); - int l_idx=generator_idx(o_idx,i_idx); + int l_idx=generator_idx(o_idx,i_idx); - const int num_rand_seed=16; - std::vector site_seeds(num_rand_seed); - for(int i=0;i site_seeds(num_rand_seed); + for(int i=0;iBroadcast(0,(void *)&site_seeds[0],sizeof(int)*site_seeds.size()); + _grid->Broadcast(0,(void *)&site_seeds[0],sizeof(int)*site_seeds.size()); - if( rank == _grid->ThisRank() ){ + if( rank == _grid->ThisRank() ){ fixedSeed ssrc(site_seeds); typename source::result_type sinit = ssrc(); _generators[l_idx] = RngEngine(sinit); - } + } } _seeded=1; } @@ -373,7 +373,7 @@ namespace Grid { for (int idx = 0; idx < words; idx++) fillScalar(pointer[idx], dist[gdx], _generators[gdx]); } - // merge into SIMD lanes + // merge into SIMD lanes, FIXME suboptimal implementation merge(l._odata[sm], buf); } } diff --git a/lib/qcd/action/fermion/FermionOperatorImpl.h b/lib/qcd/action/fermion/FermionOperatorImpl.h index 78bc52b4..f658872d 100644 --- a/lib/qcd/action/fermion/FermionOperatorImpl.h +++ b/lib/qcd/action/fermion/FermionOperatorImpl.h @@ -304,6 +304,8 @@ class DomainWallVec5dImpl : public PeriodicGaugeImpl< GaugeImplTypes< S,Nrepres // FIXME // Current implementation works, thread safe, probably suboptimal + // Passing through the local coordinate for grid transformation + // the force grid is in general very different from the Ls vectorized grid PARALLEL_FOR_LOOP for (int so = 0; so < grid->oSites(); so++) { @@ -314,9 +316,6 @@ class DomainWallVec5dImpl : public PeriodicGaugeImpl< GaugeImplTypes< S,Nrepres std::vector local_coor; std::vector icoor; grid->iCoorFromIindex(icoor,si); grid->InOutCoorToLocalCoor(ocoor, icoor, local_coor); - //for (int i = 0; i < dimU; i++) local_coor[i] = ocoor[i] + grid->_rdimensions[i]*icoor[i]; - //std::cout << "so: " << so << " si: "<< si << " local_coor: " << local_coor << std::endl; - for (int s = 0; s < LLs; s++) { std::vector slocal_coor(dimF); slocal_coor[0] = s; diff --git a/tests/hmc/Test_hmc_EODWFRatioLsVectorised_Binary.cc b/tests/hmc/Test_hmc_EODWFRatioLsVectorised_Binary.cc new file mode 100644 index 00000000..13542e04 --- /dev/null +++ b/tests/hmc/Test_hmc_EODWFRatioLsVectorised_Binary.cc @@ -0,0 +1,170 @@ +/************************************************************************************* + +Grid physics library, www.github.com/paboyle/Grid + +Source file: ./tests/Test_hmc_WilsonFermionGauge.cc + +Copyright (C) 2015 + +Author: Peter Boyle +Author: neo +Author: Guido Cossu + +This program is free software; you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation; either version 2 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License along +with this program; if not, write to the Free Software Foundation, Inc., +51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + +See the full license in the file "LICENSE" in the top level distribution +directory +*************************************************************************************/ +/* END LEGAL */ +#include + +using namespace std; +using namespace Grid; +using namespace Grid::QCD; + +namespace Grid { +namespace QCD { + + class HMCRunnerParameters : Serializable { + public: + GRID_SERIALIZABLE_CLASS_MEMBERS(HMCRunnerParameters, + double, beta, + double, mass, + int, MaxCGIterations, + double, StoppingCondition, + bool, smearedAction, + int, SaveInterval, + std::string, format, + std::string, conf_prefix, + std::string, rng_prefix, + double, rho, + int, SmearingLevels, + ); + + HMCRunnerParameters() {} + }; + +// Derive from the BinaryHmcRunner (templated for gauge fields) +class HmcRunner : public BinaryHmcRunner { +public: + void BuildTheAction(int argc, char **argv) + + { + typedef DomainWallVec5dImplR ImplPolicy; + typedef ScaledShamirFermion FermionAction; + typedef typename FermionAction::FermionField FermionField; + + const int Ls = 8; + + UGrid = SpaceTimeGrid::makeFourDimGrid(GridDefaultLatt(), GridDefaultSimd(Nd,vComplex::Nsimd()),GridDefaultMpi()); + UrbGrid = SpaceTimeGrid::makeFourDimRedBlackGrid(UGrid); + + + GridCartesian* sUGrid = SpaceTimeGrid::makeFourDimDWFGrid(GridDefaultLatt(),GridDefaultMpi()); + GridRedBlackCartesian* sUrbGrid = SpaceTimeGrid::makeFourDimRedBlackGrid(sUGrid); + + + FGrid = SpaceTimeGrid::makeFiveDimDWFGrid(Ls,UGrid); + FrbGrid = SpaceTimeGrid::makeFiveDimDWFRedBlackGrid(Ls,UGrid); + + // temporarily need a gauge field + LatticeGaugeField U(UGrid); + + // Gauge action + double beta = 4.0; + WilsonGaugeActionR Waction(beta); + + Real mass = 0.04; + Real pv = 1.0; + RealD M5 = 1.5; + RealD scale = 2.0; + FermionAction DenOp(U,*FGrid,*FrbGrid,*sUGrid,*sUrbGrid,mass,M5,scale); + FermionAction NumOp(U,*FGrid,*FrbGrid,*sUGrid,*sUrbGrid,pv,M5,scale); + + double StoppingCondition = 1.0e-8; + double MaxCGIterations = 10000; + ConjugateGradient CG(StoppingCondition,MaxCGIterations); + TwoFlavourEvenOddRatioPseudoFermionAction Nf2(NumOp, DenOp,CG,CG); + + // Set smearing (true/false), default: false + Nf2.is_smeared = true; + + // Collect actions + // here an example of 2 level integration + ActionLevel Level1(1); + Level1.push_back(&Nf2); + + // this level will integrate with a + // step that is 4 times finer + // than the previous level + ActionLevel Level2(4); + Level2.push_back(&Waction); + + TheAction.push_back(Level1); + TheAction.push_back(Level2); + + // Add observables + int SaveInterval = 1; + std::string format = std::string("IEEE64BIG"); + std::string conf_prefix = std::string("DWF_ckpoint_lat"); + std::string rng_prefix = std::string("DWF_ckpoint_rng"); + BinaryHmcCheckpointer Checkpoint(conf_prefix, rng_prefix, SaveInterval, format); + // Can implement also a specific function in the hmcrunner + // AddCheckpoint (...) that takes the same parameters + a string/tag + // defining the type of the checkpointer + // with tags can be implemented by overloading and no ifs + // Then force all checkpoint to have few common functions + // return an object that is then passed to the Run function + + PlaquetteLogger PlaqLog(std::string("Plaquette")); + ObservablesList.push_back(&PlaqLog); + ObservablesList.push_back(&Checkpoint); + + // Smearing section, omit if not needed + double rho = 0.1; // smearing parameter + int Nsmear = 2; // number of smearing levels + Smear_Stout Stout(rho); + SmearedConfiguration SmearingPolicy( + UGrid, Nsmear, Stout); + /////////////////// + + Run(argc, argv, Checkpoint, SmearingPolicy); + //Run(argc, argv, Checkpoint); // no smearing +}; +}; +} +} + +int main(int argc, char **argv) { + Grid_init(&argc, &argv); + + int threads = GridThread::GetThreads(); + std::cout << GridLogMessage << "Grid is setup to use " << threads + << " threads" << std::endl; + + HmcRunner TheHMC; + + + // Seeds for the random number generators + std::vector SerSeed({1, 2, 3, 4, 5}); + std::vector ParSeed({6, 7, 8, 9, 5}); + TheHMC.RNGSeeds(SerSeed, ParSeed); + + TheHMC.MDparameters.set(20, 1.0);// MDsteps, traj length + + TheHMC.BuildTheAction(argc, argv); + + Grid_finalize(); +} diff --git a/tests/hmc/Test_hmc_EODWFRatio_Binary.cc b/tests/hmc/Test_hmc_EODWFRatio_Binary.cc index 65833024..481ff96d 100644 --- a/tests/hmc/Test_hmc_EODWFRatio_Binary.cc +++ b/tests/hmc/Test_hmc_EODWFRatio_Binary.cc @@ -62,8 +62,7 @@ public: void BuildTheAction(int argc, char **argv) { - //typedef WilsonImplR ImplPolicy; - typedef DomainWallVec5dImplR ImplPolicy; + typedef WilsonImplR ImplPolicy; typedef ScaledShamirFermion FermionAction; typedef typename FermionAction::FermionField FermionField; @@ -71,37 +70,23 @@ public: UGrid = SpaceTimeGrid::makeFourDimGrid(GridDefaultLatt(), GridDefaultSimd(Nd,vComplex::Nsimd()),GridDefaultMpi()); UrbGrid = SpaceTimeGrid::makeFourDimRedBlackGrid(UGrid); - - - GridCartesian* sUGrid = SpaceTimeGrid::makeFourDimDWFGrid(GridDefaultLatt(),GridDefaultMpi()); - GridRedBlackCartesian* sUrbGrid = SpaceTimeGrid::makeFourDimRedBlackGrid(sUGrid); - - FGrid = SpaceTimeGrid::makeFiveDimDWFGrid(Ls,UGrid); - FrbGrid = SpaceTimeGrid::makeFiveDimDWFRedBlackGrid(Ls,UGrid); - - /* FGrid = SpaceTimeGrid::makeFiveDimGrid(Ls,UGrid); FrbGrid = SpaceTimeGrid::makeFiveDimRedBlackGrid(Ls,UGrid); - */ - // temporarily need a gauge field LatticeGaugeField U(UGrid); // Gauge action double beta = 4.0; - WilsonGaugeActionR Iaction(beta); + WilsonGaugeActionR Waction(beta); Real mass = 0.04; Real pv = 1.0; RealD M5 = 1.5; RealD scale = 2.0; - FermionAction DenOp(U,*FGrid,*FrbGrid,*sUGrid,*sUrbGrid,mass,M5,scale); - FermionAction NumOp(U,*FGrid,*FrbGrid,*sUGrid,*sUrbGrid,pv,M5,scale); - - std::cout << GridLogMessage << "Frb Osites: " << FrbGrid->oSites() << std::endl; - std::cout << GridLogMessage << "sUGrid Osites: " << sUGrid->oSites() << std::endl; + FermionAction DenOp(U,*FGrid,*FrbGrid,*UGrid,*UrbGrid,mass,M5,scale); + FermionAction NumOp(U,*FGrid,*FrbGrid,*UGrid,*UrbGrid,pv,M5,scale); double StoppingCondition = 1.0e-8; double MaxCGIterations = 10000; @@ -109,7 +94,7 @@ public: TwoFlavourEvenOddRatioPseudoFermionAction Nf2(NumOp, DenOp,CG,CG); // Set smearing (true/false), default: false - Nf2.is_smeared = false; + Nf2.is_smeared = true; // Collect actions // here an example of 2 level integration @@ -120,7 +105,7 @@ public: // step that is 4 times finer // than the previous level ActionLevel Level2(4); - Level2.push_back(&Iaction); + Level2.push_back(&Waction); TheAction.push_back(Level1); TheAction.push_back(Level2);