From 9210b0aa6ea3750e3eab7c5811d5f6b33cf19848 Mon Sep 17 00:00:00 2001 From: Nils Asmussen Date: Tue, 20 Aug 2019 15:21:23 +0100 Subject: [PATCH] remove namespace QCD from directory HMC --- HMC/Mobius2p1f.cc | 17 +- HMC/Mobius2p1fEOFA.cc | 373 +++++++++++++++++++-------------------- HMC/Mobius2p1fEOFA_F1.cc | 342 +++++++++++++++++------------------ HMC/Mobius2p1fRHMC.cc | 23 ++- 4 files changed, 376 insertions(+), 379 deletions(-) diff --git a/HMC/Mobius2p1f.cc b/HMC/Mobius2p1f.cc index fe373dcb..5f82e0e7 100644 --- a/HMC/Mobius2p1f.cc +++ b/HMC/Mobius2p1f.cc @@ -31,7 +31,6 @@ directory int main(int argc, char **argv) { using namespace Grid; - using namespace Grid::QCD; Grid_init(&argc, &argv); int threads = GridThread::GetThreads(); @@ -44,18 +43,18 @@ int main(int argc, char **argv) { typedef typename FermionAction::FermionField FermionField; typedef Grid::XmlReader Serialiser; - + //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: IntegratorParameters MD; - // typedef GenericHMCRunner HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; // MD.name = std::string("Leap Frog"); - // typedef GenericHMCRunner HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; // MD.name = std::string("Force Gradient"); - typedef GenericHMCRunner HMCWrapper; + typedef GenericHMCRunner HMCWrapper; MD.name = std::string("MinimumNorm2"); MD.MDsteps = 20; MD.trajL = 1.0; - + HMCparameters HMCparams; HMCparams.StartTrajectory = 0; HMCparams.Trajectories = 200; @@ -67,7 +66,7 @@ int main(int argc, char **argv) { // Grid from the command line arguments --grid and --mpi TheHMC.Resources.AddFourDimGrid("gauge"); // use default simd lanes decomposition - + CheckpointerParameters CPparams; CPparams.config_prefix = "ckpoint_EODWF_lat"; CPparams.rng_prefix = "ckpoint_EODWF_rng"; @@ -81,7 +80,7 @@ int main(int argc, char **argv) { TheHMC.Resources.SetRNGSeeds(RNGpar); // Construct observables - // here there is too much indirection + // here there is too much indirection typedef PlaquetteMod PlaqObs; TheHMC.Resources.AddObservable(); ////////////////////////////////////////////// @@ -118,7 +117,7 @@ int main(int argc, char **argv) { // These lines are unecessary if BC are all periodic std::vector boundary = {1,1,1,-1}; FermionAction::ImplParams Params(boundary); - + double StoppingCondition = 1e-10; double MaxCGIterations = 30000; ConjugateGradient CG(StoppingCondition,MaxCGIterations); diff --git a/HMC/Mobius2p1fEOFA.cc b/HMC/Mobius2p1fEOFA.cc index 35fb19cb..4a37bc22 100644 --- a/HMC/Mobius2p1fEOFA.cc +++ b/HMC/Mobius2p1fEOFA.cc @@ -2,7 +2,7 @@ Grid physics library, www.github.com/paboyle/Grid -Source file: +Source file: Copyright (C) 2015-2016 @@ -34,140 +34,139 @@ directory #define MIXED_PRECISION #endif -namespace Grid{ - namespace QCD{ +NAMESPACE_BEGIN(Grid); - /* - * Need a plan for gauge field update for mixed precision in HMC (2x speed up) - * -- Store the single prec action operator. - * -- Clone the gauge field from the operator function argument. - * -- Build the mixed precision operator dynamically from the passed operator and single prec clone. - */ +/* + * Need a plan for gauge field update for mixed precision in HMC (2x speed up) + * -- Store the single prec action operator. + * -- Clone the gauge field from the operator function argument. + * -- Build the mixed precision operator dynamically from the passed operator and single prec clone. + */ - template - class MixedPrecisionConjugateGradientOperatorFunction : public OperatorFunction { - public: - typedef typename FermionOperatorD::FermionField FieldD; - typedef typename FermionOperatorF::FermionField FieldF; +template +class MixedPrecisionConjugateGradientOperatorFunction : public OperatorFunction { +public: + typedef typename FermionOperatorD::FermionField FieldD; + typedef typename FermionOperatorF::FermionField FieldF; - using OperatorFunction::operator(); + using OperatorFunction::operator(); - RealD Tolerance; - RealD InnerTolerance; //Initial tolerance for inner CG. Defaults to Tolerance but can be changed - Integer MaxInnerIterations; - Integer MaxOuterIterations; - GridBase* SinglePrecGrid4; //Grid for single-precision fields - GridBase* SinglePrecGrid5; //Grid for single-precision fields - RealD OuterLoopNormMult; //Stop the outer loop and move to a final double prec solve when the residual is OuterLoopNormMult * Tolerance + RealD Tolerance; + RealD InnerTolerance; //Initial tolerance for inner CG. Defaults to Tolerance but can be changed + Integer MaxInnerIterations; + Integer MaxOuterIterations; + GridBase* SinglePrecGrid4; //Grid for single-precision fields + GridBase* SinglePrecGrid5; //Grid for single-precision fields + RealD OuterLoopNormMult; //Stop the outer loop and move to a final double prec solve when the residual is OuterLoopNormMult * Tolerance - FermionOperatorF &FermOpF; - FermionOperatorD &FermOpD;; - SchurOperatorF &LinOpF; - SchurOperatorD &LinOpD; + FermionOperatorF &FermOpF; + FermionOperatorD &FermOpD;; + SchurOperatorF &LinOpF; + SchurOperatorD &LinOpD; - Integer TotalInnerIterations; //Number of inner CG iterations - Integer TotalOuterIterations; //Number of restarts - Integer TotalFinalStepIterations; //Number of CG iterations in final patch-up step + Integer TotalInnerIterations; //Number of inner CG iterations + Integer TotalOuterIterations; //Number of restarts + Integer TotalFinalStepIterations; //Number of CG iterations in final patch-up step - MixedPrecisionConjugateGradientOperatorFunction(RealD tol, - Integer maxinnerit, - Integer maxouterit, - GridBase* _sp_grid4, - GridBase* _sp_grid5, - FermionOperatorF &_FermOpF, - FermionOperatorD &_FermOpD, - SchurOperatorF &_LinOpF, - SchurOperatorD &_LinOpD): - LinOpF(_LinOpF), - LinOpD(_LinOpD), - FermOpF(_FermOpF), - FermOpD(_FermOpD), - Tolerance(tol), - InnerTolerance(tol), - MaxInnerIterations(maxinnerit), - MaxOuterIterations(maxouterit), - SinglePrecGrid4(_sp_grid4), - SinglePrecGrid5(_sp_grid5), - OuterLoopNormMult(100.) - { - /* Debugging instances of objects; references are stored - std::cout << GridLogMessage << " Mixed precision CG wrapper LinOpF " < &LinOpU, const FieldD &src, FieldD &psi) { - - std::cout << GridLogMessage << " Mixed precision CG wrapper operator() "<(&LinOpU); - - // std::cout << GridLogMessage << " Mixed precision CG wrapper operator() FermOpU " <_Mat)<_Mat)==&(LinOpD._Mat)); - - //////////////////////////////////////////////////////////////////////////////////// - // Must snarf a single precision copy of the gauge field in Linop_d argument - //////////////////////////////////////////////////////////////////////////////////// - typedef typename FermionOperatorF::GaugeField GaugeFieldF; - typedef typename FermionOperatorF::GaugeLinkField GaugeLinkFieldF; - typedef typename FermionOperatorD::GaugeField GaugeFieldD; - typedef typename FermionOperatorD::GaugeLinkField GaugeLinkFieldD; - - GridBase * GridPtrF = SinglePrecGrid4; - GridBase * GridPtrD = FermOpD.Umu.Grid(); - GaugeFieldF U_f (GridPtrF); - GaugeLinkFieldF Umu_f(GridPtrF); - // std::cout << " Dim gauge field "<Nd()<Nd()<(FermOpD.Umu, mu); - precisionChange(Umu_f,Umu_d); - PokeIndex(FermOpF.Umu, Umu_f, mu); - } - pickCheckerboard(Even,FermOpF.UmuEven,FermOpF.Umu); - pickCheckerboard(Odd ,FermOpF.UmuOdd ,FermOpF.Umu); - - //////////////////////////////////////////////////////////////////////////////////// - // Could test to make sure that LinOpF and LinOpD agree to single prec? - //////////////////////////////////////////////////////////////////////////////////// - /* - GridBase *Fgrid = psi._grid; - FieldD tmp2(Fgrid); - FieldD tmp1(Fgrid); - LinOpU.Op(src,tmp1); - LinOpD.Op(src,tmp2); - std::cout << " Double gauge field "<< norm2(FermOpD.Umu)< MPCG(Tolerance,MaxInnerIterations,MaxOuterIterations,SinglePrecGrid5,LinOpF,LinOpD); - std::cout << GridLogMessage << "Calling mixed precision Conjugate Gradient" < &LinOpU, const FieldD &src, FieldD &psi) { + + std::cout << GridLogMessage << " Mixed precision CG wrapper operator() "<(&LinOpU); + + // std::cout << GridLogMessage << " Mixed precision CG wrapper operator() FermOpU " <_Mat)<_Mat)==&(LinOpD._Mat)); + + //////////////////////////////////////////////////////////////////////////////////// + // Must snarf a single precision copy of the gauge field in Linop_d argument + //////////////////////////////////////////////////////////////////////////////////// + typedef typename FermionOperatorF::GaugeField GaugeFieldF; + typedef typename FermionOperatorF::GaugeLinkField GaugeLinkFieldF; + typedef typename FermionOperatorD::GaugeField GaugeFieldD; + typedef typename FermionOperatorD::GaugeLinkField GaugeLinkFieldD; + + GridBase * GridPtrF = SinglePrecGrid4; + GridBase * GridPtrD = FermOpD.Umu.Grid(); + GaugeFieldF U_f (GridPtrF); + GaugeLinkFieldF Umu_f(GridPtrF); + // std::cout << " Dim gauge field "<Nd()<Nd()<(FermOpD.Umu, mu); + precisionChange(Umu_f,Umu_d); + PokeIndex(FermOpF.Umu, Umu_f, mu); + } + pickCheckerboard(Even,FermOpF.UmuEven,FermOpF.Umu); + pickCheckerboard(Odd ,FermOpF.UmuOdd ,FermOpF.Umu); + + //////////////////////////////////////////////////////////////////////////////////// + // Could test to make sure that LinOpF and LinOpD agree to single prec? + //////////////////////////////////////////////////////////////////////////////////// + /* + GridBase *Fgrid = psi._grid; + FieldD tmp2(Fgrid); + FieldD tmp1(Fgrid); + LinOpU.Op(src,tmp1); + LinOpD.Op(src,tmp2); + std::cout << " Double gauge field "<< norm2(FermOpD.Umu)< MPCG(Tolerance,MaxInnerIterations,MaxOuterIterations,SinglePrecGrid5,LinOpF,LinOpD); + std::cout << GridLogMessage << "Calling mixed precision Conjugate Gradient" < HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; // MD.name = std::string("Leap Frog"); - typedef GenericHMCRunner HMCWrapper; + typedef GenericHMCRunner HMCWrapper; MD.name = std::string("Force Gradient"); - // typedef GenericHMCRunner HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; // MD.name = std::string("MinimumNorm2"); MD.MDsteps = 6; MD.trajL = 1.0; - + HMCparameters HMCparams; HMCparams.StartTrajectory = 590; HMCparams.Trajectories = 1000; @@ -208,7 +207,7 @@ int main(int argc, char **argv) { // Grid from the command line arguments --grid and --mpi TheHMC.Resources.AddFourDimGrid("gauge"); // use default simd lanes decomposition - + CheckpointerParameters CPparams; CPparams.config_prefix = "ckpoint_EODWF_lat"; CPparams.rng_prefix = "ckpoint_EODWF_rng"; @@ -222,7 +221,7 @@ int main(int argc, char **argv) { TheHMC.Resources.SetRNGSeeds(RNGpar); // Construct observables - // here there is too much indirection + // here there is too much indirection typedef PlaquetteMod PlaqObs; TheHMC.Resources.AddObservable(); ////////////////////////////////////////////// @@ -233,7 +232,7 @@ int main(int argc, char **argv) { Real strange_mass = 0.04; Real pv_mass = 1.0; RealD M5 = 1.8; - RealD b = 1.0; + RealD b = 1.0; RealD c = 0.0; std::vector hasenbusch({ 0.1, 0.3, 0.6 }); @@ -262,7 +261,7 @@ int main(int argc, char **argv) { std::vector boundary = {1,1,1,-1}; FermionAction::ImplParams Params(boundary); FermionActionF::ImplParams ParamsF(boundary); - + double ActionStoppingCondition = 1e-10; double DerivativeStoppingCondition = 1e-6; double MaxCGIterations = 30000; @@ -293,7 +292,7 @@ int main(int argc, char **argv) { OFRp.degree = 14; OFRp.precision= 50; - + MobiusEOFAFermionR Strange_Op_L (U , *FGrid , *FrbGrid , *GridPtr , *GridRBPtr , strange_mass, strange_mass, pv_mass, 0.0, -1, M5, b, c); MobiusEOFAFermionF Strange_Op_LF(UF, *FGridF, *FrbGridF, *GridPtrF, *GridRBPtrF, strange_mass, strange_mass, pv_mass, 0.0, -1, M5, b, c); MobiusEOFAFermionR Strange_Op_R (U , *FGrid , *FrbGrid , *GridPtr , *GridRBPtr , pv_mass, strange_mass, pv_mass, -1.0, 1, M5, b, c); @@ -310,50 +309,50 @@ int main(int argc, char **argv) { LinearOperatorEOFAF Strange_LinOp_RF(Strange_Op_RF); MxPCG_EOFA ActionCGL(ActionStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_LF,Strange_Op_L, - Strange_LinOp_LF,Strange_LinOp_L); + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_LF,Strange_Op_L, + Strange_LinOp_LF,Strange_LinOp_L); MxPCG_EOFA DerivativeCGL(DerivativeStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_LF,Strange_Op_L, - Strange_LinOp_LF,Strange_LinOp_L); - - MxPCG_EOFA ActionCGR(ActionStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_RF,Strange_Op_R, - Strange_LinOp_RF,Strange_LinOp_R); - - MxPCG_EOFA DerivativeCGR(DerivativeStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_RF,Strange_Op_R, - Strange_LinOp_RF,Strange_LinOp_R); + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_LF,Strange_Op_L, + Strange_LinOp_LF,Strange_LinOp_L); - ExactOneFlavourRatioPseudoFermionAction - EOFA(Strange_Op_L, Strange_Op_R, - ActionCG, - ActionCGL, ActionCGR, - DerivativeCGL, DerivativeCGR, - OFRp, true); + MxPCG_EOFA ActionCGR(ActionStoppingCondition, + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_RF,Strange_Op_R, + Strange_LinOp_RF,Strange_LinOp_R); + + MxPCG_EOFA DerivativeCGR(DerivativeStoppingCondition, + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_RF,Strange_Op_R, + Strange_LinOp_RF,Strange_LinOp_R); + + ExactOneFlavourRatioPseudoFermionAction + EOFA(Strange_Op_L, Strange_Op_R, + ActionCG, + ActionCGL, ActionCGR, + DerivativeCGL, DerivativeCGR, + OFRp, true); #else - ExactOneFlavourRatioPseudoFermionAction - EOFA(Strange_Op_L, Strange_Op_R, - ActionCG, - ActionCG, ActionCG, - DerivativeCG, DerivativeCG, - OFRp, true); + ExactOneFlavourRatioPseudoFermionAction + EOFA(Strange_Op_L, Strange_Op_R, + ActionCG, + ActionCG, ActionCG, + DerivativeCG, DerivativeCG, + OFRp, true); #endif Level1.push_back(&EOFA); @@ -384,7 +383,7 @@ int main(int argc, char **argv) { std::vector MPCG; std::vector DenominatorsF; std::vector LinOpD; - std::vector LinOpF; + std::vector LinOpF; for(int h=0;h(*Numerators[h],*Denominators[h],*MPCG[h],*ActionMPCG[h],ActionCG)); diff --git a/HMC/Mobius2p1fEOFA_F1.cc b/HMC/Mobius2p1fEOFA_F1.cc index 3d51b16c..9d006da3 100644 --- a/HMC/Mobius2p1fEOFA_F1.cc +++ b/HMC/Mobius2p1fEOFA_F1.cc @@ -2,7 +2,7 @@ Grid physics library, www.github.com/paboyle/Grid -Source file: +Source file: Copyright (C) 2015-2016 @@ -34,123 +34,123 @@ directory #define MIXED_PRECISION #endif -namespace Grid{ - namespace QCD{ +NAMESPACE_BEGIN(Grid); - /* - * Need a plan for gauge field update for mixed precision in HMC (2x speed up) - * -- Store the single prec action operator. - * -- Clone the gauge field from the operator function argument. - * -- Build the mixed precision operator dynamically from the passed operator and single prec clone. - */ +/* + * Need a plan for gauge field update for mixed precision in HMC (2x speed up) + * -- Store the single prec action operator. + * -- Clone the gauge field from the operator function argument. + * -- Build the mixed precision operator dynamically from the passed operator and single prec clone. + */ - template - class MixedPrecisionConjugateGradientOperatorFunction : public OperatorFunction { - public: - typedef typename FermionOperatorD::FermionField FieldD; - typedef typename FermionOperatorF::FermionField FieldF; +template +class MixedPrecisionConjugateGradientOperatorFunction : public OperatorFunction { +public: + typedef typename FermionOperatorD::FermionField FieldD; + typedef typename FermionOperatorF::FermionField FieldF; - using OperatorFunction::operator(); + using OperatorFunction::operator(); - RealD Tolerance; - RealD InnerTolerance; //Initial tolerance for inner CG. Defaults to Tolerance but can be changed - Integer MaxInnerIterations; - Integer MaxOuterIterations; - GridBase* SinglePrecGrid4; //Grid for single-precision fields - GridBase* SinglePrecGrid5; //Grid for single-precision fields - RealD OuterLoopNormMult; //Stop the outer loop and move to a final double prec solve when the residual is OuterLoopNormMult * Tolerance + RealD Tolerance; + RealD InnerTolerance; //Initial tolerance for inner CG. Defaults to Tolerance but can be changed + Integer MaxInnerIterations; + Integer MaxOuterIterations; + GridBase* SinglePrecGrid4; //Grid for single-precision fields + GridBase* SinglePrecGrid5; //Grid for single-precision fields + RealD OuterLoopNormMult; //Stop the outer loop and move to a final double prec solve when the residual is OuterLoopNormMult * Tolerance - FermionOperatorF &FermOpF; - FermionOperatorD &FermOpD;; - SchurOperatorF &LinOpF; - SchurOperatorD &LinOpD; + FermionOperatorF &FermOpF; + FermionOperatorD &FermOpD;; + SchurOperatorF &LinOpF; + SchurOperatorD &LinOpD; - Integer TotalInnerIterations; //Number of inner CG iterations - Integer TotalOuterIterations; //Number of restarts - Integer TotalFinalStepIterations; //Number of CG iterations in final patch-up step + Integer TotalInnerIterations; //Number of inner CG iterations + Integer TotalOuterIterations; //Number of restarts + Integer TotalFinalStepIterations; //Number of CG iterations in final patch-up step - MixedPrecisionConjugateGradientOperatorFunction(RealD tol, - Integer maxinnerit, - Integer maxouterit, - GridBase* _sp_grid4, - GridBase* _sp_grid5, - FermionOperatorF &_FermOpF, - FermionOperatorD &_FermOpD, - SchurOperatorF &_LinOpF, - SchurOperatorD &_LinOpD): - LinOpF(_LinOpF), - LinOpD(_LinOpD), - FermOpF(_FermOpF), - FermOpD(_FermOpD), - Tolerance(tol), - InnerTolerance(tol), - MaxInnerIterations(maxinnerit), - MaxOuterIterations(maxouterit), - SinglePrecGrid4(_sp_grid4), - SinglePrecGrid5(_sp_grid5), - OuterLoopNormMult(100.) - { - /* Debugging instances of objects; references are stored - std::cout << GridLogMessage << " Mixed precision CG wrapper LinOpF " < &LinOpU, const FieldD &src, FieldD &psi) { - - std::cout << GridLogMessage << " Mixed precision CG wrapper operator() "<(&LinOpU); - - // std::cout << GridLogMessage << " Mixed precision CG wrapper operator() FermOpU " <_Mat)<_Mat)==&(LinOpD._Mat)); - - //////////////////////////////////////////////////////////////////////////////////// - // Must snarf a single precision copy of the gauge field in Linop_d argument - //////////////////////////////////////////////////////////////////////////////////// - typedef typename FermionOperatorF::GaugeField GaugeFieldF; - typedef typename FermionOperatorF::GaugeLinkField GaugeLinkFieldF; - typedef typename FermionOperatorD::GaugeField GaugeFieldD; - typedef typename FermionOperatorD::GaugeLinkField GaugeLinkFieldD; - - GridBase * GridPtrF = SinglePrecGrid4; - GridBase * GridPtrD = FermOpD.Umu.Grid(); - GaugeFieldF U_f (GridPtrF); - GaugeLinkFieldF Umu_f(GridPtrF); - // std::cout << " Dim gauge field "<Nd()<Nd()<(FermOpD.Umu, mu); - precisionChange(Umu_f,Umu_d); - PokeIndex(FermOpF.Umu, Umu_f, mu); - } - pickCheckerboard(Even,FermOpF.UmuEven,FermOpF.Umu); - pickCheckerboard(Odd ,FermOpF.UmuOdd ,FermOpF.Umu); - - //////////////////////////////////////////////////////////////////////////////////// - // Make a mixed precision conjugate gradient - //////////////////////////////////////////////////////////////////////////////////// - MixedPrecisionConjugateGradient MPCG(Tolerance,MaxInnerIterations,MaxOuterIterations,SinglePrecGrid5,LinOpF,LinOpD); - std::cout << GridLogMessage << "Calling mixed precision Conjugate Gradient" < &LinOpU, const FieldD &src, FieldD &psi) { + + std::cout << GridLogMessage << " Mixed precision CG wrapper operator() "<(&LinOpU); + + // std::cout << GridLogMessage << " Mixed precision CG wrapper operator() FermOpU " <_Mat)<_Mat)==&(LinOpD._Mat)); + + //////////////////////////////////////////////////////////////////////////////////// + // Must snarf a single precision copy of the gauge field in Linop_d argument + //////////////////////////////////////////////////////////////////////////////////// + typedef typename FermionOperatorF::GaugeField GaugeFieldF; + typedef typename FermionOperatorF::GaugeLinkField GaugeLinkFieldF; + typedef typename FermionOperatorD::GaugeField GaugeFieldD; + typedef typename FermionOperatorD::GaugeLinkField GaugeLinkFieldD; + + GridBase * GridPtrF = SinglePrecGrid4; + GridBase * GridPtrD = FermOpD.Umu.Grid(); + GaugeFieldF U_f (GridPtrF); + GaugeLinkFieldF Umu_f(GridPtrF); + // std::cout << " Dim gauge field "<Nd()<Nd()<(FermOpD.Umu, mu); + precisionChange(Umu_f,Umu_d); + PokeIndex(FermOpF.Umu, Umu_f, mu); + } + pickCheckerboard(Even,FermOpF.UmuEven,FermOpF.Umu); + pickCheckerboard(Odd ,FermOpF.UmuOdd ,FermOpF.Umu); + + //////////////////////////////////////////////////////////////////////////////////// + // Make a mixed precision conjugate gradient + //////////////////////////////////////////////////////////////////////////////////// + MixedPrecisionConjugateGradient MPCG(Tolerance,MaxInnerIterations,MaxOuterIterations,SinglePrecGrid5,LinOpF,LinOpD); + std::cout << GridLogMessage << "Calling mixed precision Conjugate Gradient" < HMCWrapper; - // typedef GenericHMCRunner HMCWrapper; - typedef GenericHMCRunner HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; + typedef GenericHMCRunner HMCWrapper; HMCparameters HMCparams; { @@ -184,7 +184,7 @@ int main(int argc, char **argv) { // Grid from the command line arguments --grid and --mpi TheHMC.Resources.AddFourDimGrid("gauge"); // use default simd lanes decomposition - + CheckpointerParameters CPparams; CPparams.config_prefix = "ckpoint_EODWF_lat"; CPparams.rng_prefix = "ckpoint_EODWF_rng"; @@ -198,7 +198,7 @@ int main(int argc, char **argv) { TheHMC.Resources.SetRNGSeeds(RNGpar); // Construct observables - // here there is too much indirection + // here there is too much indirection typedef PlaquetteMod PlaqObs; TheHMC.Resources.AddObservable(); ////////////////////////////////////////////// @@ -209,7 +209,7 @@ int main(int argc, char **argv) { Real strange_mass = 0.02144; Real pv_mass = 1.0; RealD M5 = 1.8; - RealD b = 1.5; + RealD b = 1.5; RealD c = 0.5; // Copied from paper @@ -222,7 +222,7 @@ int main(int argc, char **argv) { /////////////////////////////////////////////////////////////////////////////////////////////// //Bad choices with large dH. Equalising force L2 norm was not wise. /////////////////////////////////////////////////////////////////////////////////////////////// - //std::vector hasenbusch({ 0.03, 0.2, 0.3, 0.5, 0.8 }); + //std::vector hasenbusch({ 0.03, 0.2, 0.3, 0.5, 0.8 }); //std::vector hasenbusch({ 0.05, 0.2, 0.4, 0.6, 0.8 }); auto GridPtr = TheHMC.Resources.GetCartesian(); @@ -249,7 +249,7 @@ int main(int argc, char **argv) { std::vector boundary = {1,1,1,-1}; FermionAction::ImplParams Params(boundary); FermionActionF::ImplParams ParamsF(boundary); - + double ActionStoppingCondition = 1e-10; double DerivativeStoppingCondition = 1e-7; double MaxCGIterations = 30000; @@ -280,7 +280,7 @@ int main(int argc, char **argv) { OFRp.degree = 12; OFRp.precision= 50; - + MobiusEOFAFermionR Strange_Op_L (U , *FGrid , *FrbGrid , *GridPtr , *GridRBPtr , strange_mass, strange_mass, pv_mass, 0.0, -1, M5, b, c); MobiusEOFAFermionF Strange_Op_LF(UF, *FGridF, *FrbGridF, *GridPtrF, *GridRBPtrF, strange_mass, strange_mass, pv_mass, 0.0, -1, M5, b, c); MobiusEOFAFermionR Strange_Op_R (U , *FGrid , *FrbGrid , *GridPtr , *GridRBPtr , pv_mass, strange_mass, pv_mass, -1.0, 1, M5, b, c); @@ -298,51 +298,51 @@ int main(int argc, char **argv) { LinearOperatorEOFAF Strange_LinOp_RF(Strange_Op_RF); MxPCG_EOFA ActionCGL(ActionStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_LF,Strange_Op_L, - Strange_LinOp_LF,Strange_LinOp_L); + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_LF,Strange_Op_L, + Strange_LinOp_LF,Strange_LinOp_L); MxPCG_EOFA DerivativeCGL(DerivativeStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_LF,Strange_Op_L, - Strange_LinOp_LF,Strange_LinOp_L); - - MxPCG_EOFA ActionCGR(ActionStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_RF,Strange_Op_R, - Strange_LinOp_RF,Strange_LinOp_R); - - MxPCG_EOFA DerivativeCGR(DerivativeStoppingCondition, - MX_inner, - MaxCGIterations, - GridPtrF, - FrbGridF, - Strange_Op_RF,Strange_Op_R, - Strange_LinOp_RF,Strange_LinOp_R); + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_LF,Strange_Op_L, + Strange_LinOp_LF,Strange_LinOp_L); - ExactOneFlavourRatioPseudoFermionAction - EOFA(Strange_Op_L, Strange_Op_R, - ActionCG, - ActionCGL, ActionCGR, - DerivativeCGL, DerivativeCGR, - OFRp, true); + MxPCG_EOFA ActionCGR(ActionStoppingCondition, + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_RF,Strange_Op_R, + Strange_LinOp_RF,Strange_LinOp_R); + + MxPCG_EOFA DerivativeCGR(DerivativeStoppingCondition, + MX_inner, + MaxCGIterations, + GridPtrF, + FrbGridF, + Strange_Op_RF,Strange_Op_R, + Strange_LinOp_RF,Strange_LinOp_R); + + ExactOneFlavourRatioPseudoFermionAction + EOFA(Strange_Op_L, Strange_Op_R, + ActionCG, + ActionCGL, ActionCGR, + DerivativeCGL, DerivativeCGR, + OFRp, true); #else - ExactOneFlavourRatioPseudoFermionAction - EOFA(Strange_Op_L, Strange_Op_R, - ActionCG, - ActionCG, ActionCG, - ActionCG, ActionCG, - // DerivativeCG, DerivativeCG, - OFRp, true); + ExactOneFlavourRatioPseudoFermionAction + EOFA(Strange_Op_L, Strange_Op_R, + ActionCG, + ActionCG, ActionCG, + ActionCG, ActionCG, + // DerivativeCG, DerivativeCG, + OFRp, true); #endif Level1.push_back(&EOFA); @@ -373,7 +373,7 @@ int main(int argc, char **argv) { std::vector MPCG; std::vector DenominatorsF; std::vector LinOpD; - std::vector LinOpF; + std::vector LinOpF; for(int h=0;h(*Numerators[h],*Denominators[h],*MPCG[h],*ActionMPCG[h],ActionCG)); diff --git a/HMC/Mobius2p1fRHMC.cc b/HMC/Mobius2p1fRHMC.cc index 04fb0ee5..82ca4d37 100644 --- a/HMC/Mobius2p1fRHMC.cc +++ b/HMC/Mobius2p1fRHMC.cc @@ -31,7 +31,6 @@ directory int main(int argc, char **argv) { using namespace Grid; - using namespace Grid::QCD; Grid_init(&argc, &argv); int threads = GridThread::GetThreads(); @@ -44,18 +43,18 @@ int main(int argc, char **argv) { typedef typename FermionAction::FermionField FermionField; typedef Grid::XmlReader Serialiser; - + //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: IntegratorParameters MD; - // typedef GenericHMCRunner HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; // MD.name = std::string("Leap Frog"); - // typedef GenericHMCRunner HMCWrapper; + // typedef GenericHMCRunner HMCWrapper; // MD.name = std::string("Force Gradient"); - typedef GenericHMCRunner HMCWrapper; + typedef GenericHMCRunner HMCWrapper; MD.name = std::string("MinimumNorm2"); MD.MDsteps = 20; MD.trajL = 1.0; - + HMCparameters HMCparams; HMCparams.StartTrajectory = 30; HMCparams.Trajectories = 200; @@ -68,7 +67,7 @@ int main(int argc, char **argv) { // Grid from the command line arguments --grid and --mpi TheHMC.Resources.AddFourDimGrid("gauge"); // use default simd lanes decomposition - + CheckpointerParameters CPparams; CPparams.config_prefix = "ckpoint_EODWF_lat"; CPparams.rng_prefix = "ckpoint_EODWF_rng"; @@ -82,7 +81,7 @@ int main(int argc, char **argv) { TheHMC.Resources.SetRNGSeeds(RNGpar); // Construct observables - // here there is too much indirection + // here there is too much indirection typedef PlaquetteMod PlaqObs; TheHMC.Resources.AddObservable(); ////////////////////////////////////////////// @@ -93,11 +92,11 @@ int main(int argc, char **argv) { Real strange_mass = 0.04; Real pv_mass = 1.0; RealD M5 = 1.8; - RealD b = 1.0; + RealD b = 1.0; RealD c = 0.0; - + // FIXME: - // Same in MC and MD + // Same in MC and MD // Need to mix precision too OneFlavourRationalParams OFRp; OFRp.lo = 4.0e-3; @@ -122,7 +121,7 @@ int main(int argc, char **argv) { // These lines are unecessary if BC are all periodic std::vector boundary = {1,1,1,-1}; FermionAction::ImplParams Params(boundary); - + double StoppingCondition = 1e-10; double MaxCGIterations = 30000; ConjugateGradient CG(StoppingCondition,MaxCGIterations);