hmc_Leapfrog.h

Go to the documentation of this file.

#ifndef HMC_LEAPFROG_INCLUDED
#define HMC_LEAPFROG_INCLUDED

#include "action.h"
#include "director.h"
#include "integrator.h"
#include "randomNumbers.h"
#include "staples.h"

#include "bridgeIO.h"
using Bridge::vout;


//- parameters class
class Parameters_HMC_Leapfrog : virtual public Parameters
{
 public:
  Parameters_HMC_Leapfrog();
};
//- end

class HMC_Leapfrog
{
 private:
  int    m_Nmdc;
  int    m_Nprec;
  int    m_Mtpl_test;
  double m_Estep;

  Staples *m_staple;

  std::valarray<Action *>   m_action;
  std::valarray<Director *> m_director;
  RandomNumbers             *m_rand;
  Bridge::VerboseLevel      m_vl;

 public:
  HMC_Leapfrog(std::valarray<Action *> action,
               RandomNumbers *rand)
    : m_vl(CommonParameters::Vlevel())
  {
    m_action.resize(action.size());
    for (int i = 0; i < action.size(); ++i) {
      m_action[i] = action[i];
    }
    m_staple    = new Staples;
    m_rand      = rand;
    m_Estep     = 0.0;
    m_Nmdc      = 0;
    m_Nprec     = 0;
    m_Mtpl_test = 0;
  }

  HMC_Leapfrog(std::valarray<Action *> action,
               std::valarray<Director *> director,
               RandomNumbers *rand)
    : m_vl(CommonParameters::Vlevel())
  {
    m_action.resize(action.size());
    for (int i = 0; i < action.size(); ++i) {
      m_action[i] = action[i];
    }
    m_director.resize(director.size());
    for (int i = 0; i < director.size(); ++i) {
      m_director[i] = director[i];
    }
    m_staple    = new Staples;
    m_rand      = rand;
    m_Estep     = 0.0;
    m_Nmdc      = 0;
    m_Nprec     = 0;
    m_Mtpl_test = 0;
  }

  ~HMC_Leapfrog()
  {
    delete m_staple;
  }

  void set_parameters(const Parameters& params);
  void set_parameters(double Estep, int Nmdc, int Nprec, int Mtpl_test);

  void set_parameter_verboselevel(const Bridge::VerboseLevel vl) { m_vl = vl; }

  //  Mtpl_test=0: no test, !=0: test
  double update(Field_G&);

  double langevin(Field_G& iP, Field_G& U);
  double langevin_P(Field_G& iP);

  double calc_Hamiltonian(Field_G& iP, Field_G& U);
  double calcH_P(Field_G& iP);

  void integrate(Field_G& iP, Field_G& U);

  void update_U(double estep, Field_G& iP, Field_G& U);
  void update_P(double estep, Field_G& iP, Field_G& U);
};
#endif