portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2025-04-09 21:50:45 +01:00
Code Releases Activity

Namespace

Browse Source
This commit is contained in:
paboyle 2018-01-14 22:52:22 +00:00
parent 7828887604
commit c2b0e0269a
1 changed files with 138 additions and 138 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

276
lib/qcd/action/gauge/Photon.h
View File

@ -23,93 +23,93 @@
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
See the full license in the file "LICENSE" in the top level distribution directory See the full license in the file "LICENSE" in the top level distribution directory
*************************************************************************************/ *************************************************************************************/
/* END LEGAL */ /* END LEGAL */
#ifndef QCD_PHOTON_ACTION_H #ifndef QCD_PHOTON_ACTION_H
#define QCD_PHOTON_ACTION_H #define QCD_PHOTON_ACTION_H
namespace Grid{ NAMESPACE_BEGIN(Grid);
namespace QCD{
template <class S>
class QedGimpl
{
public:
typedef S Simd;
template <typename vtype>
using iImplGaugeLink = iScalar<iScalar<iScalar<vtype>>>;
template <typename vtype>
using iImplGaugeField = iVector<iScalar<iScalar<vtype>>, Nd>;
typedef iImplGaugeLink<Simd> SiteLink;
typedef iImplGaugeField<Simd> SiteField;
typedef SiteField SiteComplex;
typedef Lattice<SiteLink> LinkField;
typedef Lattice<SiteField> Field;
typedef Field ComplexField;
};
typedef QedGimpl<vComplex> QedGimplR;
template<class Gimpl>
class Photon
{
public:
INHERIT_GIMPL_TYPES(Gimpl);
GRID_SERIALIZABLE_ENUM(Gauge, undef, feynman, 1, coulomb, 2, landau, 3);
GRID_SERIALIZABLE_ENUM(ZmScheme, undef, qedL, 1, qedTL, 2);
public:
Photon(Gauge gauge, ZmScheme zmScheme);
virtual ~Photon(void) = default;
void FreePropagator(const GaugeField &in, GaugeField &out);
void MomentumSpacePropagator(const GaugeField &in, GaugeField &out);
void StochasticWeight(GaugeLinkField &weight);
void StochasticField(GaugeField &out, GridParallelRNG &rng);
void StochasticField(GaugeField &out, GridParallelRNG &rng,
const GaugeLinkField &weight);
private:
void invKHatSquared(GaugeLinkField &out);
void zmSub(GaugeLinkField &out);
private:
Gauge gauge_;
ZmScheme zmScheme_;
};
typedef Photon<QedGimplR> PhotonR; template <class S>
class QedGimpl
{
public:
typedef S Simd;
template <typename vtype>
using iImplGaugeLink = iScalar<iScalar<iScalar<vtype>>>;
template <typename vtype>
using iImplGaugeField = iVector<iScalar<iScalar<vtype>>, Nd>;
typedef iImplGaugeLink<Simd> SiteLink;
typedef iImplGaugeField<Simd> SiteField;
typedef SiteField SiteComplex;
typedef Lattice<SiteLink> LinkField;
typedef Lattice<SiteField> Field;
typedef Field ComplexField;
};
template<class Gimpl> typedef QedGimpl<vComplex> QedGimplR;
Photon<Gimpl>::Photon(Gauge gauge, ZmScheme zmScheme)
template<class Gimpl>
class Photon
{
public:
INHERIT_GIMPL_TYPES(Gimpl);
GRID_SERIALIZABLE_ENUM(Gauge, undef, feynman, 1, coulomb, 2, landau, 3);
GRID_SERIALIZABLE_ENUM(ZmScheme, undef, qedL, 1, qedTL, 2);
public:
Photon(Gauge gauge, ZmScheme zmScheme);
virtual ~Photon(void) = default;
void FreePropagator(const GaugeField &in, GaugeField &out);
void MomentumSpacePropagator(const GaugeField &in, GaugeField &out);
void StochasticWeight(GaugeLinkField &weight);
void StochasticField(GaugeField &out, GridParallelRNG &rng);
void StochasticField(GaugeField &out, GridParallelRNG &rng,
const GaugeLinkField &weight);
private:
void invKHatSquared(GaugeLinkField &out);
void zmSub(GaugeLinkField &out);
private:
Gauge gauge_;
ZmScheme zmScheme_;
};
typedef Photon<QedGimplR> PhotonR;
template<class Gimpl>
Photon<Gimpl>::Photon(Gauge gauge, ZmScheme zmScheme)
: gauge_(gauge), zmScheme_(zmScheme) : gauge_(gauge), zmScheme_(zmScheme)
{} {}
template<class Gimpl> template<class Gimpl>
void Photon<Gimpl>::FreePropagator (const GaugeField &in,GaugeField &out) void Photon<Gimpl>::FreePropagator (const GaugeField &in,GaugeField &out)
{ {
FFT theFFT(in._grid); FFT theFFT(in._grid);
GaugeField in_k(in._grid); GaugeField in_k(in._grid);
GaugeField prop_k(in._grid); GaugeField prop_k(in._grid);
theFFT.FFT_all_dim(in_k,in,FFT::forward); theFFT.FFT_all_dim(in_k,in,FFT::forward);
MomentumSpacePropagator(prop_k,in_k); MomentumSpacePropagator(prop_k,in_k);
theFFT.FFT_all_dim(out,prop_k,FFT::backward); theFFT.FFT_all_dim(out,prop_k,FFT::backward);
} }
template<class Gimpl> template<class Gimpl>
void Photon<Gimpl>::invKHatSquared(GaugeLinkField &out) void Photon<Gimpl>::invKHatSquared(GaugeLinkField &out)
{ {
GridBase *grid = out._grid; GridBase *grid = out._grid;
GaugeLinkField kmu(grid), one(grid); GaugeLinkField kmu(grid), one(grid);
const unsigned int nd = grid->_ndimension; const unsigned int nd = grid->_ndimension;
std::vector<int> &l = grid->_fdimensions; std::vector<int> &l = grid->_fdimensions;
std::vector<int> zm(nd,0); std::vector<int> zm(nd,0);
TComplex Tone = Complex(1.0,0.0); TComplex Tone = Complex(1.0,0.0);
TComplex Tzero= Complex(0.0,0.0); TComplex Tzero= Complex(0.0,0.0);
one = Complex(1.0,0.0); one = Complex(1.0,0.0);
out = zero; out = zero;
for(int mu = 0; mu < nd; mu++) for(int mu = 0; mu < nd; mu++)
{ {
Real twoPiL = M_PI*2./l[mu]; Real twoPiL = M_PI*2./l[mu];
@ -117,20 +117,20 @@ namespace QCD{
kmu = 2.*sin(.5*twoPiL*kmu); kmu = 2.*sin(.5*twoPiL*kmu);
out = out + kmu*kmu; out = out + kmu*kmu;
} }
pokeSite(Tone, out, zm); pokeSite(Tone, out, zm);
out = one/out; out = one/out;
pokeSite(Tzero, out, zm); pokeSite(Tzero, out, zm);
} }
template<class Gimpl> template<class Gimpl>
void Photon<Gimpl>::zmSub(GaugeLinkField &out) void Photon<Gimpl>::zmSub(GaugeLinkField &out)
{ {
GridBase *grid = out._grid; GridBase *grid = out._grid;
const unsigned int nd = grid->_ndimension; const unsigned int nd = grid->_ndimension;
switch (zmScheme_) switch (zmScheme_)
{ {
case ZmScheme::qedTL: case ZmScheme::qedTL:
{ {
std::vector<int> zm(nd,0); std::vector<int> zm(nd,0);
TComplex Tzero = Complex(0.0,0.0); TComplex Tzero = Complex(0.0,0.0);
@ -139,86 +139,86 @@ namespace QCD{
break; break;
} }
case ZmScheme::qedL: case ZmScheme::qedL:
{ {
LatticeInteger spNrm(grid), coor(grid); LatticeInteger spNrm(grid), coor(grid);
GaugeLinkField z(grid); GaugeLinkField z(grid);
spNrm = zero; spNrm = zero;
for(int d = 0; d < grid->_ndimension - 1; d++) for(int d = 0; d < grid->_ndimension - 1; d++)
{ {
LatticeCoordinate(coor,d); LatticeCoordinate(coor,d);
spNrm = spNrm + coor*coor; spNrm = spNrm + coor*coor;
} }
out = where(spNrm == Integer(0), 0.*out, out); out = where(spNrm == Integer(0), 0.*out, out);
break; break;
} }
default: default:
break; break;
} }
} }
template<class Gimpl> template<class Gimpl>
void Photon<Gimpl>::MomentumSpacePropagator(const GaugeField &in, void Photon<Gimpl>::MomentumSpacePropagator(const GaugeField &in,
GaugeField &out) GaugeField &out)
{ {
GridBase *grid = out._grid; GridBase *grid = out._grid;
LatticeComplex k2Inv(grid); LatticeComplex k2Inv(grid);
invKHatSquared(k2Inv); invKHatSquared(k2Inv);
zmSub(k2Inv); zmSub(k2Inv);
out = in*k2Inv; out = in*k2Inv;
} }
template<class Gimpl> template<class Gimpl>
void Photon<Gimpl>::StochasticWeight(GaugeLinkField &weight) void Photon<Gimpl>::StochasticWeight(GaugeLinkField &weight)
{ {
auto *grid = dynamic_cast<GridCartesian *>(weight._grid); auto *grid = dynamic_cast<GridCartesian *>(weight._grid);
const unsigned int nd = grid->_ndimension; const unsigned int nd = grid->_ndimension;
std::vector<int> latt_size = grid->_fdimensions; std::vector<int> latt_size = grid->_fdimensions;
Integer vol = 1; Integer vol = 1;
for(int d = 0; d < nd; d++) for(int d = 0; d < nd; d++)
{ {
vol = vol * latt_size[d]; vol = vol * latt_size[d];
} }
invKHatSquared(weight); invKHatSquared(weight);
weight = sqrt(vol*real(weight)); weight = sqrt(vol*real(weight));
zmSub(weight); zmSub(weight);
} }
template<class Gimpl> template<class Gimpl>
void Photon<Gimpl>::StochasticField(GaugeField &out, GridParallelRNG &rng) void Photon<Gimpl>::StochasticField(GaugeField &out, GridParallelRNG &rng)
{ {
auto *grid = dynamic_cast<GridCartesian *>(out._grid); auto *grid = dynamic_cast<GridCartesian *>(out._grid);
GaugeLinkField weight(grid); GaugeLinkField weight(grid);
StochasticWeight(weight); StochasticWeight(weight);
StochasticField(out, rng, weight); StochasticField(out, rng, weight);
} }
template<class Gimpl> template<class Gimpl>
void Photon<Gimpl>::StochasticField(GaugeField &out, GridParallelRNG &rng, void Photon<Gimpl>::StochasticField(GaugeField &out, GridParallelRNG &rng,
const GaugeLinkField &weight) const GaugeLinkField &weight)
{ {
auto *grid = dynamic_cast<GridCartesian *>(out._grid); auto *grid = dynamic_cast<GridCartesian *>(out._grid);
const unsigned int nd = grid->_ndimension; const unsigned int nd = grid->_ndimension;
GaugeLinkField r(grid); GaugeLinkField r(grid);
GaugeField aTilde(grid); GaugeField aTilde(grid);
FFT fft(grid); FFT fft(grid);
for(int mu = 0; mu < nd; mu++) for(int mu = 0; mu < nd; mu++)
{ {
gaussian(rng, r); gaussian(rng, r);
r = weight*r; r = weight*r;
pokeLorentz(aTilde, r, mu); pokeLorentz(aTilde, r, mu);
} }
fft.FFT_all_dim(out, aTilde, FFT::backward); fft.FFT_all_dim(out, aTilde, FFT::backward);
out = real(out); out = real(out);
} }
// template<class Gimpl> // template<class Gimpl>
// void Photon<Gimpl>::FeynmanGaugeMomentumSpacePropagator_L(GaugeField &out, // void Photon<Gimpl>::FeynmanGaugeMomentumSpacePropagator_L(GaugeField &out,
// const GaugeField &in) // const GaugeField &in)
@ -282,5 +282,5 @@ namespace QCD{
// out = in*denom; // out = in*denom;
// }; // };
}} NAMESPACE_END(Grid);
#endif #endif