1
0
mirror of https://github.com/paboyle/Grid.git synced 2024-11-14 17:55:38 +00:00
Grid/lib/qcd/utils/LinalgUtils.h

137 lines
3.7 KiB
C
Raw Normal View History

#ifndef GRID_QCD_LINALG_UTILS_H
#define GRID_QCD_LINALG_UTILS_H
namespace Grid{
namespace QCD{
////////////////////////////////////////////////////////////////////////
//This file brings additional linear combination assist that is helpful
//to QCD such as chiral projectors and spin matrices applied to one of the inputs.
//These routines support five-D chiral fermions and contain s-subslice indexing
//on the 5d (rb4d) checkerboarded lattices
////////////////////////////////////////////////////////////////////////
template<class vobj>
void axpby_ssp(Lattice<vobj> &z, RealD a,const Lattice<vobj> &x,RealD b,const Lattice<vobj> &y,int s,int sp)
{
z.checkerboard = x.checkerboard;
conformable(x,y);
conformable(x,z);
GridBase *grid=x._grid;
int Ls = grid->_rdimensions[0];
PARALLEL_FOR_LOOP
for(int ss=0;ss<grid->oSites();ss+=Ls){ // adds Ls
vobj tmp = a*x._odata[ss+s]+b*y._odata[ss+sp];
vstream(z._odata[ss+s],tmp);
}
}
template<class vobj>
void ag5xpby_ssp(Lattice<vobj> &z,RealD a,const Lattice<vobj> &x,RealD b,const Lattice<vobj> &y,int s,int sp)
{
z.checkerboard = x.checkerboard;
conformable(x,y);
conformable(x,z);
GridBase *grid=x._grid;
int Ls = grid->_rdimensions[0];
2015-08-13 10:51:29 +01:00
Gamma G5(Gamma::Gamma5);
PARALLEL_FOR_LOOP
for(int ss=0;ss<grid->oSites();ss+=Ls){ // adds Ls
vobj tmp;
2015-08-13 10:51:29 +01:00
tmp = G5*x._odata[ss+s]*a;
tmp = tmp + b*y._odata[ss+sp];
vstream(z._odata[ss+s],tmp);
}
}
template<class vobj>
void axpbg5y_ssp(Lattice<vobj> &z,RealD a,const Lattice<vobj> &x,RealD b,const Lattice<vobj> &y,int s,int sp)
{
z.checkerboard = x.checkerboard;
conformable(x,y);
conformable(x,z);
GridBase *grid=x._grid;
int Ls = grid->_rdimensions[0];
2015-08-13 10:51:29 +01:00
Gamma G5(Gamma::Gamma5);
PARALLEL_FOR_LOOP
for(int ss=0;ss<grid->oSites();ss+=Ls){ // adds Ls
vobj tmp;
2015-08-13 10:51:29 +01:00
tmp = G5*y._odata[ss+sp]*b;
tmp = tmp + a*x._odata[ss+s];
vstream(z._odata[ss+s],tmp);
}
}
template<class vobj>
void ag5xpbg5y_ssp(Lattice<vobj> &z,RealD a,const Lattice<vobj> &x,RealD b,const Lattice<vobj> &y,int s,int sp)
{
z.checkerboard = x.checkerboard;
conformable(x,y);
conformable(x,z);
GridBase *grid=x._grid;
int Ls = grid->_rdimensions[0];
2015-08-13 10:51:29 +01:00
Gamma G5(Gamma::Gamma5);
PARALLEL_FOR_LOOP
for(int ss=0;ss<grid->oSites();ss+=Ls){ // adds Ls
vobj tmp1;
vobj tmp2;
tmp1 = a*x._odata[ss+s]+b*y._odata[ss+sp];
2015-08-13 10:51:29 +01:00
tmp2 = G5*tmp1;
vstream(z._odata[ss+s],tmp2);
}
}
template<class vobj>
void axpby_ssp_pminus(Lattice<vobj> &z,RealD a,const Lattice<vobj> &x,RealD b,const Lattice<vobj> &y,int s,int sp)
{
z.checkerboard = x.checkerboard;
conformable(x,y);
conformable(x,z);
GridBase *grid=x._grid;
int Ls = grid->_rdimensions[0];
PARALLEL_FOR_LOOP
for(int ss=0;ss<grid->oSites();ss+=Ls){ // adds Ls
vobj tmp;
spProj5m(tmp,y._odata[ss+sp]);
tmp = a*x._odata[ss+s]+b*tmp;
vstream(z._odata[ss+s],tmp);
}
}
template<class vobj>
void axpby_ssp_pplus(Lattice<vobj> &z,RealD a,const Lattice<vobj> &x,RealD b,const Lattice<vobj> &y,int s,int sp)
{
z.checkerboard = x.checkerboard;
conformable(x,y);
conformable(x,z);
GridBase *grid=x._grid;
int Ls = grid->_rdimensions[0];
PARALLEL_FOR_LOOP
for(int ss=0;ss<grid->oSites();ss+=Ls){ // adds Ls
vobj tmp;
spProj5p(tmp,y._odata[ss+sp]);
tmp = a*x._odata[ss+s]+b*tmp;
vstream(z._odata[ss+s],tmp);
}
}
2015-06-05 18:16:25 +01:00
template<class vobj>
void G5R5(Lattice<vobj> &z,const Lattice<vobj> &x)
{
GridBase *grid=x._grid;
z.checkerboard = x.checkerboard;
conformable(x,z);
int Ls = grid->_rdimensions[0];
2015-08-13 10:51:29 +01:00
Gamma G5(Gamma::Gamma5);
2015-06-05 18:16:25 +01:00
PARALLEL_FOR_LOOP
for(int ss=0;ss<grid->oSites();ss+=Ls){ // adds Ls
vobj tmp;
for(int s=0;s<Ls;s++){
int sp = Ls-1-s;
2015-08-13 10:51:29 +01:00
tmp = G5*x._odata[ss+s];
2015-06-05 18:16:25 +01:00
vstream(z._odata[ss+sp],tmp);
}
}
}
}}
#endif