fopr_CloverTerm_eo.h

Go to the documentation of this file.

#ifndef FOPR_CLOVERTERM_EO_IMP_INCLUDED
#define FOPR_CLOVERTERM_EO_IMP_INCLUDED

#include <vector>

#include "Fopr/fopr.h"
#include "Tools/gammaMatrixSet.h"
#include "Field/index_eo.h"
#include "Field/shiftField_eo.h"

#include "Tools/mat_SU_N.h"
#include "Tools/vec_SU_N.h"

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


namespace Imp {
  class Fopr_CloverTerm_eo : public Fopr
  {
   public:
    static const std::string class_name;

// This class return D_ee = 1-f_ee or D_oo = 1-f_oo
   private:
    double           m_kappa;
    double           m_cSW;
    std::vector<int> m_boundary;
    std::string      m_repr;
    std::string      m_mode;

    int m_Nvol, m_Nvol2;
    int m_Ndim;
    int m_Nc, m_Nd, m_Ndm2;
    int m_NinF;

    std::vector<GammaMatrix> m_GM, m_SG;

    void (Fopr_CloverTerm_eo::*m_mult)(Field&, const Field&);

    Index_eo idx;

    const Field_G *m_Ueo;

    Index_eo      m_idx;
    ShiftField_eo m_shift_eo;

    Field_F  *m_fee_inv;
    Field_F  *m_foo_inv;
    Field_F  m_vf, m_ff;
    Vec_SU_N v1, v2;

    Field_G m_T;

    std::vector<Field_G> m_T2;

   public:
    Fopr_CloverTerm_eo(std::string repr)
      : m_Nvol(CommonParameters::Nvol()),
        m_Nvol2(m_Nvol / 2),
        m_Ndim(CommonParameters::Ndim()),
        m_Nc(CommonParameters::Nc()),
        m_Nd(CommonParameters::Nd()),
        m_Ndm2(m_Nd * m_Nd / 2),
        m_T(m_Nvol, m_Ndm2)
    {
      init(repr);
    }

    ~Fopr_CloverTerm_eo()
    {
      delete m_fee_inv;
      delete m_foo_inv;
    }

    void set_parameters(const Parameters& params);
    void set_parameters(const double kappa, const double cSW,
                        const std::vector<int> bc);

    void set_config(Field *Ueo);

    void set_config(unique_ptr<Field_G>& Ueo)
    {
      set_config(Ueo.get());
    }

    void set_mode(std::string mode)
    {
      m_mode = mode;
    }

    std::string get_mode() const
    {
      return m_mode;
    }

    const Field mult(const Field& f)
    {
      Field v(f.nin(), f.nvol(), f.nex());

      mult(v, f);
      return v;
    }

    const Field mult_dag(const Field& f)
    {
      Field v(f.nin(), f.nvol(), f.nex());

      mult_dag(v, f);
      return v;
    }

    void mult(Field& v, const Field& f)
    {
      if (m_mode == "even") {
        D(v, f, 0);
      } else if (m_mode == "odd") {
        D(v, f, 1);
      } else {
        vout.crucial("Error at %s: undefined mode = %s\n", class_name.c_str(), m_mode.c_str());
        exit(EXIT_FAILURE);
      }
    }

    void mult_dag(Field& v, const Field& f)
    {
      mult(v, f);
    }

    void mult_isigma(Field_F&, const Field_F&,
                     const int mu, const int nu);

    void D(Field& v, const Field& f, const int ieo);

    void D_dirac(Field& v, const Field& f, const int ieo);

    void D_chiral(Field& v, const Field& f, const int ieo);

    //const Field_F mult_csw_inv(const Field_F&, const int ieo);
    //const Field_G trSigmaInv(const int mu, const int nu);
    void trSigmaInv(Field_G&, const int mu, const int nu);

    void mult_csw_inv(Field&, const Field&, const int ieo);

    void mult_csw_inv_dirac(Field&, const Field&, const int ieo);

    void mult_csw_inv_chiral(Field&, const Field&, const int ieo);

    std::vector<double> csmatrix(const int&);

    int field_nvol() { return m_Nvol2; }
    int field_nin() { return 2 * m_Nc * m_Nd; }
    int field_nex() { return 1; }

    double flop_count();

   private:
    void init(std::string repr);
    void tidyup();

    void solve_csw_inv();

    void set_csw();

    void set_csw_dirac();

    void set_csw_chiral();

    void mult_csw(Field_F&, const Field_F&, const int ieo);
    void set_fieldstrength(Field_G&, const int, const int);

    int sg_index(int mu, int nu) { return mu * m_Ndim + nu; }
  };
}
#endif