mirror of
https://github.com/paboyle/Grid.git
synced 2024-11-10 15:55:37 +00:00
QedFVol: calculate plaquette and 2x2 Wilson loop of stochastic QED field
This commit is contained in:
James Harrison
parent
26d124283e
commit
3ab4c8c0bb
167
programs/qed-fvol/WilsonLoops.h
Normal file
167
programs/qed-fvol/WilsonLoops.h
Normal file
@ -0,0 +1,167 @@
|
||||
#ifndef QEDFVOL_WILSONLOOPS_H
|
||||
#define QEDFVOL_WILSONLOOPS_H
|
||||
|
||||
#include <Global.hpp>
|
||||
|
||||
BEGIN_QEDFVOL_NAMESPACE
|
||||
|
||||
template <class Gimpl> class WilsonLoops : public Gimpl {
|
||||
public:
|
||||
INHERIT_GIMPL_TYPES(Gimpl);
|
||||
|
||||
typedef typename Gimpl::GaugeLinkField GaugeMat;
|
||||
typedef typename Gimpl::GaugeField GaugeLorentz;
|
||||
|
||||
//////////////////////////////////////////////////
|
||||
// directed plaquette oriented in mu,nu plane
|
||||
//////////////////////////////////////////////////
|
||||
static void dirPlaquette(GaugeMat &plaq, const std::vector<GaugeMat> &U,
|
||||
const int mu, const int nu) {
|
||||
// Annoyingly, must use either scope resolution to find dependent base
|
||||
// class,
|
||||
// or this-> ; there is no "this" in a static method. This forces explicit
|
||||
// Gimpl scope
|
||||
// resolution throughout the usage in this file, and rather defeats the
|
||||
// purpose of deriving
|
||||
// from Gimpl.
|
||||
plaq = Gimpl::CovShiftBackward(
|
||||
U[mu], mu, Gimpl::CovShiftBackward(
|
||||
U[nu], nu, Gimpl::CovShiftForward(U[mu], mu, U[nu])));
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// trace of directed plaquette oriented in mu,nu plane
|
||||
//////////////////////////////////////////////////
|
||||
static void traceDirPlaquette(LatticeComplex &plaq,
|
||||
const std::vector<GaugeMat> &U, const int mu,
|
||||
const int nu) {
|
||||
GaugeMat sp(U[0]._grid);
|
||||
dirPlaquette(sp, U, mu, nu);
|
||||
plaq = trace(sp);
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// sum over all planes of plaquette
|
||||
//////////////////////////////////////////////////
|
||||
static void sitePlaquette(LatticeComplex &Plaq,
|
||||
const std::vector<GaugeMat> &U) {
|
||||
LatticeComplex sitePlaq(U[0]._grid);
|
||||
Plaq = zero;
|
||||
for (int mu = 1; mu < Nd; mu++) {
|
||||
for (int nu = 0; nu < mu; nu++) {
|
||||
traceDirPlaquette(sitePlaq, U, mu, nu);
|
||||
Plaq = Plaq + sitePlaq;
|
||||
}
|
||||
}
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// sum over all x,y,z,t and over all planes of plaquette
|
||||
//////////////////////////////////////////////////
|
||||
static RealD sumPlaquette(const GaugeLorentz &Umu) {
|
||||
std::vector<GaugeMat> U(4, Umu._grid);
|
||||
|
||||
for (int mu = 0; mu < Nd; mu++) {
|
||||
U[mu] = PeekIndex<LorentzIndex>(Umu, mu);
|
||||
}
|
||||
|
||||
LatticeComplex Plaq(Umu._grid);
|
||||
|
||||
sitePlaquette(Plaq, U);
|
||||
|
||||
TComplex Tp = sum(Plaq);
|
||||
Complex p = TensorRemove(Tp);
|
||||
return p.real();
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// average over all x,y,z,t and over all planes of plaquette
|
||||
//////////////////////////////////////////////////
|
||||
static RealD avgPlaquette(const GaugeLorentz &Umu) {
|
||||
RealD sumplaq = sumPlaquette(Umu);
|
||||
double vol = Umu._grid->gSites();
|
||||
double faces = (1.0 * Nd * (Nd - 1)) / 2.0;
|
||||
return sumplaq / vol / faces / Nc; // Nd , Nc dependent... FIXME
|
||||
}
|
||||
|
||||
//////////////////////////////////////////////////
|
||||
// Wilson loop of size (R1, R2), oriented in mu,nu plane
|
||||
//////////////////////////////////////////////////
|
||||
static void wilsonLoop(GaugeMat &wl, const std::vector<GaugeMat> &U,
|
||||
const int Rmu, const int Rnu,
|
||||
const int mu, const int nu) {
|
||||
wl = U[nu];
|
||||
|
||||
for(int i = 0; i < Rnu-1; i++){
|
||||
wl = Gimpl::CovShiftForward(U[nu], nu, wl);
|
||||
}
|
||||
|
||||
for(int i = 0; i < Rmu; i++){
|
||||
wl = Gimpl::CovShiftForward(U[mu], mu, wl);
|
||||
}
|
||||
|
||||
for(int i = 0; i < Rnu; i++){
|
||||
wl = Gimpl::CovShiftBackward(U[nu], nu, wl);
|
||||
}
|
||||
|
||||
for(int i = 0; i < Rmu; i++){
|
||||
wl = Gimpl::CovShiftBackward(U[mu], mu, wl);
|
||||
}
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// trace of Wilson Loop oriented in mu,nu plane
|
||||
//////////////////////////////////////////////////
|
||||
static void traceWilsonLoop(LatticeComplex &wl,
|
||||
const std::vector<GaugeMat> &U,
|
||||
const int Rmu, const int Rnu,
|
||||
const int mu, const int nu) {
|
||||
GaugeMat sp(U[0]._grid);
|
||||
WilsonLoop(sp, U, Rmu, Rnu, mu, nu);
|
||||
wl = trace(sp);
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// sum over all planes of Wilson loop
|
||||
//////////////////////////////////////////////////
|
||||
static void siteWilsonLoop(LatticeComplex &Wl,
|
||||
const std::vector<GaugeMat> &U
|
||||
const int R1, const int R2) {
|
||||
LatticeComplex siteWl(U[0]._grid);
|
||||
Wl = zero;
|
||||
for (int mu = 1; mu < Nd; mu++) {
|
||||
for (int nu = 0; nu < mu; nu++) {
|
||||
traceWilsonLoop(siteWl, U, R1, R2, mu, nu);
|
||||
Wl = Wl + siteWl;
|
||||
traceWilsonLoop(siteWl, U, R2, R1, mu, nu);
|
||||
Wl = Wl + siteWl;
|
||||
}
|
||||
}
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// sum over all x,y,z,t and over all planes of Wilson loop
|
||||
//////////////////////////////////////////////////
|
||||
static RealD sumWilsonLoop(const GaugeLorentz &Umu,
|
||||
const int R1, const int R2) {
|
||||
std::vector<GaugeMat> U(4, Umu._grid);
|
||||
|
||||
for (int mu = 0; mu < Nd; mu++) {
|
||||
U[mu] = PeekIndex<LorentzIndex>(Umu, mu);
|
||||
}
|
||||
|
||||
LatticeComplex Wl(Umu._grid);
|
||||
|
||||
siteWilsonLoop(Wl, U, R1, R2);
|
||||
|
||||
TComplex Tp = sum(Wl);
|
||||
Complex p = TensorRemove(Tp);
|
||||
return p.real();
|
||||
}
|
||||
//////////////////////////////////////////////////
|
||||
// average over all x,y,z,t and over all planes of Wilson loop
|
||||
//////////////////////////////////////////////////
|
||||
static RealD avgPlaquette(const GaugeLorentz &Umu,
|
||||
const int R1, const int R2) {
|
||||
RealD sumWl = sumWilsonLoop(Umu);
|
||||
double vol = Umu._grid->gSites();
|
||||
double faces = 1.0 * Nd * (Nd - 1);
|
||||
return sumWl / vol / faces / Nc; // Nd , Nc dependent... FIXME
|
||||
}
|
||||
|
||||
END_QEDFVOL_NAMESPACE
|
||||
|
||||
#endif // QEDFVOL_WILSONLOOPS_H
|
||||
@ -64,6 +64,51 @@ int main(int argc, char *argv[])
|
||||
pRNG.SeedRandomDevice();
|
||||
photon.StochasticField(a, pRNG);
|
||||
|
||||
// Calculate log of plaquette
|
||||
EmComp plaqA(&grid);
|
||||
EmComp wlA(&grid);
|
||||
EmComp tmp(&grid);
|
||||
std::vector<EmComp> a_comp(4, &grid);
|
||||
|
||||
for (int dir = 0; dir < Nd; dir++) {
|
||||
a_comp[dir] = PeekIndex<LorentzIndex>(a, dir);
|
||||
}
|
||||
|
||||
plaqA = zero;
|
||||
wlA = zero;
|
||||
|
||||
for(int mu = 1; mu < Nd; mu++) {
|
||||
for(int nu = 0; nu < mu; nu++) {
|
||||
tmp = a_comp[mu] + Cshift(a_comp[nu], mu, 1) - Cshift(a_comp[mu], nu, 1) - a_comp[nu];
|
||||
plaqA = plaqA + cos(tmp);
|
||||
|
||||
tmp = a_comp[mu] + Cshift(a_comp[mu], mu, 1)
|
||||
+ Cshift(a_comp[nu], mu, 2) + Cshift(Cshift(a_comp[nu], mu, 2), nu, 1)
|
||||
- Cshift(Cshift(a_comp[mu], nu, 2), mu, 1) - Cshift(a_comp[mu], nu, 2)
|
||||
- Cshift(a_comp[nu], nu, 1) - a_comp[nu];
|
||||
wlA = wlA + cos(tmp);
|
||||
}
|
||||
}
|
||||
|
||||
double vol = grid.gSites();
|
||||
double faces = (1.0 * Nd * (Nd - 1)) / 2.0;
|
||||
|
||||
Complex avgPlaqA = sum(trace(plaqA));
|
||||
avgPlaqA = avgPlaqA / vol / faces;
|
||||
|
||||
Complex avgWlA = sum(trace(wlA));
|
||||
avgWlA = avgWlA / vol / faces;
|
||||
|
||||
TComplex tplaqsite;
|
||||
LatticeComplex plaqtrace = trace(plaqA);
|
||||
std::vector<int> site0 = {0,0,0,0};
|
||||
peekSite(tplaqsite, plaqtrace, site0);
|
||||
Complex plaqsite = TensorRemove(tplaqsite);
|
||||
|
||||
LOG(Message) << "Plaquette average: " << avgPlaqA << std::endl;
|
||||
LOG(Message) << "2x2 Wilson Loop average: " << avgWlA << std::endl;
|
||||
LOG(Message) << "Plaquette (one site): " << plaqsite / faces << std::endl;
|
||||
|
||||
// epilogue
|
||||
LOG(Message) << "Grid is finalizing now" << std::endl;
|
||||
Grid_finalize();
|
||||
|
||||
Loading…
Reference in New Issue
Block a user