mirror of
https://github.com/paboyle/Grid.git
synced 2024-11-10 07:55:35 +00:00
Simplifying HMC syntax for the final user
This commit is contained in:
parent
510f55ba30
commit
1e9317e5cf
@ -20,10 +20,10 @@
|
|||||||
#define GRID_COMMS_NONE 1
|
#define GRID_COMMS_NONE 1
|
||||||
|
|
||||||
/* GRID_DEFAULT_PRECISION is DOUBLE */
|
/* GRID_DEFAULT_PRECISION is DOUBLE */
|
||||||
/* #undef GRID_DEFAULT_PRECISION_DOUBLE */
|
#define GRID_DEFAULT_PRECISION_DOUBLE 1
|
||||||
|
|
||||||
/* GRID_DEFAULT_PRECISION is SINGLE */
|
/* GRID_DEFAULT_PRECISION is SINGLE */
|
||||||
#define GRID_DEFAULT_PRECISION_SINGLE 1
|
/* #undef GRID_DEFAULT_PRECISION_SINGLE */
|
||||||
|
|
||||||
/* Support Altivec instructions */
|
/* Support Altivec instructions */
|
||||||
/* #undef HAVE_ALTIVEC */
|
/* #undef HAVE_ALTIVEC */
|
||||||
|
@ -17,7 +17,15 @@ namespace Grid{
|
|||||||
namespace QCD{
|
namespace QCD{
|
||||||
|
|
||||||
typedef Action<LatticeLorentzColourMatrix>* ActPtr; // now force the same colours as the rest of the code
|
typedef Action<LatticeLorentzColourMatrix>* ActPtr; // now force the same colours as the rest of the code
|
||||||
typedef std::vector<ActPtr> ActionLevel;
|
struct ActionLevel{
|
||||||
|
int multiplier;
|
||||||
|
public:
|
||||||
|
std::vector<ActPtr> actions;
|
||||||
|
explicit ActionLevel(int mul = 1):multiplier(mul){assert (mul > 0);};
|
||||||
|
void push_back(ActPtr ptr){
|
||||||
|
actions.push_back(ptr);
|
||||||
|
}
|
||||||
|
};
|
||||||
typedef std::vector<ActionLevel> ActionSet;
|
typedef std::vector<ActionLevel> ActionSet;
|
||||||
typedef std::vector<Observer*> ObserverList;
|
typedef std::vector<Observer*> ObserverList;
|
||||||
|
|
||||||
@ -45,7 +53,6 @@ namespace Grid{
|
|||||||
private:
|
private:
|
||||||
IntegratorParameters Params;
|
IntegratorParameters Params;
|
||||||
const ActionSet as;
|
const ActionSet as;
|
||||||
const std::vector<int> Nrel; //relative step size per level
|
|
||||||
std::unique_ptr<LatticeLorentzColourMatrix> P;
|
std::unique_ptr<LatticeLorentzColourMatrix> P;
|
||||||
GridParallelRNG pRNG;
|
GridParallelRNG pRNG;
|
||||||
//ObserverList observers; // not yet
|
//ObserverList observers; // not yet
|
||||||
@ -56,9 +63,9 @@ namespace Grid{
|
|||||||
void notify_observers();
|
void notify_observers();
|
||||||
|
|
||||||
void update_P(LatticeLorentzColourMatrix&U, int level,double ep){
|
void update_P(LatticeLorentzColourMatrix&U, int level,double ep){
|
||||||
for(int a=0; a<as[level].size(); ++a){
|
for(int a=0; a<as[level].actions.size(); ++a){
|
||||||
LatticeLorentzColourMatrix force(U._grid);
|
LatticeLorentzColourMatrix force(U._grid);
|
||||||
as[level].at(a)->deriv(U,force);
|
as[level].actions.at(a)->deriv(U,force);
|
||||||
*P -= force*ep;
|
*P -= force*ep;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@ -84,9 +91,8 @@ namespace Grid{
|
|||||||
Integrator<IntegratorAlgorithm>* Integ);
|
Integrator<IntegratorAlgorithm>* Integ);
|
||||||
public:
|
public:
|
||||||
Integrator(GridBase* grid, IntegratorParameters Par,
|
Integrator(GridBase* grid, IntegratorParameters Par,
|
||||||
ActionSet& Aset, std::vector<int> Nrel_):
|
ActionSet& Aset):
|
||||||
Params(Par),as(Aset),Nrel(Nrel_),P(new LatticeLorentzColourMatrix(grid)),pRNG(grid){
|
Params(Par),as(Aset),P(new LatticeLorentzColourMatrix(grid)),pRNG(grid){
|
||||||
assert(as.size() == Nrel.size());
|
|
||||||
pRNG.SeedRandomDevice();
|
pRNG.SeedRandomDevice();
|
||||||
};
|
};
|
||||||
|
|
||||||
@ -99,8 +105,8 @@ namespace Grid{
|
|||||||
|
|
||||||
MDutils::generate_momenta(*P,pRNG);
|
MDutils::generate_momenta(*P,pRNG);
|
||||||
for(int level=0; level< as.size(); ++level){
|
for(int level=0; level< as.size(); ++level){
|
||||||
for(int actionID=0; actionID<as.at(level).size(); ++actionID){
|
for(int actionID=0; actionID<as[level].actions.size(); ++actionID){
|
||||||
as[level].at(actionID)->init(U, pRNG);
|
as[level].actions.at(actionID)->init(U, pRNG);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@ -123,8 +129,8 @@ namespace Grid{
|
|||||||
|
|
||||||
// Actions
|
// Actions
|
||||||
for(int level=0; level<as.size(); ++level)
|
for(int level=0; level<as.size(); ++level)
|
||||||
for(int actionID=0; actionID<as.at(level).size(); ++actionID)
|
for(int actionID=0; actionID<as[level].actions.size(); ++actionID)
|
||||||
H += as[level].at(actionID)->S(U);
|
H += as[level].actions.at(actionID)->S(U);
|
||||||
|
|
||||||
return H;
|
return H;
|
||||||
}
|
}
|
||||||
|
@ -27,14 +27,14 @@ namespace Grid{
|
|||||||
int fl = Integ->as.size() -1;
|
int fl = Integ->as.size() -1;
|
||||||
double eps = Integ->Params.stepsize;
|
double eps = Integ->Params.stepsize;
|
||||||
|
|
||||||
for(int l=0; l<=level; ++l) eps/= 2.0*Integ->Nrel[l];
|
for(int l=0; l<=level; ++l) eps/= 2.0*Integ->as[l].multiplier;
|
||||||
|
|
||||||
int fin = Integ->Nrel[0];
|
int fin = Integ->as[0].multiplier;
|
||||||
for(int l=1; l<=level; ++l) fin*= 2.0*Integ->Nrel[l];
|
for(int l=1; l<=level; ++l) fin*= 2.0*Integ->as[l].multiplier;
|
||||||
fin = 3*Integ->Params.MDsteps*fin -1;
|
fin = 3*Integ->Params.MDsteps*fin -1;
|
||||||
|
|
||||||
|
|
||||||
for(int e=0; e<Integ->Nrel[level]; ++e){
|
for(int e=0; e<Integ->as[level].multiplier; ++e){
|
||||||
|
|
||||||
if(clock[level] == 0){ // initial half step
|
if(clock[level] == 0){ // initial half step
|
||||||
Integ->update_P(U,level,lambda*eps);
|
Integ->update_P(U,level,lambda*eps);
|
||||||
@ -101,13 +101,13 @@ namespace Grid{
|
|||||||
double eps = Integ->Params.stepsize;
|
double eps = Integ->Params.stepsize;
|
||||||
|
|
||||||
// Get current level step size
|
// Get current level step size
|
||||||
for(int l=0; l<=level; ++l) eps/= Integ->Nrel[l];
|
for(int l=0; l<=level; ++l) eps/= Integ->as[l].multiplier;
|
||||||
|
|
||||||
int fin = 1;
|
int fin = 1;
|
||||||
for(int l=0; l<=level; ++l) fin*= Integ->Nrel[l];
|
for(int l=0; l<=level; ++l) fin*= Integ->as[l].multiplier;
|
||||||
fin = 2*Integ->Params.MDsteps*fin - 1;
|
fin = 2*Integ->Params.MDsteps*fin - 1;
|
||||||
|
|
||||||
for(int e=0; e<Integ->Nrel[level]; ++e){
|
for(int e=0; e<Integ->as[level].multiplier; ++e){
|
||||||
|
|
||||||
if(clock[level] == 0){ // initial half step
|
if(clock[level] == 0){ // initial half step
|
||||||
Integ->update_P(U, level,eps/2.0);
|
Integ->update_P(U, level,eps/2.0);
|
||||||
|
@ -39,10 +39,9 @@ int main (int argc, char ** argv)
|
|||||||
FullSet.push_back(Level1);
|
FullSet.push_back(Level1);
|
||||||
|
|
||||||
// Create integrator
|
// Create integrator
|
||||||
typedef MinimumNorm2 IntegratorAlgorithm;// change here to change the algorithm
|
typedef MinimumNorm2 IntegratorAlgorithm;// change here to modify the algorithm
|
||||||
IntegratorParameters MDpar(12,5,1.0);
|
IntegratorParameters MDpar(12,5,1.0);
|
||||||
std::vector<int> rel ={1};
|
Integrator<IntegratorAlgorithm> MDynamics(&Fine,MDpar, FullSet);
|
||||||
Integrator<IntegratorAlgorithm> MDynamics(&Fine,MDpar, FullSet,rel);
|
|
||||||
|
|
||||||
// Create HMC
|
// Create HMC
|
||||||
HMCparameters HMCpar;
|
HMCparameters HMCpar;
|
||||||
|
Loading…
Reference in New Issue
Block a user