From 94ea84d83f132b4c3f0f3c4fe93bdde4ad1bc84a Mon Sep 17 00:00:00 2001 From: Azusa Yamaguchi Date: Fri, 5 Jun 2015 10:02:36 +0100 Subject: [PATCH] Adding some wilson loop support --- lib/qcd/action/gauge/GaugeActionBase.h | 20 ++++ lib/qcd/action/gauge/WilsonGaugeAction.h | 24 +++++ lib/qcd/utils/CovariantCshift.h | 22 ++++ lib/qcd/utils/WilsonLoops.h | 130 +++++++++++++++++++++++ 4 files changed, 196 insertions(+) create mode 100644 lib/qcd/action/gauge/GaugeActionBase.h create mode 100644 lib/qcd/action/gauge/WilsonGaugeAction.h create mode 100644 lib/qcd/utils/CovariantCshift.h create mode 100644 lib/qcd/utils/WilsonLoops.h diff --git a/lib/qcd/action/gauge/GaugeActionBase.h b/lib/qcd/action/gauge/GaugeActionBase.h new file mode 100644 index 00000000..55141965 --- /dev/null +++ b/lib/qcd/action/gauge/GaugeActionBase.h @@ -0,0 +1,20 @@ +#ifndef QCD_GAUGE_ACTION_BASE +#define QCD_GAUGE_ACTION_BASE +namespace Grid { +namespace QCD{ + +template +class GaugeActionBase { // derive this from TermInAction? + + public: + virtual RealD S(const GaugeField &U) = 0; // evaluate the action + virtual void deriv(const GaugeField &U,GaugeField & dSdU ) = 0; // evaluate the action derivative + virtual void staple(const GaugeField &U,GaugeField & dSdU ) = 0; // evaluate the action derivative + virtual void refresh(const GaugeField & ) {}; + // Boundary conditions? + // Heatbath? + virtual ~GaugeActionBase() {}; +}; + +}} +#endif diff --git a/lib/qcd/action/gauge/WilsonGaugeAction.h b/lib/qcd/action/gauge/WilsonGaugeAction.h new file mode 100644 index 00000000..8ecb1f81 --- /dev/null +++ b/lib/qcd/action/gauge/WilsonGaugeAction.h @@ -0,0 +1,24 @@ +#ifndef QCD_WILSON_GAUGE_ACTION_H +#define QCD_WILSON_GAUGE_ACTION_H + +//////////////////////////////////////////////////////////////////////// +// Wilson Gauge Action .. should I template the Nc etc.. +//////////////////////////////////////////////////////////////////////// +template +class WilsonGaugeAction : public GaugeActionBase { + public: + virtual RealD S(const GaugeField &U) { + return WilsonLoops::sumPlaquette(U); + }; + virtual RealD deriv(GaugeField &U,GaugeField & dSdU ) { + WilsonLoops::Staple(dSdU,U,mu); + }; + virtual void staple(const MatrixField &stap,GaugeField & U,int mu ) { + WilsonLoops::Staple(stap,U,mu); + }; +}; + + +#endif + +#endif diff --git a/lib/qcd/utils/CovariantCshift.h b/lib/qcd/utils/CovariantCshift.h new file mode 100644 index 00000000..0c621069 --- /dev/null +++ b/lib/qcd/utils/CovariantCshift.h @@ -0,0 +1,22 @@ +#ifndef QCD_UTILS_COVARIANT_CSHIFT_H +#define QCD_UTILS_COVARIANT_CSHIFT_H +namespace Grid { +namespace QCD { +//////////////////////////////////////////////////////////////////////// +// Low performance implementation of CovariantCshift API +//////////////////////////////////////////////////////////////////////// +template Lattice CovShiftForward(const Lattice &Link, + int mu, + const Lattice &field) +{ + return Link*Cshift(field,mu,1);// moves towards negative mu +} +template Lattice CovShiftBackward(const Lattice &Link, + int mu, + const Lattice &field) +{ + Lattice tmp(field._grid); + tmp = adj(Link)*field; + return Cshift(tmp,mu,-1);// moves towards positive mu +} +}} diff --git a/lib/qcd/utils/WilsonLoops.h b/lib/qcd/utils/WilsonLoops.h new file mode 100644 index 00000000..f2d6d4de --- /dev/null +++ b/lib/qcd/utils/WilsonLoops.h @@ -0,0 +1,130 @@ +#ifndef QCD_UTILS_WILSON_LOOPS_H +#define QCD_UTILS_WILSON_LOOPS_H +namespace Grid { +namespace QCD { + +// Common wilson loop observables +template +class WilsonLoops { +public: + + ////////////////////////////////////////////////// + // directed plaquette oriented in mu,nu plane + ////////////////////////////////////////////////// + static void dirPlaquette(GaugeMat &plaq,const std::vector &U, const int mu, const int nu) + { + plaq=CovShiftForward(U[mu],mu,U[nu])*adj(CovShiftForward(U[nu],nu,U[mu])); + } + ////////////////////////////////////////////////// + // trace of directed plaquette oriented in mu,nu plane + ////////////////////////////////////////////////// + static void traceDirPlaquette(LatticeComplex &plaq, const std::vector &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 &U) + { + LatticeComplex sitePlaq(U[0]._grid); + for(int mu=1;mu U(4,Umu._grid); + for(int mu=0;mu(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 + } + ////////////////////////////////////////////////// + // the sum over all staples on each site + ////////////////////////////////////////////////// + static void Staple(GaugeMat &staple,GaugeLorentz &Umu,int mu){ + + std::vector U(4,Umu._grid); + for(int d=0;d(Umu,d); + } + + staple = zero; + GaugeMat tmp(U[0]._grid); + + for(int nu=0;nu &U, const int mu, const int nu){ +RealD avgRectangle(const std::vector &U){} +RealD avgRectangle(const std::vector &U, const int mu, const int nu){} +void traceRectangle(LatticeComplex &plaq,const std::vector &U, const int mu, const int nu){} +void siteRectangle(GaugeMat &plaq,const std::vector &U, const int mu, const int nu){} + */ + +}; +}} + + + + + + +#endif