diff --git a/tests/forces/Test_fthmc.cc b/tests/forces/Test_fthmc.cc index 0a64cdb9..50f60f99 100644 --- a/tests/forces/Test_fthmc.cc +++ b/tests/forces/Test_fthmc.cc @@ -32,9 +32,12 @@ Author: Peter Boyle using namespace std; using namespace Grid; +typedef MobiusFermionD FermionAction; +typedef WilsonImplD FimplD; +typedef WilsonImplD FermionImplPolicy; template -void ForceTest(Action &action,SmearedConfigurationMasked & smU,MomentumFilterBase &Filter) +void ForceTest(Action &action,ConfigurationBase & smU,MomentumFilterBase &Filter) { LatticeGaugeField U = smU.get_U(false); // unsmeared config GridBase *UGrid = U.Grid(); @@ -51,14 +54,14 @@ void ForceTest(Action &action,SmearedConfigurationMasked &action,SmearedConfigurationMasked &action,SmearedConfigurationMasked(UdSdU,mu); Pmu= PeekIndex(P,mu); - dS = dS - trace(Pmu*UdSdUmu)*eps*2.0*2.0; + dS = dS - trace(Pmu*UdSdUmu)*eps*2.0*HMC_MOMENTUM_DENOMINATOR; } ComplexD dSpred = sum(dS); RealD diff = S2-S1-dSpred.real(); @@ -125,8 +128,11 @@ int main (int argc, char ** argv) const int Ls=12; const int Nt = latt_size[3]; GridCartesian * UGrid = SpaceTimeGrid::makeFourDimGrid(GridDefaultLatt(), GridDefaultSimd(Nd,vComplex::Nsimd()),GridDefaultMpi()); + GridRedBlackCartesian * UrbGrid = SpaceTimeGrid::makeFourDimRedBlackGrid(UGrid); GridCartesian * FGrid = SpaceTimeGrid::makeFiveDimGrid(Ls,UGrid); + GridRedBlackCartesian * FrbGrid = SpaceTimeGrid::makeFiveDimRedBlackGrid(Ls,UGrid); + ///////////////////// Gauge Field and Gauge Forces //////////////////////////// LatticeGaugeField U(UGrid); @@ -141,17 +147,40 @@ int main (int argc, char ** argv) #endif - RealD beta=6.0; - WilsonGaugeActionR PlaqAction(beta); - IwasakiGaugeActionR RectAction(beta); + WilsonGaugeActionR PlaqAction(6.0); + IwasakiGaugeActionR RectAction(2.13); + PlaqAction.is_smeared = true; + RectAction.is_smeared = true; + //////////////////////////////////// + // Fermion Action + //////////////////////////////////// + RealD mass=0.01; + RealD pvmass=1.0; + RealD M5=1.8; + RealD b=1.5; + RealD c=0.5; + + // Double versions + std::vector boundary = {1,1,1,-1}; + FermionAction::ImplParams Params(boundary); + FermionAction DdwfPeriodic(U,*FGrid,*FrbGrid,*UGrid,*UrbGrid,mass,M5,b,c,Params); + FermionAction PVPeriodic (U,*FGrid,*FrbGrid,*UGrid,*UrbGrid,pvmass,M5,b,c,Params); + double StoppingCondition = 1.0e-8; + double MaxCGIterations = 50000; + ConjugateGradient CG(StoppingCondition,MaxCGIterations); + + TwoFlavourRatioPseudoFermionAction Nf2(PVPeriodic, DdwfPeriodic,CG,CG); + Nf2.is_smeared = true; + //////////////////////////////////////////////// // Plaquette only FTHMC smearer //////////////////////////////////////////////// double rho = 0.1; Smear_Stout Smearer(rho); - SmearedConfigurationMasked SmartConfig(UGrid,2*Nd,Smearer,true); + SmearedConfigurationMasked SmartConfig(UGrid,2*Nd,Smearer); + SmearedConfiguration StoutConfig(UGrid,1,Smearer); JacobianAction Jacobian(&SmartConfig); @@ -159,12 +188,32 @@ int main (int argc, char ** argv) // Run some tests //////////////////////////////////////////////// MomentumFilterNone FilterNone; + + std::cout << " ********* FIELD TRANSFORM SMEARING ***** "<(PlaqAction,SmartConfig,FilterNone); + SmartConfig.set_Field(U); ForceTest(RectAction,SmartConfig,FilterNone); + SmartConfig.set_Field(U); ForceTest(Jacobian,SmartConfig,FilterNone); + SmartConfig.set_Field(U); + ForceTest(Nf2,SmartConfig,FilterNone); + + std::cout << " ********* STOUT SMEARING ***** "<(PlaqAction,StoutConfig,FilterNone); + + StoutConfig.set_Field(U); + ForceTest(RectAction,StoutConfig,FilterNone); + + StoutConfig.set_Field(U); + ForceTest(Nf2,StoutConfig,FilterNone); + + Grid_finalize(); }