fopr_Overlap_5d.h

Go to the documentation of this file.

#ifndef FOPR_OVERLAP_5D_INCLUDED
#define FOPR_OVERLAP_5D_INCLUDED

#include "fopr_Wilson_eo.h"
#include "math_Sign_Zolotarev.h"
#include "bridge_complex.h"

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


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

class Fopr_Overlap_5d : public Fopr
{
 private:
  //- parameters common to overlap fermion
  double             m_mq;             // quark mass
  double             m_M0;             // domain-wall height
  int                m_Np;             // number of poles in rational approx.
  double             m_x_min, m_x_max; // valid range of approximate sign function
  int                m_Niter_ms;       // max iteration of shiftsolver (dummy)
  double             m_Stop_cond_ms;   // stopping condition of shift solver (dummy)
  std::valarray<int> m_boundary;

  //- low-mode subtraction
  int m_Nsbt;
  std::valarray<double> m_ev;
  std::valarray<Field>  m_vk;

  //- following parameters are used pnly in this class
  std::valarray<double> m_cl;
  std::valarray<double> m_bl;
  std::valarray<double> m_sigma;

  std::valarray<double> m_p_sqrt;
  std::valarray<double> m_q_sqrt;
  double                m_p0_parameter, m_R_parameter, m_h, m_u0;
  std::valarray<double> m_rl;
  std::valarray<double> m_sl;

  std::valarray<double>   m_prf;
  std::valarray<dcomplex> m_u0c_e, m_u0cinv_e;
  std::valarray<dcomplex> m_u0c_o, m_u0cinv_o;

  Fopr_Wilson_eo *m_fopr_w;
  Index_eo       m_index_eo;
  std::string    m_mode;

 public:
  Fopr_Overlap_5d(Fopr_Wilson_eo *fopr) :
    Fopr(), m_fopr_w(fopr)
  {
    m_Nsbt = 0;
    m_mode = "DdagD_eo";
  }

  void set_parameters(const Parameters& params);
  void set_parameters(const double mq, const double M0, const int Np,
                      const double x_min, const double x_max,
                      const int Niter, const double Stop_cond,
                      const std::valarray<int> bc);

  void set_config(Field *U)
  {
    m_fopr_w->set_config(U);
  }

  void set_lowmodes(int Nsbt, std::valarray<double> *, std::valarray<Field> *);

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

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

  const Field mult(const Field& f)
  {
    if (m_mode == "D_eo") {
      return DD_5d_eo(f, 1);
    } else if (m_mode == "DdagD_eo") {
      return DdagD_eo(f);
    } else if (m_mode == "Ddag_eo") {
      return DD_5d_eo(f, -1);
    } else {
      vout.crucial(m_vl, "Fopr_Overlap: mode undefined\n");
      abort();
    }
  }

  const Field mult_dag(const Field& f)
  {
    if (m_mode == "D_eo") {
      return DD_5d_eo(f, -1);
    } else if (m_mode == "DdagD_eo") {
      return DdagD_eo(f);
    } else if (m_mode == "Ddag_eo") {
      return DD_5d_eo(f, 1);
    } else {
      vout.crucial(m_vl, "Fopr_Overlap: mode undefined\n");
      abort();
    }
  }

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

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

  const Field DdagD_eo(const Field&);

  const Field DD_5d_eo(const Field& w, const int jd);

  void Mopr_5d_eo(Field& v, const Field& w, const int ieo);

  void LUprecond(Field& v, const Field& w, const int ieo);

  void Proj_H_eo(const int ieo1, const int ieo2,
                 Field& v1, const Field& w1);

  void Proj_L_mult_eo(const int ieo1, const int ieo2,
                      Field& v1, const Field& w1);

  void mult_u0inv(Field& v1, const Field& w1, const int ieo);

  void Calc_Coeff_u0inv();
  void Solv_Coeff_u0inv(std::valarray<dcomplex>&,
                        const std::valarray<dcomplex>&,
                        const std::valarray<dcomplex>&);

  void innerprd_c(double& prd_r, double& prd_i,
                  const Field& v, const Field& w);

  void add_c(Field& v, const Field& w, const double a_r, const double a_i);

  double norm_c(const std::valarray<dcomplex>&);

  double innerprd_c(const std::valarray<dcomplex>&,
                    const std::valarray<dcomplex>&);

  void mult_WdagW(std::valarray<dcomplex>&,
                  const std::valarray<dcomplex>&,
                  const std::valarray<dcomplex>&);

  int field_nvol() { return m_fopr_w->field_nvol(); }
  int field_nin() { return m_fopr_w->field_nin();  }
  int field_nex() { return 2 * m_Np + 1; }
};
#endif