alg_int.h

Go to the documentation of this file.

#include<config.h>
CPS_START_NAMESPACE
//------------------------------------------------------------------
//------------------------------------------------------------------


#ifndef INCLUDED_ALG_INT_H
#define INCLUDED_ALG_INT_H

CPS_END_NAMESPACE
#include <util/lattice.h>
#include <util/smalloc.h>
#include <util/pmalloc.h>
#include <alg/alg_base.h>
#include <alg/alg_meas.h>
#include <alg/cg_stats.h>
#include <util/checksum.h>
#include <alg/int_arg.h>
#include <alg/remez_arg.h>
#include <util/qcdio.h>
#include <alg/alg_eig.h> 
CPS_START_NAMESPACE

class AlgInt {

 private:
  char *cname;

 protected:
  int traj;
  IntegratorType int_type;

 public:
  
  AlgInt();
  virtual ~AlgInt();

  virtual void heatbath() = 0;

  virtual void evolve (Float dt, int steps) = 0;

  virtual Float energy() = 0;

  virtual void cost(CgStats*) = 0;

  virtual void reverse() = 0;

  virtual void init() = 0;

};

class AlgIntAB : public AlgInt {

 private:
  char *cname;

 protected:
  AlgInt *A;
  AlgInt *B;
  int A_steps, B_steps;
  int A_calls, B_calls;
  
  IntegratorLevel level; 
  unsigned long step_cnt;
  IntABArg *ab_arg;

 public:
  AlgIntAB(AlgInt &A, AlgInt &B, IntABArg &);
  virtual ~AlgIntAB();

  void heatbath();

  Float energy();

  virtual void evolve(Float dt, int steps) = 0;

  void cost(CgStats*);
  void reverse();

  void init();

  static AlgIntAB& Create(AlgInt &A, AlgInt &B, IntABArg &ab_arg);
  static void Destroy(AlgIntAB&);

};

class AlgIntLeap : public AlgIntAB {

private:
  char *cname;

public:
  AlgIntLeap(AlgInt &A, AlgInt &B, IntABArg &);
  virtual ~AlgIntLeap();

  void evolve(Float dt, int steps);
  
};

class AlgIntOmelyan : public AlgIntAB {

private:
  char *cname;
  Float lambda;

public:
  AlgIntOmelyan(AlgInt &A, AlgInt &B, IntABArg &);
  virtual ~AlgIntOmelyan();

  void evolve(Float dt, int steps);
  
};

class AlgIntCampostrini : public AlgIntAB {

private:
  char *cname;
  Float sigma;
  Float epsilon;

public:
  AlgIntCampostrini(AlgInt &A, AlgInt &B, IntABArg &);
  virtual ~AlgIntCampostrini();

  void evolve(Float dt, int steps);
  
};

class AlgIntOmelyan44 : public AlgIntAB {

private:
  char *cname;
  Float rho;
  Float theta;
  Float lambda;

public:
  AlgIntOmelyan44(AlgInt &A, AlgInt &B, IntABArg &);
  virtual ~AlgIntOmelyan44();

  void evolve(Float dt, int steps);
  
};

class AlgIntOmelyan45 : public AlgIntAB {

private:
  char *cname;
  Float theta;
  Float rho;
  float lambda;
  Float mu;

public:
  AlgIntOmelyan45(AlgInt &A, AlgInt &B, IntABArg &);
  virtual ~AlgIntOmelyan45();

  void evolve(Float dt, int steps);
  
};

class AlgIntSum : public AlgIntAB {

private:
  char *cname;

public:
  AlgIntSum(AlgInt &A, AlgInt &B, IntABArg &);
  virtual ~AlgIntSum();

  void evolve(Float dt, int steps);

};

class AlgHamiltonian : public AlgInt {

 protected:
  int g_size;

 private:
  char *cname;

 public:
  AlgHamiltonian();
  virtual ~AlgHamiltonian();

  virtual void heatbath() = 0;
  virtual Float energy() = 0;
  virtual void evolve(Float dt, int steps) = 0;
  virtual void cost(CgStats*) = 0;

};

class AlgMomentum : public AlgHamiltonian {

 private:
  char *cname;

  char *md_time_str;
  Matrix *mom;

 public:
  AlgMomentum();
  virtual ~AlgMomentum();

  void heatbath();

  Float energy();

  void evolve(Float dt, int steps);

  void cost(CgStats*);

  Matrix* getMom();

  void reverse();

  void init();
};

class AlgAction : public AlgHamiltonian {

 private:
  char *cname;

 protected:
  Matrix *mom;
  ForceMeasure force_measure;
  char *force_label;
  ForceArg Fdt;

 public:
  AlgAction();
  AlgAction(AlgMomentum &mom, ActionArg &action_arg);
  ~AlgAction();

  virtual void heatbath() = 0;
  virtual Float energy() = 0;
  virtual void evolve(Float dt, int steps) = 0;
  virtual void cost(CgStats*) = 0;
  void reverse();

};

class AlgActionBilinear : public AlgAction {

 private:
  char *cname;

 protected:
  ActionBilinearArg *bi_arg;

  int n_masses;
  FclassType fermion;

  Float *mass;

  int *max_num_iter;

  int f_sites;
  int f_vec_count;
  int f_size;
  int Ncb;

  CgStats cg_stats;
  int cg_iter;

  Vector **phi;

  int md_steps;

 public:
  AlgActionBilinear();
  AlgActionBilinear(AlgMomentum &, ActionBilinearArg &);
  virtual ~AlgActionBilinear();
  void cost(CgStats*);
  void updateCgStats(CgArg*);

  int getNmass();
  Float getMass(int);
  FclassType getFermion();

  virtual void heatbath() = 0;
  virtual Float energy() = 0;
  virtual void evolve(Float dt, int steps) = 0;

  void init();

};

class AlgActionRational : public AlgActionBilinear {

 private:
  char *cname;
  int **fractionSplit;
  int **splitCheck;
  ActionRationalArg *rat_arg;

  RemezArg *remez_arg_md;
  RemezArg *remez_arg_mc;

 protected:

  int evolved;

  int heatbathEval;

  int energyEval;

  Float h_init;

  CgArg ***frm_cg_arg_md;
  CgArg ***frm_cg_arg_mc;

  Vector** frmn;

  Vector** frmn_d;

  Vector** frmn_tmp;

  int total_size;

  int max_size;
  
  Float *all_res;

  EigArg eig_arg;

  char eig_file[256];

  CommonArg ca_eig;

  Float **lambda_low;
  Float **lambda_high;

  void generateApprox(Float *mass, RemezArg **remez_arg_md, 
                      RemezArg **remez_arg_mc, RationalDescr *rat);
  
  //<! Free approximations
  void destroyApprox(RemezArg *remez_arg_md, RemezArg *remez_arg_mc);
  
  void generateCgArg(Float *mass, CgArg ****cg_arg_md, 
                     CgArg ****cg_arg_mc, char *label, 
                     RationalDescr *rat_des);

  //<! Free cg args
  void destroyCgArg(CgArg ***cg_arg_md, CgArg ***cg_arg_mc,
                    char *label, RemezArg *remez_arg_md,
                    RemezArg *remez_arg_mc);

  void generateEigArg(EigenDescr eigen);

  void destroyEigArg();

  void checkApprox(Float *mass, RemezArg *remez_arg, EigenDescr eigen);

 public:

  AlgActionRational();
  AlgActionRational(AlgMomentum &mom, ActionRationalArg &rat_arg, int traj_num=0);
  AlgActionRational(AlgMomentum &mom, ActionBilinearArg &bi_arg);
  virtual ~AlgActionRational();

  void heatbath();
  Float energy();

  void evolve(Float dt, int steps);
  void evolve(Float dt, int steps, int **fractionSplit);

  static int compareApprox(RemezArg &, RemezArg &);

  void init(int traj_num);

  void setSplit(int i, int j);
  void checkSplit();
  
};

class AlgActionRationalSplit : public AlgActionRational {

 private:
  char *cname;
  AlgActionRational *rat;
  ActionRationalSplitArg *rat_split_arg;
  int **fractionSplit;

 public:
  AlgActionRationalSplit(AlgActionRational &rat, 
                         ActionRationalSplitArg &rat_split_arg); 
                         
  virtual ~AlgActionRationalSplit();

  void heatbath();
  Float energy();

  void evolve(Float dt, int steps);
  
  int getNmass();
  Float getMass(int);

  void cost(CgStats*);
};

class AlgActionBoson : public AlgActionBilinear {

 private:
  char *cname;

  ActionBosonArg *bsn_arg;
  CgArg **bsn_cg_arg;   

 public:

  AlgActionBoson(AlgMomentum &mom, ActionBosonArg &boson_arg);
  virtual ~AlgActionBoson();
  
  void heatbath();

  Float energy();

  void evolve(Float dt, int steps);

};

class AlgActionFermion : public AlgActionBilinear {

 private:
  char *cname;

  ActionFermionArg *frm_arg;

  CgArg **frm_cg_arg_md;   
  CgArg **frm_cg_arg_mc;   

  int evolved;
  Float h_init;

  Vector ***v;
  Vector ***cg_sol_old;
  Vector *cg_sol;

  Vector ***vm;

  int *chrono;

 public:

  AlgActionFermion(AlgMomentum &mom, ActionFermionArg &frm_arg);
  virtual ~AlgActionFermion();
  
  void heatbath();

  Float energy();

  void evolve(Float dt, int steps);

  void init();
};

class AlgActionQuotient : public AlgActionBilinear {

 private:
  char *cname;

  ActionQuotientArg *quo_arg;

  CgArg **bsn_cg_arg;   
  CgArg **frm_cg_arg_md;   
  CgArg **frm_cg_arg_mc;   

  Float *bsn_mass;  
      // ~~added for twisted mass Wilson fermions
  Float *bsn_mass_epsilon;  
  Float *frm_mass;  
      // ~~added for twisted mass Wilson fermions
  Float *frm_mass_epsilon;  

  int evolved;
  Float h_init;

  Vector ***v;
  Vector ***cg_sol_old;
  Vector *cg_sol;
  Vector *tmp1;
  Vector *tmp2;

  Vector ***vm;

  int *chrono;

 public:

  AlgActionQuotient(AlgMomentum &mom, ActionQuotientArg &frm_arg);
  virtual ~AlgActionQuotient();
  
  void reweight(Float *rw_fac,Float *norm);
  void heatbath();

  Float energy();

  void evolve(Float dt, int steps);

  void init();
};

class AlgActionRationalQuotient : public AlgActionRational {

 private:
  char *cname;
  int **fractionSplit;
  ActionRationalQuotientArg *rat_quo_arg;

 protected:

  Float *bsn_mass;  
  Float *frm_mass;  

  RemezArg *frm_remez_arg_md;
  RemezArg *frm_remez_arg_mc;
  RemezArg *bsn_remez_arg_md;
  RemezArg *bsn_remez_arg_mc;

  CgArg ***bsn_cg_arg_md;
  CgArg ***bsn_cg_arg_mc;

  Vector **eta; 

 public:

  AlgActionRationalQuotient();
  AlgActionRationalQuotient(AlgMomentum &mom, 
                            ActionRationalQuotientArg &rat_quo_arg, int traj_num=0);
  virtual ~AlgActionRationalQuotient();

  void reweight(Float *rw_fac,Float *norm);
  void heatbath();
  Float energy();

  void evolve(Float dt, int steps);
  void evolve(Float dt, int steps, int **fractionSplit);

};

class AlgActionGauge : public AlgAction {

 private:
  char *cname;
  ActionGaugeArg *gauge_arg;
  GclassType gluon;

 public:

  AlgActionGauge(AlgMomentum &mom, ActionGaugeArg &gauge_arg);
  virtual ~AlgActionGauge();

  void heatbath();

  Float energy();

  void evolve(Float dt, int steps);

  void cost(CgStats*);

  void init();

};

#endif

CPS_END_NAMESPACE

---