alg_nuc3pt.h

Go to the documentation of this file.

//------------------------------------------------------------------
//
// alg_nuc3pt.h
//
// Header file for all alg classes relevant to Wilson-type fermion
// nuc3pt' spectrum. The type of glue or fermion is given as
// an argument of type Lattice& to the constructor. If the 
// fermion type is not F_CLASS_WILSON or F_CLASS_CLOVER or 
// F_CLASS_DWF the QPropW constructors exit with a general error.
//
//------------------------------------------------------------------


#ifndef INCLUDED_ALG_NUC3PT_H
#define INCLUDED_ALG_NUC3PT_H

#include <util/lattice.h>
#include <alg/alg_base.h>

#include <alg/qpropw.h>

#include <alg/corrfunc.h>
#include <util/momentum.h>
#include <alg/nuc3pt_arg.h>
#include <alg/nuc2pt.h>
#include <alg/nuc3pt.h>
#include <alg/nuc3pt_gammar.h>

#if TARGET == QCDOC
#include <comms/sysfunc_cps.h>
#endif


CPS_START_NAMESPACE


enum DIR {
 X =  0, 
 Y =  1,
 Z =  2,
 T =  3,
 G5 = -5
} ;


//------------------------------------------------------------------
//
// AlgNuc3pt is derived from Alg and is relevant to  
// meson three point functions with Wilson type 
// fermions (i.e. Wilson, Clover, Dwf). 
// The type of fermion is determined by the argument to the 
// constructor. If the fermion type is not F_CLASS_WILSON or 
// F_CLASS_CLOVER or F_CLASS_DWF the constructors exit with a 
// general error.
//
//------------------------------------------------------------------
class AlgNuc3pt : public Alg
{
 private:
    char* cname;

    Nuc3ptArg* Nuc3pt_arg ;
        // The argument structure for the
        // three point calculation

    QPropW* q_prop ; // Stores the quark propagator
    QPropWSeqBar* u_s_prop ; // Stores the up   quark sequential propagator
    QPropWSeqBar* d_s_prop ; // Stores the down quark sequential propagator

    FILE *fp ; // The I/O file pointer

    void GetThePropagator(int,Float) ;
    void GetTheSeqPropagator(int, Float, SourceType, int*, ProjectType) ;
    
    //void QIO_SaveProp(QPropW& prop, char*, char*, int) ;
    //void QIO_ReadProp(QPropW& prop, char*) ;

    void OpenFile();
    void CloseFile() 
      {
        if(fp != NULL)
          Fclose(fp);
      }

   
    void calc_Scalar()
      {     
        Gamma S;
        Nuc3ptGamma Scal(S) ;
        Scal.Calc3pt(*u_s_prop,*q_prop);
        Scal.Calc3pt(*d_s_prop,*q_prop);
        OpenFile(); 
        Scal.Print(fp) ;
        CloseFile();
      }
    


#ifdef _NEDM

    void calc_NEDM(){ 
      Gamma Gt(T);
      Nuc3ptGammaR NEDM(Gt,Z) ;
      NEDM.Calc3pt(*u_s_prop,*q_prop);
      NEDM.Calc3pt(*d_s_prop,*q_prop);
      OpenFile(); 
      NEDM.Print(fp) ;
      CloseFile();
    }

    void calc_NMM(){
      Gamma Gx(X);
      Nuc3ptGammaR NMM(Gx,Y) ;
      NMM.Calc3pt(*u_s_prop,*q_prop);
      NMM.Calc3pt(*d_s_prop,*q_prop);
      OpenFile(); 
      NMM.Print(fp) ;
      CloseFile();
    }

#endif

    void calc_Vector(const ThreeMom& q) ;

    void calc_Axial(const ThreeMom& q) ;


    void calc_PScalar(const ThreeMom& q) ;

    void calc_EnergyMomentum(const ThreeMom& mom) ;
    
    //UnPolarized
    void calc_Vector()
      {     
        Gamma Gt(T);
        Nuc3ptGamma VectCurr(Gt) ;
        VectCurr.Calc3pt(*u_s_prop,*q_prop);
        VectCurr.Calc3pt(*d_s_prop,*q_prop);
        OpenFile(); 
        VectCurr.Print(fp) ;
        CloseFile();
      }

    void calc_X_q_a()
      {
        Gamma Gx(X);
        Gamma Gt(T);
        Derivative Der_t(T);
        Derivative Der_x(T);

        Nuc3ptStru Xq_xt(Gx, Der_t);
        Xq_xt.Calc3pt(*u_s_prop, *q_prop);
        Xq_xt.Calc3pt(*d_s_prop, *q_prop);

        Nuc3ptStru Xq_tx(Gt,Der_x);
        Xq_tx.Calc3pt(*u_s_prop, *q_prop);
        Xq_tx.Calc3pt(*d_s_prop, *q_prop);
        
        Xq_xt += Xq_tx ;

        OpenFile(); 
        Fprintf(fp,"The next is: 14 + 41\n");
        Xq_xt.Print(fp) ; 
        CloseFile();
      }

    void calc_X_q_b()
      {
        Gamma Gt(T) ;
        Derivative Der_t(T) ;
        Nuc3ptStru Xq_tt(Gt,Der_t) ;
        Xq_tt.Calc3pt(*u_s_prop,*q_prop);
        Xq_tt.Calc3pt(*d_s_prop,*q_prop);
         
        for(int k(X);k<T;k++)
          {
            Gamma Gk(k) ;
            Derivative Der_k(k) ;
            Nuc3ptStru tmp(Complex(-1.0/3.0,0.0), Gk, Der_k) ;
            tmp.Calc3pt(*u_s_prop,*q_prop);
            tmp.Calc3pt(*d_s_prop,*q_prop);
            Xq_tt+=tmp ;
          }
        OpenFile(); 
        Fprintf(fp,"The next is: 44 - 1/3 (11 + 22 + 33)\n");
        Xq_tt.Print(fp) ; 
        CloseFile();
      }

    void calc_X2_q()
      {
        Gamma Gt(T);
        Derivative D_xx(X,X);
        Nuc3ptStru X2q_txx(Gt,D_xx) ;
        X2q_txx.Calc3pt(*u_s_prop,*q_prop);
        X2q_txx.Calc3pt(*d_s_prop,*q_prop);
         
        for(int k(Y);k<T;k++)
          {
            Derivative D_kk(k,k);
            Nuc3ptStru tmp(Complex(-0.5,0.0),Gt,D_kk) ;
            tmp.Calc3pt(*u_s_prop,*q_prop);
            tmp.Calc3pt(*d_s_prop,*q_prop);
            X2q_txx+=tmp ;
          }
        OpenFile(); 
        Fprintf(fp,"The next is: 411 - 1/2(422+433)\n");
        X2q_txx.Print(fp) ; 
        CloseFile();    
      }

    void calc_X3_q()
      {
        Gamma Gx(X);
        Gamma Gy(Y);

        Derivative D_xtt(X,T,T);
        Nuc3ptStru X3q_xxtt(Gx,D_xtt) ;
        X3q_xxtt.Calc3pt(*u_s_prop,*q_prop);
        X3q_xxtt.Calc3pt(*d_s_prop,*q_prop);
        
        {
          Derivative D_yzz(Y,Z,Z) ;
          Nuc3ptStru tmp(Gy,D_yzz) ;
          tmp.Calc3pt(*u_s_prop,*q_prop);
          tmp.Calc3pt(*d_s_prop,*q_prop);
          X3q_xxtt += tmp ;
        }
        
        {
          Derivative D_xzz(X,Z,Z) ;
          Nuc3ptStru tmp(Complex(-1.0,0.0),Gx,D_xzz) ;
          tmp.Calc3pt(*u_s_prop,*q_prop);
          tmp.Calc3pt(*d_s_prop,*q_prop);
          X3q_xxtt += tmp ;
        }
        
        {
          Derivative D_ytt(Y,T,T) ;
          Nuc3ptStru tmp(Complex(-1.0,0.0),Gy, D_ytt) ;
          tmp.Calc3pt(*u_s_prop,*q_prop);
          tmp.Calc3pt(*d_s_prop,*q_prop);
          X3q_xxtt += tmp ;
        }

        
        OpenFile(); 
        Fprintf(fp,"The next is: 1144 + 2233 - 1133 - 2244\n");
        X3q_xxtt.Print(fp) ; 
        CloseFile();
      }

    //Polarized

    void calc_Axial()
      {
        Gamma G5z(G5,Z);
        Nuc3ptGamma AxialCurr(Complex(0.0,1.0),G5z) ;
        AxialCurr.Calc3pt(*u_s_prop,*q_prop);
        AxialCurr.Calc3pt(*d_s_prop,*q_prop);
        
        OpenFile();
        AxialCurr.Print(fp) ;
        CloseFile();
      }


    void calc_X_Dq_a()
      {
        Nuc3ptStru XDq_xz(Complex(0.0,1.0),Gamma(G5,X),Derivative(Z)) ;
        XDq_xz.Calc3pt(*u_s_prop,*q_prop);
        XDq_xz.Calc3pt(*d_s_prop,*q_prop);
         
        Nuc3ptStru XDq_zx(Complex(0.0,1.0),Gamma(G5,Z),Derivative(X)) ;
        XDq_zx.Calc3pt(*u_s_prop,*q_prop);
        XDq_zx.Calc3pt(*d_s_prop,*q_prop);

        XDq_xz+=XDq_xz ;

        OpenFile();
        Fprintf(fp,"The next is: 13 + 31\n");
        XDq_zx.Print(fp);
        CloseFile();
      }

    void calc_X_Dq_b()
      {
        Nuc3ptStru XDq_zt(Complex(0.0,1.0),Gamma(G5,Z),Derivative(T)) ;
        XDq_zt.Calc3pt(*u_s_prop,*q_prop);
        XDq_zt.Calc3pt(*d_s_prop,*q_prop);
         
        Nuc3ptStru XDq_tz(Complex(0.0,1.0),Gamma(G5,T),Derivative(Z)) ;
        XDq_tz.Calc3pt(*u_s_prop,*q_prop);
        XDq_tz.Calc3pt(*d_s_prop,*q_prop);

        XDq_zt+=XDq_tz ;

        XDq_zt.setTag("p") ;
        OpenFile();
        Fprintf(fp,"The next is: 34 + 43\n");
        XDq_zt.Print(fp);
        CloseFile();

      }


    void calc_X2_Dq()
      {
      }

    void calc_d1() // Could be compined with X_Dq_a
      {
        Nuc3ptStru d1_zt(Complex(0.0,1.0),Gamma(G5,Z),Derivative(T)) ;
        d1_zt.Calc3pt(*u_s_prop,*q_prop);
        d1_zt.Calc3pt(*d_s_prop,*q_prop);
         
        Nuc3ptStru XDq_tz(Complex(0.0,-1.0),Gamma(G5,T),Derivative(Z)) ;
        XDq_tz.Calc3pt(*u_s_prop,*q_prop);
        XDq_tz.Calc3pt(*d_s_prop,*q_prop);

        d1_zt+=XDq_tz ;
        d1_zt.setTag("m") ;

        OpenFile();
        Fprintf(fp,"The next is: 34 - 43\n");
        d1_zt.Print(fp);
        CloseFile();
      }
    void calc_d2()
      {
      }

    
    //Transversity

    void calc_Tensor()
      {
        Nuc3ptGamma Tensor(Complex(0.0,1.0),Gamma(G5,Z,T)) ;
        Tensor.Calc3pt(*u_s_prop,*q_prop);
        Tensor.Calc3pt(*d_s_prop,*q_prop);

        OpenFile();
        Tensor.Print(fp) ;
        CloseFile();
      }
    
    void calc_X_dq()
      {
        Nuc3ptStru Xdq_ztx(Complex(0.0,1.0),Gamma(G5,Z,T),Derivative(X)) ;
        Xdq_ztx.Calc3pt(*u_s_prop, *q_prop);
        Xdq_ztx.Calc3pt(*d_s_prop, *q_prop);
         
        Nuc3ptStru Xdq_zxt(Complex(0.0,1.0),Gamma(G5,Z,X),Derivative(T)) ;
        Xdq_zxt.Calc3pt(*u_s_prop, *q_prop);
        Xdq_zxt.Calc3pt(*d_s_prop, *q_prop);
        
        Xdq_ztx += Xdq_zxt ;

        OpenFile();
        Fprintf(fp,"The next is: 341 + 314\n");
        Xdq_ztx.Print(fp) ;
        CloseFile();
      }

    void calc_Cons_Vector(int i, ThreeMom*) ;

    void calc_Cons_Axial_Vector(int i, ThreeMom*) ;


 public:
    AlgNuc3pt(Lattice & latt, CommonArg* c_arg, Nuc3ptArg* arg);

    virtual ~AlgNuc3pt();

    void run();

    void test_run();

    void run_hbd();

    void run_2pt();
};

CPS_END_NAMESPACE

#endif

---