mirror of
https://github.com/paboyle/Grid.git
synced 2024-11-14 17:55:38 +00:00
0b603225d1
Still to check plaq and <e(-dH)>, but nevertheless this is progress
179 lines
6.9 KiB
C++
179 lines
6.9 KiB
C++
#ifndef GRID_QCD_WILSON_FERMION_H
|
|
#define GRID_QCD_WILSON_FERMION_H
|
|
|
|
namespace Grid {
|
|
|
|
namespace QCD {
|
|
|
|
class WilsonFermion : public FermionOperator<LatticeFermion,LatticeGaugeField>
|
|
{
|
|
public:
|
|
|
|
///////////////////////////////////////////////////////////////
|
|
// Implement the abstract base
|
|
///////////////////////////////////////////////////////////////
|
|
GridBase *GaugeGrid(void) { return _grid ;}
|
|
GridBase *GaugeRedBlackGrid(void) { return _cbgrid ;}
|
|
GridBase *FermionGrid(void) { return _grid;}
|
|
GridBase *FermionRedBlackGrid(void) { return _cbgrid;}
|
|
|
|
// override multiply
|
|
virtual RealD M (const LatticeFermion &in, LatticeFermion &out);
|
|
virtual RealD Mdag (const LatticeFermion &in, LatticeFermion &out);
|
|
|
|
// half checkerboard operaions
|
|
void Meooe (const LatticeFermion &in, LatticeFermion &out);
|
|
void MeooeDag (const LatticeFermion &in, LatticeFermion &out);
|
|
|
|
virtual void Mooee (const LatticeFermion &in, LatticeFermion &out); // remain virtual so we
|
|
virtual void MooeeDag (const LatticeFermion &in, LatticeFermion &out); // can derive Clover
|
|
virtual void MooeeInv (const LatticeFermion &in, LatticeFermion &out); // from Wilson base
|
|
virtual void MooeeInvDag (const LatticeFermion &in, LatticeFermion &out);
|
|
|
|
////////////////////////
|
|
//
|
|
// Force term: d/dtau S = 0
|
|
//
|
|
// It is simplest to consider the two flavour force term
|
|
//
|
|
// S[U,phi] = phidag (MdagM)^-1 phi
|
|
//
|
|
// But simplify even this to
|
|
//
|
|
// S[U,phi] = phidag MdagM phi
|
|
//
|
|
// (other options exist depending on nature of action fragment.)
|
|
//
|
|
// Require momentum be traceless anti-hermitian to move within group manifold [ P = i P^a T^a ]
|
|
//
|
|
// Define the HMC hamiltonian
|
|
//
|
|
// H = 1/2 Tr P^2 + S(U,phi)
|
|
//
|
|
// .
|
|
// U = P U (lorentz & color indices multiplied)
|
|
//
|
|
// Hence
|
|
//
|
|
// .c c c c
|
|
// U = U P = - U P (c == dagger)
|
|
//
|
|
// So, taking some liberty with implicit indices
|
|
// . . .c c
|
|
// dH/dt = 0 = Tr P P +Tr[ U dS/dU + U dS/dU ]
|
|
//
|
|
// . c c
|
|
// = Tr P P + i Tr[ P U dS/dU - U P dS/dU ]
|
|
//
|
|
// . c c
|
|
// = Tr P (P + i ( U dS/dU - P dS/dU U ]
|
|
//
|
|
// . c c
|
|
// => P = -i [ U dS/dU - dS/dU U ] generates HMC EoM
|
|
//
|
|
// Simple case work this out using S = phi^dag MdagM phi for wilson:
|
|
// c c
|
|
// dSdt = dU_xdt dSdUx + dUxdt dSdUx
|
|
//
|
|
// = Tr i P U_x [ (\phi^\dag)_x (1+g) (M \phi)_x+\mu +(\phi^\dag M^\dag)_x (1-g) \phi_{x+\mu} ]
|
|
// c
|
|
// - i U_x P [ (\phi^\dag)_x+mu (1-g) (M \phi)_x +(\phi^\dag M^\dag)_(x+\mu) (1+g) \phi_{x} ]
|
|
//
|
|
// = i [(\phi^\dag)_x ]_j P_jk [U_x(1+g) (M \phi)_x+\mu]_k (1)
|
|
// + i [(\phi^\dagM^\dag)_x]_j P_jk [U_x(1-g) (\phi)_x+\mu]_k (2)
|
|
// - i [(\phi^\dag)_x+mu (1-g) U^dag_x]_j P_jk [(M \phi)_xk (3)
|
|
// - i [(\phi^\dagM^\dag)_x+mu (1+g) U^dag_x]_j P_jk [ \phi]_xk (4)
|
|
//
|
|
// Observe that (1)* = (4)
|
|
// (2)* = (3)
|
|
//
|
|
// Write as .
|
|
// P_{kj} = - i ( [U_x(1+g) (M \phi)_x+\mu] (x) [(\phi^\dag)_x] + [U_x(1-g) (\phi)_x+\mu] (x) [(\phi^\dagM^\dag)_x] - h.c )
|
|
//
|
|
// where (x) denotes outer product in colour and spins are traced.
|
|
//
|
|
// Need only evaluate (1) and (2) [Chroma] or (2) and (4) [IroIro] and take the
|
|
// traceless anti hermitian part (of term in brackets w/o the "i")
|
|
//
|
|
// Generalisation to S=phi^dag (MdagM)^{-1} phi is simple:
|
|
//
|
|
// For more complicated DWF etc... apply product rule in differentiation
|
|
//
|
|
////////////////////////
|
|
void DhopDeriv (LatticeGaugeField &mat,const LatticeFermion &U,const LatticeFermion &V,int dag);
|
|
void DhopDerivEO(LatticeGaugeField &mat,const LatticeFermion &U,const LatticeFermion &V,int dag);
|
|
void DhopDerivOE(LatticeGaugeField &mat,const LatticeFermion &U,const LatticeFermion &V,int dag);
|
|
|
|
// Extra support internal
|
|
void DerivInternal(CartesianStencil & st,
|
|
LatticeDoubledGaugeField & U,
|
|
LatticeGaugeField &mat,
|
|
const LatticeFermion &A,
|
|
const LatticeFermion &B,
|
|
int dag);
|
|
|
|
|
|
// non-hermitian hopping term; half cb or both
|
|
void Dhop (const LatticeFermion &in, LatticeFermion &out,int dag);
|
|
void DhopOE(const LatticeFermion &in, LatticeFermion &out,int dag);
|
|
void DhopEO(const LatticeFermion &in, LatticeFermion &out,int dag);
|
|
|
|
// Multigrid assistance
|
|
void Mdir (const LatticeFermion &in, LatticeFermion &out,int dir,int disp);
|
|
void DhopDir(const LatticeFermion &in, LatticeFermion &out,int dir,int disp);
|
|
void DhopDirDisp(const LatticeFermion &in, LatticeFermion &out,int dirdisp,int gamma,int dag);
|
|
|
|
///////////////////////////////////////////////////////////////
|
|
// Extra methods added by derived
|
|
///////////////////////////////////////////////////////////////
|
|
void DhopInternal(CartesianStencil & st,
|
|
LatticeDoubledGaugeField &U,
|
|
const LatticeFermion &in,
|
|
LatticeFermion &out,
|
|
int dag);
|
|
|
|
// Constructor
|
|
WilsonFermion(LatticeGaugeField &_Umu,GridCartesian &Fgrid,GridRedBlackCartesian &Hgrid,RealD _mass);
|
|
|
|
// DoubleStore
|
|
virtual void ImportGauge(const LatticeGaugeField &_Umu);
|
|
void DoubleStore(LatticeDoubledGaugeField &Uds,const LatticeGaugeField &Umu);
|
|
|
|
///////////////////////////////////////////////////////////////
|
|
// Data members require to support the functionality
|
|
///////////////////////////////////////////////////////////////
|
|
static int HandOptDslash; // these are a temporary hack
|
|
static int MortonOrder;
|
|
|
|
// protected:
|
|
public:
|
|
|
|
RealD mass;
|
|
|
|
GridBase * _grid;
|
|
GridBase * _cbgrid;
|
|
|
|
static const int npoint=8;
|
|
static const std::vector<int> directions ;
|
|
static const std::vector<int> displacements;
|
|
|
|
//Defines the stencils for even and odd
|
|
CartesianStencil Stencil;
|
|
CartesianStencil StencilEven;
|
|
CartesianStencil StencilOdd;
|
|
|
|
// Copy of the gauge field , with even and odd subsets
|
|
LatticeDoubledGaugeField Umu;
|
|
LatticeDoubledGaugeField UmuEven;
|
|
LatticeDoubledGaugeField UmuOdd;
|
|
|
|
// Comms buffer
|
|
std::vector<vHalfSpinColourVector,alignedAllocator<vHalfSpinColourVector> > comm_buf;
|
|
|
|
|
|
};
|
|
|
|
}
|
|
}
|
|
#endif
|