portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2025-04-24 12:45:56 +01:00
Code Releases Activity

Several updates

Browse Source
This commit is contained in:
Peter Boyle 2022-02-14 17:29:41 +01:00
parent 3355ceea9f
commit 53e1b00cde
4 changed files with 52 additions and 42 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
Grid/qcd/hmc/integrators/Integrator.h
View File

@ -138,9 +138,11 @@ protected:
std::cout << GridLogIntegrator << "Smearing (on/off): " << as[level].actions.at(a)->is_smeared << std::endl; std::cout << GridLogIntegrator << "Smearing (on/off): " << as[level].actions.at(a)->is_smeared << std::endl;
auto name = as[level].actions.at(a)->action_name(); auto name = as[level].actions.at(a)->action_name();
if (as[level].actions.at(a)->is_smeared) Smearer.smeared_force(force); if (as[level].actions.at(a)->is_smeared) Smearer.smeared_force(force);
DumpSliceNorm("force before Ta",force,Nd-1);
force = FieldImplementation::projectForce(force); // Ta for gauge fields force = FieldImplementation::projectForce(force); // Ta for gauge fields
double end_force = usecond(); double end_force = usecond();
DumpSliceNorm("force before filter",force,Nd-1);
MomFilter->applyFilter(force); MomFilter->applyFilter(force);
Real force_abs = std::sqrt(norm2(force)/U.Grid()->gSites()); //average per-site norm. nb. norm2(latt) = \sum_x norm2(latt[x]) Real force_abs = std::sqrt(norm2(force)/U.Grid()->gSites()); //average per-site norm. nb. norm2(latt) = \sum_x norm2(latt[x])
@ -162,20 +164,12 @@ protected:
double time_force = (end_force - start_force) / 1e3; double time_force = (end_force - start_force) / 1e3;
std::cout << GridLogMessage << "["<<level<<"]["<<a<<"] P update elapsed time: " << time_full << " ms (force: " << time_force << " ms)" << std::endl; std::cout << GridLogMessage << "["<<level<<"]["<<a<<"] P update elapsed time: " << time_full << " ms (force: " << time_force << " ms)" << std::endl;
DumpSliceNorm("force",force,Nd-1); DumpSliceNorm("force after filter",force,Nd-1);
/*
auto pol = PeekIndex<LorentzIndex>(force,Nd-1);
DumpSliceNorm("force_t",pol);
pol=Zero();
PokeIndex<LorentzIndex>(force,pol,Nd-1);
DumpSliceNorm("force_xyz",force);
*/
} }
// Force from the other representations // Force from the other representations
as[level].apply(update_P_hireps, Representations, Mom, U, ep); as[level].apply(update_P_hireps, Representations, Mom, U, ep);
MomFilter->applyFilter(Mom);
} }
void update_U(Field& U, double ep) void update_U(Field& U, double ep)
@ -189,8 +183,12 @@ protected:
void update_U(MomentaField& Mom, Field& U, double ep) void update_U(MomentaField& Mom, Field& U, double ep)
{ {
MomentaField MomFiltered(Mom.Grid());
MomFiltered = Mom;
MomFilter->applyFilter(MomFiltered);
// exponential of Mom*U in the gauge fields case // exponential of Mom*U in the gauge fields case
FieldImplementation::update_field(Mom, U, ep); FieldImplementation::update_field(MomFiltered, U, ep);
// Update the smeared fields, can be implemented as observer // Update the smeared fields, can be implemented as observer
Smearer.set_Field(U); Smearer.set_Field(U);
@ -366,7 +364,6 @@ public:
as[level].apply(refresh_hireps, Representations, sRNG, pRNG); as[level].apply(refresh_hireps, Representations, sRNG, pRNG);
} }
MomFilter->applyFilter(P);
} }
// to be used by the actionlevel class to iterate // to be used by the actionlevel class to iterate

10
HMC/Mobius2f.cc
View File

@ -56,12 +56,12 @@ int main(int argc, char **argv) {
MD.trajL = 1.0; MD.trajL = 1.0;
HMCparameters HMCparams; HMCparameters HMCparams;
HMCparams.StartTrajectory = 0; HMCparams.StartTrajectory = 17;
HMCparams.Trajectories = 200; HMCparams.Trajectories = 200;
HMCparams.NoMetropolisUntil= 20; HMCparams.NoMetropolisUntil= 0;
// "[HotStart, ColdStart, TepidStart, CheckpointStart]\n"; // "[HotStart, ColdStart, TepidStart, CheckpointStart]\n";
HMCparams.StartingType =std::string("ColdStart"); // HMCparams.StartingType =std::string("ColdStart");
// HMCparams.StartingType =std::string("CheckpointStart"); HMCparams.StartingType =std::string("CheckpointStart");
HMCparams.MD = MD; HMCparams.MD = MD;
HMCWrapper TheHMC(HMCparams); HMCWrapper TheHMC(HMCparams);
@ -94,7 +94,7 @@ int main(int argc, char **argv) {
RealD b = 1.0; RealD b = 1.0;
RealD c = 0.0; RealD c = 0.0;
std::vector<Real> hasenbusch({ 0.1 }); std::vector<Real> hasenbusch({ 0.1, 0.4, 0.7 });
auto GridPtr = TheHMC.Resources.GetCartesian(); auto GridPtr = TheHMC.Resources.GetCartesian();
auto GridRBPtr = TheHMC.Resources.GetRBCartesian(); auto GridRBPtr = TheHMC.Resources.GetRBCartesian();

33
HMC/Mobius2f_DDHMC_mixed.cc
View File

@ -51,7 +51,7 @@ public:
{}; {};
virtual void refreshRestrict(FermionField &eta) virtual void refreshRestrict(FermionField &eta)
{ {
Domains.ProjectDomain(eta,1); Domains.ProjectDomain(eta,0);
DumpSliceNorm("refresh Restrict eta",eta); DumpSliceNorm("refresh Restrict eta",eta);
}; };
}; };
@ -97,20 +97,29 @@ int main(int argc, char **argv)
// typedef GenericHMCRunner<ForceGradient> HMCWrapper; // typedef GenericHMCRunner<ForceGradient> HMCWrapper;
// MD.name = std::string("Force Gradient"); // MD.name = std::string("Force Gradient");
typedef GenericHMCRunner<MinimumNorm2> HMCWrapper; typedef GenericHMCRunner<MinimumNorm2> HMCWrapper;
/*
MD.name = std::string("MinimumNorm2"); MD.name = std::string("MinimumNorm2");
MD.MDsteps = 4; // dH = 0.08 MD.MDsteps = 4; // dH = 0.08
// MD.MDsteps = 3; // dH = 0.8 // MD.MDsteps = 3; // dH = 0.8
MD.trajL = 1.0; MD.trajL = 1.0;
*/
HMCparameters HMCparams; HMCparameters HMCparams;
HMCparams.StartTrajectory = 48; {
HMCparams.Trajectories = 20; XmlReader HMCrd("HMCparameters.xml");
read(HMCrd,"HMCparameters",HMCparams);
std::cout << GridLogMessage<< HMCparams <<std::endl;
}
HMCWrapper TheHMC(HMCparams);
/*
HMCparams.StartTrajectory = 66;
HMCparams.Trajectories = 200;
HMCparams.NoMetropolisUntil= 0; HMCparams.NoMetropolisUntil= 0;
// "[HotStart, ColdStart, TepidStart, CheckpointStart]\n"; // "[HotStart, ColdStart, TepidStart, CheckpointStart]\n";
// HMCparams.StartingType =std::string("ColdStart"); // HMCparams.StartingType =std::string("ColdStart");
HMCparams.StartingType =std::string("CheckpointStart"); HMCparams.StartingType =std::string("CheckpointStart");
HMCparams.MD = MD; HMCparams.MD = MD;
HMCWrapper TheHMC(HMCparams); */
// Grid from the command line arguments --grid and --mpi // Grid from the command line arguments --grid and --mpi
TheHMC.Resources.AddFourDimGrid("gauge"); // use default simd lanes decomposition TheHMC.Resources.AddFourDimGrid("gauge"); // use default simd lanes decomposition
@ -177,8 +186,10 @@ int main(int argc, char **argv)
double ActionStoppingCondition = 1e-10; double ActionStoppingCondition = 1e-10;
double DerivativeStoppingCondition = 1e-10; double DerivativeStoppingCondition = 1e-10;
// double BoundaryDerivativeStoppingCondition = 1e-6; // double BoundaryDerivativeStoppingCondition = 1e-10; decent acceptance
double BoundaryDerivativeStoppingCondition = 1e-10; double BoundaryDerivativeStoppingCondition = 1e-7; // decent acceptance
// double BoundaryDerivativeStoppingCondition = 1e-6; // bit bigger not huge
// double BoundaryDerivativeStoppingCondition = 1e-5; // Large dH poor acceptance
double MaxCGIterations = 30000; double MaxCGIterations = 30000;
//////////////////////////////////// ////////////////////////////////////
@ -239,8 +250,8 @@ int main(int argc, char **argv)
std::vector<LinearOperatorD *> LinOpD; std::vector<LinearOperatorD *> LinOpD;
std::vector<LinearOperatorF *> LinOpF; std::vector<LinearOperatorF *> LinOpF;
const int MX_inner = 1000; int MX_inner = 1000;
const RealD MX_tol = 1.0e-8; RealD MX_tol = 1.0e-5;
for(int h=0;h<n_hasenbusch+1;h++){ for(int h=0;h<n_hasenbusch+1;h++){
@ -335,7 +346,9 @@ int main(int argc, char **argv)
DirichletHasenD,DirichletHasenF, DirichletHasenD,DirichletHasenF,
Block); Block);
#if 1
std::cout << GridLogMessage << " Boundary NO ratio "<< std::endl; std::cout << GridLogMessage << " Boundary NO ratio "<< std::endl;
MX_tol = 1.0e-5;
Level1.push_back(new Level1.push_back(new
DomainDecomposedBoundaryTwoFlavourPseudoFermion<FimplD,FimplF> DomainDecomposedBoundaryTwoFlavourPseudoFermion<FimplD,FimplF>
(BoundaryDenominator, (BoundaryDenominator,
@ -344,13 +357,13 @@ int main(int argc, char **argv)
DomainDecomposedBoundaryTwoFlavourBosonPseudoFermion<FimplD,FimplF> DomainDecomposedBoundaryTwoFlavourBosonPseudoFermion<FimplD,FimplF>
(BoundaryNumerator, (BoundaryNumerator,
BoundaryDerivativeStoppingCondition,ActionStoppingCondition,MX_tol)); BoundaryDerivativeStoppingCondition,ActionStoppingCondition,MX_tol));
/* #else
Level1.push_back(new Level1.push_back(new
DomainDecomposedBoundaryTwoFlavourRatioPseudoFermion<FimplD,FimplF> DomainDecomposedBoundaryTwoFlavourRatioPseudoFermion<FimplD,FimplF>
(BoundaryNumerator, (BoundaryNumerator,
BoundaryDenominator, BoundaryDenominator,
BoundaryDerivativeStoppingCondition,ActionStoppingCondition)); BoundaryDerivativeStoppingCondition,ActionStoppingCondition));
*/ #endif
///////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////
// Gauge action // Gauge action

30
HMC/Mobius2p1fEOFA.cc
View File

@ -59,20 +59,20 @@ int main(int argc, char **argv) {
IntegratorParameters MD; IntegratorParameters MD;
// typedef GenericHMCRunner<LeapFrog> HMCWrapper; // typedef GenericHMCRunner<LeapFrog> HMCWrapper;
// MD.name = std::string("Leap Frog"); // MD.name = std::string("Leap Frog");
// typedef GenericHMCRunner<ForceGradient> HMCWrapper; typedef GenericHMCRunner<ForceGradient> HMCWrapper;
// MD.name = std::string("Force Gradient"); MD.name = std::string("Force Gradient");
typedef GenericHMCRunner<MinimumNorm2> HMCWrapper; //typedef GenericHMCRunner<MinimumNorm2> HMCWrapper;
MD.name = std::string("MinimumNorm2"); //MD.name = std::string("MinimumNorm2");
MD.MDsteps = 12; MD.MDsteps = 15;
MD.trajL = 1.0; MD.trajL = 1.0;
HMCparameters HMCparams; HMCparameters HMCparams;
HMCparams.StartTrajectory = 211; HMCparams.StartTrajectory = 0;
HMCparams.Trajectories = 1000; HMCparams.Trajectories = 1000;
HMCparams.NoMetropolisUntil= 0; HMCparams.NoMetropolisUntil= 10;
// "[HotStart, ColdStart, TepidStart, CheckpointStart]\n"; // "[HotStart, ColdStart, TepidStart, CheckpointStart]\n";
// HMCparams.StartingType =std::string("ColdStart"); HMCparams.StartingType =std::string("ColdStart");
HMCparams.StartingType =std::string("CheckpointStart"); //HMCparams.StartingType =std::string("CheckpointStart");
HMCparams.MD = MD; HMCparams.MD = MD;
HMCWrapper TheHMC(HMCparams); HMCWrapper TheHMC(HMCparams);
@ -97,16 +97,16 @@ int main(int argc, char **argv) {
TheHMC.Resources.AddObservable<PlaqObs>(); TheHMC.Resources.AddObservable<PlaqObs>();
////////////////////////////////////////////// //////////////////////////////////////////////
const int Ls = 16; const int Ls = 24;
Real beta = 2.13; Real beta = 2.13;
Real light_mass = 0.01; Real light_mass = 0.005;
Real strange_mass = 0.04; Real strange_mass = 0.0362;
Real pv_mass = 1.0; Real pv_mass = 1.0;
RealD M5 = 1.8; RealD M5 = 1.8;
RealD b = 1.0; RealD b = 1.5;
RealD c = 0.0; RealD c = 0.5;
std::vector<Real> hasenbusch({ 0.1, 0.3, 0.6 }); std::vector<Real> hasenbusch({ 0.02, 0.2, 0.6 });
auto GridPtr = TheHMC.Resources.GetCartesian(); auto GridPtr = TheHMC.Resources.GetCartesian();
auto GridRBPtr = TheHMC.Resources.GetRBCartesian(); auto GridRBPtr = TheHMC.Resources.GetRBCartesian();