mirror of
https://github.com/paboyle/Grid.git
synced 2025-04-10 14:10:46 +01:00
First conserved current implementation for Wilson fermions only. Not implemented for Gparity or 5D-vectorised Wilson fermions.
This commit is contained in:
Lanny91
parent
1425afc72f
commit
44260643f6
@ -489,6 +489,14 @@ namespace QCD {
|
||||
return traceIndex<ColourIndex>(lhs);
|
||||
}
|
||||
|
||||
//////////////////////////////////////////
|
||||
// Current types
|
||||
//////////////////////////////////////////
|
||||
GRID_SERIALIZABLE_ENUM(Current, undef,
|
||||
Vector, 0,
|
||||
Axial, 1,
|
||||
Tadpole, 2);
|
||||
|
||||
} //namespace QCD
|
||||
} // Grid
|
||||
|
||||
|
||||
@ -112,6 +112,21 @@ namespace Grid {
|
||||
///////////////////////////////////////////////
|
||||
virtual void ImportGauge(const GaugeField & _U)=0;
|
||||
|
||||
//////////////////////////////////////////////////////////////////////
|
||||
// Conserved currents, either contract at sink or insert sequentially.
|
||||
//////////////////////////////////////////////////////////////////////
|
||||
virtual void ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu)=0;
|
||||
virtual void SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax)=0;
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
@ -395,6 +395,31 @@ void ImprovedStaggeredFermion<Impl>::DhopInternal(StencilImpl &st, LebesgueOrder
|
||||
}
|
||||
};
|
||||
|
||||
////////////////////////////////////////////////////////
|
||||
// Conserved current - not yet implemented.
|
||||
////////////////////////////////////////////////////////
|
||||
template <class Impl>
|
||||
void ImprovedStaggeredFermion<Impl>::ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu)
|
||||
{
|
||||
assert(0);
|
||||
}
|
||||
|
||||
template <class Impl>
|
||||
void ImprovedStaggeredFermion<Impl>::SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax)
|
||||
{
|
||||
assert(0);
|
||||
}
|
||||
|
||||
FermOpStaggeredTemplateInstantiate(ImprovedStaggeredFermion);
|
||||
|
||||
//AdjointFermOpTemplateInstantiate(ImprovedStaggeredFermion);
|
||||
|
||||
@ -157,6 +157,22 @@ class ImprovedStaggeredFermion : public StaggeredKernels<Impl>, public ImprovedS
|
||||
|
||||
LebesgueOrder Lebesgue;
|
||||
LebesgueOrder LebesgueEvenOdd;
|
||||
|
||||
///////////////////////////////////////////////////////////////
|
||||
// Conserved current utilities
|
||||
///////////////////////////////////////////////////////////////
|
||||
void ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu);
|
||||
void SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax);
|
||||
};
|
||||
|
||||
typedef ImprovedStaggeredFermion<StaggeredImplF> ImprovedStaggeredFermionF;
|
||||
|
||||
@ -345,6 +345,30 @@ void ImprovedStaggeredFermion5D<Impl>::MooeeInvDag(const FermionField &in,
|
||||
MooeeInv(in, out);
|
||||
}
|
||||
|
||||
////////////////////////////////////////////////////////
|
||||
// Conserved current - not yet implemented.
|
||||
////////////////////////////////////////////////////////
|
||||
template <class Impl>
|
||||
void ImprovedStaggeredFermion5D<Impl>::ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu)
|
||||
{
|
||||
assert(0);
|
||||
}
|
||||
|
||||
template <class Impl>
|
||||
void ImprovedStaggeredFermion5D<Impl>::SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax)
|
||||
{
|
||||
assert(0);
|
||||
}
|
||||
|
||||
FermOpStaggeredTemplateInstantiate(ImprovedStaggeredFermion5D);
|
||||
FermOpStaggeredVec5dTemplateInstantiate(ImprovedStaggeredFermion5D);
|
||||
|
||||
@ -160,6 +160,21 @@ namespace QCD {
|
||||
// Comms buffer
|
||||
std::vector<SiteHalfSpinor,alignedAllocator<SiteHalfSpinor> > comm_buf;
|
||||
|
||||
///////////////////////////////////////////////////////////////
|
||||
// Conserved current utilities
|
||||
///////////////////////////////////////////////////////////////
|
||||
void ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu);
|
||||
void SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax);
|
||||
};
|
||||
|
||||
}}
|
||||
|
||||
@ -347,6 +347,53 @@ void WilsonFermion<Impl>::DhopInternal(StencilImpl &st, LebesgueOrder &lo,
|
||||
}
|
||||
};
|
||||
|
||||
/*******************************************************************************
|
||||
* Conserved current utilities for Wilson fermions, for contracting propagators
|
||||
* to make a conserved current sink or inserting the conserved current
|
||||
* sequentially.
|
||||
******************************************************************************/
|
||||
template <class Impl>
|
||||
void WilsonFermion<Impl>::ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu)
|
||||
{
|
||||
Gamma g5(Gamma::Algebra::Gamma5);
|
||||
conformable(_grid, q_in_1._grid);
|
||||
conformable(_grid, q_in_2._grid);
|
||||
conformable(_grid, q_out._grid);
|
||||
Kernels::ContractConservedCurrentInternal(q_in_1, q_in_2, q_out,
|
||||
Umu, curr_type, mu);
|
||||
}
|
||||
|
||||
template <class Impl>
|
||||
void WilsonFermion<Impl>::SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax)
|
||||
{
|
||||
conformable(_grid, q_in._grid);
|
||||
conformable(_grid, q_out._grid);
|
||||
Lattice<iSinglet<Simd>> ph(_grid), coor(_grid);
|
||||
Complex i(0.0,1.0);
|
||||
|
||||
// Momentum projection
|
||||
ph = zero;
|
||||
for(unsigned int mu = 0; mu < Nd - 1; mu++)
|
||||
{
|
||||
LatticeCoordinate(coor, mu);
|
||||
ph = ph + mom[mu]*coor*((1./(_grid->_fdimensions[mu])));
|
||||
}
|
||||
ph = exp((Real)(2*M_PI)*i*ph);
|
||||
|
||||
Kernels::SeqConservedCurrentInternal(q_in, q_out, Umu, curr_type, mu, ph,
|
||||
tmin, tmax);
|
||||
}
|
||||
|
||||
FermOpTemplateInstantiate(WilsonFermion);
|
||||
AdjointFermOpTemplateInstantiate(WilsonFermion);
|
||||
TwoIndexFermOpTemplateInstantiate(WilsonFermion);
|
||||
|
||||
@ -146,6 +146,22 @@ class WilsonFermion : public WilsonKernels<Impl>, public WilsonFermionStatic {
|
||||
|
||||
LebesgueOrder Lebesgue;
|
||||
LebesgueOrder LebesgueEvenOdd;
|
||||
|
||||
///////////////////////////////////////////////////////////////
|
||||
// Conserved current utilities
|
||||
///////////////////////////////////////////////////////////////
|
||||
void ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu);
|
||||
void SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax);
|
||||
};
|
||||
|
||||
typedef WilsonFermion<WilsonImplF> WilsonFermionF;
|
||||
|
||||
@ -679,6 +679,83 @@ void WilsonFermion5D<Impl>::MomentumSpacePropagatorHw(FermionField &out,const Fe
|
||||
|
||||
}
|
||||
|
||||
/*******************************************************************************
|
||||
* Conserved current utilities for Wilson fermions, for contracting propagators
|
||||
* to make a conserved current sink or inserting the conserved current
|
||||
* sequentially.
|
||||
******************************************************************************/
|
||||
template <class Impl>
|
||||
void WilsonFermion5D<Impl>::ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu)
|
||||
{
|
||||
conformable(q_in_1._grid, FermionGrid());
|
||||
conformable(q_in_1._grid, q_in_2._grid);
|
||||
conformable(_FourDimGrid, q_out._grid);
|
||||
|
||||
PropagatorField q1_s(_FourDimGrid);
|
||||
PropagatorField q2_s(_FourDimGrid);
|
||||
PropagatorField tmp(_FourDimGrid);
|
||||
|
||||
// Contract across 5th dimension.
|
||||
q_out = zero;
|
||||
for (int s = 0; s < Ls; ++s)
|
||||
{
|
||||
ExtractSlice(q1_s, q_in_1, 0, s);
|
||||
ExtractSlice(q2_s, q_in_2, 0, Ls - s - 1);
|
||||
Kernels::ContractConservedCurrentInternal(q1_s, q2_s, tmp, Umu, curr_type, mu);
|
||||
|
||||
// Axial current sign
|
||||
Real G_s = (curr_type == Current::Axial) ? ((s < Ls/2) ? -1. : 1.) : 1.;
|
||||
q_out += G_s*tmp;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
template <class Impl>
|
||||
void WilsonFermion5D<Impl>::SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax)
|
||||
{
|
||||
conformable(q_in._grid, FermionGrid());
|
||||
conformable(q_in._grid, q_out._grid);
|
||||
Lattice<iSinglet<Simd>> ph(_FourDimGrid), coor(_FourDimGrid);
|
||||
Complex i(0.0, 1.0);
|
||||
|
||||
// Momentum projection
|
||||
ph = zero;
|
||||
for(unsigned int nu = 0; nu < Nd - 1; nu++)
|
||||
{
|
||||
LatticeCoordinate(coor, nu);
|
||||
ph = ph + mom[nu]*coor*((1./(_FourDimGrid->_fdimensions[nu])));
|
||||
}
|
||||
ph = exp((Real)(2*M_PI)*i*ph);
|
||||
|
||||
// Sequential insertion
|
||||
Kernels::SeqConservedCurrentInternal(q_in, q_out, Umu, curr_type,
|
||||
mu, ph, tmin, tmax);
|
||||
|
||||
// Axial current sign.
|
||||
if (curr_type == Current::Axial)
|
||||
{
|
||||
SitePropagator result;
|
||||
parallel_for(int sU = 0; sU < Umu._grid->oSites(); sU++)
|
||||
{
|
||||
int sF = sU * Ls;
|
||||
for (int s = 0; s < Ls/2; s++)
|
||||
{
|
||||
vstream(q_out._odata[sF], -q_out._odata[sF]);
|
||||
sF++;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
FermOpTemplateInstantiate(WilsonFermion5D);
|
||||
GparityFermOpTemplateInstantiate(WilsonFermion5D);
|
||||
|
||||
@ -214,6 +214,21 @@ namespace QCD {
|
||||
// Comms buffer
|
||||
std::vector<SiteHalfSpinor,alignedAllocator<SiteHalfSpinor> > comm_buf;
|
||||
|
||||
///////////////////////////////////////////////////////////////
|
||||
// Conserved current utilities
|
||||
///////////////////////////////////////////////////////////////
|
||||
void ContractConservedCurrent(PropagatorField &q_in_1,
|
||||
PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu);
|
||||
void SeqConservedCurrent(PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
std::vector<Real> mom,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax);
|
||||
};
|
||||
|
||||
}}
|
||||
|
||||
@ -567,6 +567,95 @@ void WilsonKernels<Impl>::DhopDir( StencilImpl &st, DoubledGaugeField &U,SiteHal
|
||||
vstream(out._odata[sF], result);
|
||||
}
|
||||
|
||||
/*******************************************************************************
|
||||
* Conserved current utilities for Wilson fermions, for contracting propagators
|
||||
* to make a conserved current sink or inserting the conserved current
|
||||
* sequentially. Common to both 4D and 5D.
|
||||
******************************************************************************/
|
||||
#define WilsonCurrentFwd(expr, mu) (0.5*(Gamma::gmu[mu]*expr - expr))
|
||||
#define WilsonCurrentBwd(expr, mu) (0.5*(Gamma::gmu[mu]*expr + expr))
|
||||
|
||||
template<class Impl>
|
||||
void WilsonKernels<Impl>::ContractConservedCurrentInternal(const PropagatorField &q_in_1,
|
||||
const PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
DoubledGaugeField &U,
|
||||
Current curr_type,
|
||||
unsigned int mu)
|
||||
{
|
||||
Gamma g5(Gamma::Algebra::Gamma5);
|
||||
PropagatorField tmp(q_out._grid);
|
||||
GaugeLinkField Umu(U._grid);
|
||||
Umu = PeekIndex<LorentzIndex>(U, mu);
|
||||
|
||||
tmp = this->CovShiftForward(Umu, mu, q_in_1);
|
||||
q_out = (g5*adj(q_in_2)*g5)*WilsonCurrentFwd(tmp, mu);
|
||||
|
||||
tmp = adj(Umu)*q_in_1;
|
||||
q_out += (g5*adj(this->CovShiftForward(Umu, mu, q_in_2))*g5)*WilsonCurrentBwd(q_in_1, mu);
|
||||
}
|
||||
|
||||
|
||||
template <class Impl>
|
||||
void WilsonKernels<Impl>::SeqConservedCurrentInternal(const PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
DoubledGaugeField &U,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
Lattice<iSinglet<Simd>> &ph,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax)
|
||||
{
|
||||
int tshift = (mu == Nd - 1) ? 1 : 0;
|
||||
Real G_T = (curr_type == Current::Tadpole) ? -1. : 1.;
|
||||
PropagatorField tmp(q_in._grid);
|
||||
GaugeLinkField Umu(U._grid);
|
||||
Umu = PeekIndex<LorentzIndex>(U, mu);
|
||||
Lattice<iScalar<vInteger>> t(q_in._grid);
|
||||
|
||||
tmp = this->CovShiftForward(Umu, mu, q_in)*ph;
|
||||
where((t >= tmin) and (t <= tmax), tmp, 0.*tmp);
|
||||
q_out = G_T*WilsonCurrentFwd(tmp, mu);
|
||||
|
||||
tmp = q_in*ph;
|
||||
tmp = this->CovShiftBackward(Umu, mu, tmp);
|
||||
where((t >= tmin + tshift) and (t <= tmax + tshift), tmp, 0.*tmp);
|
||||
q_out += WilsonCurrentBwd(tmp, mu);
|
||||
}
|
||||
|
||||
|
||||
// GParity, (Z)DomainWallVec5D -> require special implementation
|
||||
#define NO_CURR(Impl) \
|
||||
template <> void \
|
||||
WilsonKernels<Impl>::ContractConservedCurrentInternal(const PropagatorField &q_in_1, \
|
||||
const PropagatorField &q_in_2, \
|
||||
PropagatorField &q_out, \
|
||||
DoubledGaugeField &U, \
|
||||
Current curr_type, \
|
||||
unsigned int mu) \
|
||||
{ \
|
||||
assert(0); \
|
||||
} \
|
||||
template <> void \
|
||||
WilsonKernels<Impl>::SeqConservedCurrentInternal(const PropagatorField &q_in, \
|
||||
PropagatorField &q_out, \
|
||||
DoubledGaugeField &U, \
|
||||
Current curr_type, \
|
||||
unsigned int mu, \
|
||||
Lattice<iSinglet<Simd>> &ph, \
|
||||
unsigned int tmin, \
|
||||
unsigned int tmax) \
|
||||
{ \
|
||||
assert(0); \
|
||||
}
|
||||
|
||||
NO_CURR(GparityWilsonImplF);
|
||||
NO_CURR(GparityWilsonImplD);
|
||||
NO_CURR(DomainWallVec5dImplF);
|
||||
NO_CURR(DomainWallVec5dImplD);
|
||||
NO_CURR(ZDomainWallVec5dImplF);
|
||||
NO_CURR(ZDomainWallVec5dImplD);
|
||||
|
||||
FermOpTemplateInstantiate(WilsonKernels);
|
||||
AdjointFermOpTemplateInstantiate(WilsonKernels);
|
||||
TwoIndexFermOpTemplateInstantiate(WilsonKernels);
|
||||
|
||||
@ -166,6 +166,24 @@ public:
|
||||
void DhopDir(StencilImpl &st, DoubledGaugeField &U,SiteHalfSpinor * buf,
|
||||
int sF, int sU, const FermionField &in, FermionField &out, int dirdisp, int gamma);
|
||||
|
||||
//////////////////////////////////////////////////////////////////////////////
|
||||
// Utilities for inserting Wilson conserved current.
|
||||
//////////////////////////////////////////////////////////////////////////////
|
||||
void ContractConservedCurrentInternal(const PropagatorField &q_in_1,
|
||||
const PropagatorField &q_in_2,
|
||||
PropagatorField &q_out,
|
||||
DoubledGaugeField &U,
|
||||
Current curr_type,
|
||||
unsigned int mu);
|
||||
void SeqConservedCurrentInternal(const PropagatorField &q_in,
|
||||
PropagatorField &q_out,
|
||||
DoubledGaugeField &U,
|
||||
Current curr_type,
|
||||
unsigned int mu,
|
||||
Lattice<iSinglet<Simd>> &ph,
|
||||
unsigned int tmin,
|
||||
unsigned int tmax);
|
||||
|
||||
private:
|
||||
// Specialised variants
|
||||
void GenericDhopSite(StencilImpl &st, LebesgueOrder &lo, DoubledGaugeField &U, SiteHalfSpinor * buf,
|
||||
|
||||
Loading…
x
Reference in New Issue
Block a user