mirror of
				https://github.com/paboyle/Grid.git
				synced 2025-11-04 14:04:32 +00:00 
			
		
		
		
	Simplifying HMC syntax for the final user
This commit is contained in:
		@@ -20,10 +20,10 @@
 | 
			
		||||
#define GRID_COMMS_NONE 1
 | 
			
		||||
 | 
			
		||||
/* GRID_DEFAULT_PRECISION is DOUBLE */
 | 
			
		||||
/* #undef GRID_DEFAULT_PRECISION_DOUBLE */
 | 
			
		||||
#define GRID_DEFAULT_PRECISION_DOUBLE 1
 | 
			
		||||
 | 
			
		||||
/* GRID_DEFAULT_PRECISION is SINGLE */
 | 
			
		||||
#define GRID_DEFAULT_PRECISION_SINGLE 1
 | 
			
		||||
/* #undef GRID_DEFAULT_PRECISION_SINGLE */
 | 
			
		||||
 | 
			
		||||
/* Support Altivec instructions */
 | 
			
		||||
/* #undef HAVE_ALTIVEC */
 | 
			
		||||
 
 | 
			
		||||
@@ -17,7 +17,15 @@ namespace Grid{
 | 
			
		||||
  namespace QCD{
 | 
			
		||||
 | 
			
		||||
    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<Observer*> ObserverList;
 | 
			
		||||
    
 | 
			
		||||
@@ -45,7 +53,6 @@ namespace Grid{
 | 
			
		||||
    private:
 | 
			
		||||
      IntegratorParameters Params;
 | 
			
		||||
      const ActionSet as;
 | 
			
		||||
      const std::vector<int> Nrel; //relative step size per level
 | 
			
		||||
      std::unique_ptr<LatticeLorentzColourMatrix> P;
 | 
			
		||||
      GridParallelRNG pRNG;
 | 
			
		||||
      //ObserverList observers; // not yet
 | 
			
		||||
@@ -56,9 +63,9 @@ namespace Grid{
 | 
			
		||||
      void notify_observers();
 | 
			
		||||
 | 
			
		||||
      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);
 | 
			
		||||
	  as[level].at(a)->deriv(U,force);
 | 
			
		||||
	  as[level].actions.at(a)->deriv(U,force);
 | 
			
		||||
	  *P -= force*ep;
 | 
			
		||||
	}
 | 
			
		||||
      }
 | 
			
		||||
@@ -84,9 +91,8 @@ namespace Grid{
 | 
			
		||||
					     Integrator<IntegratorAlgorithm>* Integ);
 | 
			
		||||
    public:
 | 
			
		||||
    Integrator(GridBase* grid, IntegratorParameters Par,
 | 
			
		||||
		 ActionSet& Aset, std::vector<int> Nrel_):
 | 
			
		||||
      Params(Par),as(Aset),Nrel(Nrel_),P(new LatticeLorentzColourMatrix(grid)),pRNG(grid){
 | 
			
		||||
	assert(as.size() == Nrel.size());
 | 
			
		||||
		 ActionSet& Aset):
 | 
			
		||||
      Params(Par),as(Aset),P(new LatticeLorentzColourMatrix(grid)),pRNG(grid){
 | 
			
		||||
	pRNG.SeedRandomDevice();
 | 
			
		||||
      };
 | 
			
		||||
      
 | 
			
		||||
@@ -99,8 +105,8 @@ namespace Grid{
 | 
			
		||||
 | 
			
		||||
	MDutils::generate_momenta(*P,pRNG);
 | 
			
		||||
	for(int level=0; level< as.size(); ++level){
 | 
			
		||||
	  for(int actionID=0; actionID<as.at(level).size(); ++actionID){
 | 
			
		||||
	    as[level].at(actionID)->init(U, pRNG);
 | 
			
		||||
	  for(int actionID=0; actionID<as[level].actions.size(); ++actionID){
 | 
			
		||||
	    as[level].actions.at(actionID)->init(U, pRNG);
 | 
			
		||||
	  }
 | 
			
		||||
	}
 | 
			
		||||
      }
 | 
			
		||||
@@ -123,8 +129,8 @@ namespace Grid{
 | 
			
		||||
 | 
			
		||||
	// Actions
 | 
			
		||||
	for(int level=0; level<as.size(); ++level)
 | 
			
		||||
	  for(int actionID=0; actionID<as.at(level).size(); ++actionID)
 | 
			
		||||
	    H += as[level].at(actionID)->S(U);
 | 
			
		||||
	  for(int actionID=0; actionID<as[level].actions.size(); ++actionID)
 | 
			
		||||
	    H += as[level].actions.at(actionID)->S(U);
 | 
			
		||||
	
 | 
			
		||||
	return H;
 | 
			
		||||
      }
 | 
			
		||||
 
 | 
			
		||||
@@ -27,14 +27,14 @@ namespace Grid{
 | 
			
		||||
	int fl = Integ->as.size() -1;
 | 
			
		||||
	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];
 | 
			
		||||
	for(int l=1; l<=level; ++l) fin*= 2.0*Integ->Nrel[l];
 | 
			
		||||
	int fin = Integ->as[0].multiplier;
 | 
			
		||||
	for(int l=1; l<=level; ++l) fin*= 2.0*Integ->as[l].multiplier;
 | 
			
		||||
	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 
 | 
			
		||||
	    Integ->update_P(U,level,lambda*eps);
 | 
			
		||||
@@ -101,13 +101,13 @@ namespace Grid{
 | 
			
		||||
	double eps = Integ->Params.stepsize;
 | 
			
		||||
	
 | 
			
		||||
	// 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;
 | 
			
		||||
	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;
 | 
			
		||||
	
 | 
			
		||||
	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
 | 
			
		||||
	    Integ->update_P(U, level,eps/2.0);
 | 
			
		||||
 
 | 
			
		||||
@@ -39,10 +39,9 @@ int main (int argc, char ** argv)
 | 
			
		||||
  FullSet.push_back(Level1);
 | 
			
		||||
 | 
			
		||||
  // 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);
 | 
			
		||||
  std::vector<int> rel ={1};
 | 
			
		||||
  Integrator<IntegratorAlgorithm> MDynamics(&Fine,MDpar, FullSet,rel);
 | 
			
		||||
  Integrator<IntegratorAlgorithm> MDynamics(&Fine,MDpar, FullSet);
 | 
			
		||||
 | 
			
		||||
  // Create HMC
 | 
			
		||||
  HMCparameters HMCpar;
 | 
			
		||||
 
 | 
			
		||||
		Reference in New Issue
	
	Block a user