portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2024-09-19 16:55:37 +01:00
Code Releases Activity

Qslash term added

Browse Source
This commit is contained in:
Peter Boyle 2023-09-14 16:14:03 -04:00
parent ad2b699d2b
commit e29b97b3ea
4 changed files with 74 additions and 5 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

5
Grid/qcd/action/fermion/CayleyFermion5D.h
View File

@ -124,6 +124,11 @@ public:
RealD _b;
RealD _c;
// possible boost
std::vector<ComplexD> qmu;
void set_qmu(std::vector<ComplexD> _qmu) { qmu=_qmu; assert(qmu.size()==Nd);};
void addQmu(const FermionField &in, FermionField &out, int dag);
// Cayley form Moebius (tanh and zolotarev)
Vector<Coeff_t> omega;
Vector<Coeff_t> bs; // S dependent coeffs

40
Grid/qcd/action/fermion/implementation/CayleyFermion5DImplementation.h
View File

@ -48,7 +48,8 @@ CayleyFermion5D<Impl>::CayleyFermion5D(GaugeField &_Umu,
FourDimGrid,
FourDimRedBlackGrid,_M5,p),
mass_plus(_mass), mass_minus(_mass)
{
{
// qmu defaults to zero size;
}
///////////////////////////////////////////////////////////////
@ -270,6 +271,34 @@ void CayleyFermion5D<Impl>::MeooeDag5D (const FermionField &psi, FermionField
M5Ddag(psi,psi,Din,lower,diag,upper);
}
template<class Impl>
void CayleyFermion5D<Impl>::addQmu(const FermionField &psi,FermionField &chi, int dag)
{
if ( qmu.size() ) {
Gamma::Algebra Gmu [] = {
Gamma::Algebra::GammaX,
Gamma::Algebra::GammaY,
Gamma::Algebra::GammaZ,
Gamma::Algebra::GammaT
};
std::vector<ComplexD> coeff(Nd);
ComplexD ci(0,1);
assert(qmu.size()==Nd);
for(int mu=0;mu<Nd;mu++){
coeff[mu] = ci*qmu[mu];
if ( dag ) coeff[mu] = conjugate(coeff[mu]);
}
chi = chi + Gamma(Gmu[0])*psi*coeff[0];
for(int mu=1;mu<Nd;mu++){
chi = chi + Gamma(Gmu[mu])*psi*coeff[mu];
}
}
}
template<class Impl>
void CayleyFermion5D<Impl>::M (const FermionField &psi, FermionField &chi)
{
@ -277,8 +306,12 @@ void CayleyFermion5D<Impl>::M (const FermionField &psi, FermionField &chi)
// Assemble Din
Meooe5D(psi,Din);
this->DW(Din,chi,DaggerNo);
// add i q_mu gamma_mu here
addQmu(Din,chi,DaggerNo);
// ((b D_W + D_w hop terms +1) on s-diag
axpby(chi,1.0,1.0,chi,psi);
@ -295,6 +328,9 @@ void CayleyFermion5D<Impl>::Mdag (const FermionField &psi, FermionField &chi)
FermionField Din(psi.Grid());
// Apply Dw
this->DW(psi,Din,DaggerYes);
// add -i conj(q_mu) gamma_mu here ... if qmu is real, gammm_5 hermitian, otherwise not.
addQmu(psi,Din,DaggerYes);
MeooeDag5D(Din,chi);

8
tests/solver/Test_dwf_cg_prec.cc
View File

@ -1,4 +1,4 @@
/*************************************************************************************
*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid
@ -67,7 +67,13 @@ int main(int argc, char** argv) {
result = Zero();
LatticeGaugeField Umu(UGrid);
#if 0
FieldMetaData header;
std::string file("ckpoint_lat.4000");
NerscIO::readConfiguration(Umu,header,file);
#else
SU<Nc>::HotConfiguration(RNG4, Umu);
#endif
std::cout << GridLogMessage << "Lattice dimensions: " << GridDefaultLatt()
<< " Ls: " << Ls << std::endl;

26
tests/solver/Test_dwf_cg_unprec.cc
View File

@ -54,15 +54,30 @@ int main (int argc, char ** argv)
GridCartesian * FGrid = SpaceTimeGrid::makeFiveDimGrid(Ls,UGrid);
GridRedBlackCartesian * FrbGrid = SpaceTimeGrid::makeFiveDimRedBlackGrid(Ls,UGrid);
std::vector<ComplexD> qmu;
qmu.push_back(ComplexD(0.1,0.0));
qmu.push_back(ComplexD(0.0,0.0));
qmu.push_back(ComplexD(0.0,0.0));
qmu.push_back(ComplexD(0.0,0.01));
std::vector<int> seeds4({1,2,3,4});
std::vector<int> seeds5({5,6,7,8});
GridParallelRNG RNG5(FGrid); RNG5.SeedFixedIntegers(seeds5);
GridParallelRNG RNG4(UGrid); RNG4.SeedFixedIntegers(seeds4);
LatticeFermion tmp(FGrid);
LatticeFermion src(FGrid); random(RNG5,src);
LatticeFermion result(FGrid); result=Zero();
LatticeGaugeField Umu(UGrid); SU<Nc>::HotConfiguration(RNG4,Umu);
LatticeGaugeField Umu(UGrid);
#if 0
FieldMetaData header;
std::string file("ckpoint_lat.4000");
NerscIO::readConfiguration(Umu,header,file);
#else
SU<Nc>::HotConfiguration(RNG4,Umu);
#endif
std::vector<LatticeColourMatrix> U(4,UGrid);
for(int mu=0;mu<Nd;mu++){
U[mu] = PeekIndex<LorentzIndex>(Umu,mu);
@ -71,8 +86,15 @@ int main (int argc, char ** argv)
RealD mass=0.1;
RealD M5=1.8;
DomainWallFermionD Ddwf(Umu,*FGrid,*FrbGrid,*UGrid,*UrbGrid,mass,M5);
Ddwf.qmu = qmu;
Ddwf.M(src,tmp);
std::cout << " |M src|^2 "<<norm2(tmp)<<std::endl;
MdagMLinearOperator<DomainWallFermionD,LatticeFermion> HermOp(Ddwf);
HermOp.HermOp(src,tmp);
std::cout << " <src|MdagM| src> "<<innerProduct(src,tmp)<<std::endl;
ConjugateGradient<LatticeFermion> CG(1.0e-6,10000);
CG(HermOp,src,result);