portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2024-11-09 23:45:36 +00:00
Code Releases Activity

minor improvements

Browse Source
This commit is contained in:
Alessandro Lupo 2023-06-23 10:49:41 +01:00
parent 2372275b2c
commit de30c4e22a
8 changed files with 53 additions and 107 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
Grid/qcd/action/fermion/WilsonFermion.h
View File

@ -194,7 +194,5 @@ public:
typedef WilsonFermion<WilsonImplF> WilsonFermionF;
typedef WilsonFermion<WilsonImplD> WilsonFermionD;
//typedef WilsonFermion<SpWilsonImplF> SpWilsonFermionF;
//typedef WilsonFermion<SpWilsonImplD> SpWilsonFermionD;
NAMESPACE_END(Grid);

30
Grid/qcd/action/gauge/GaugeImplTypes.h
View File

@ -61,7 +61,6 @@ NAMESPACE_BEGIN(Grid);
typedef typename Impl::Field Field;
// hardcodes the exponential approximation in the template
//template <class S, int Nrepresentation = Nc, int Nexp = 12 > class GaugeImplTypes {
template <class S, int Nrepresentation = Nc, int Nexp = 12, class Group = SU<Nc> > class GaugeImplTypes {
public:
typedef S Simd;
@ -71,7 +70,6 @@ public:
template <typename vtype> using iImplGaugeLink = iScalar<iScalar<iMatrix<vtype, Nrepresentation> > >;
template <typename vtype> using iImplGaugeField = iVector<iScalar<iMatrix<vtype, Nrepresentation> >, Nd>;
typedef iImplScalar<Simd> SiteComplex;
typedef iImplGaugeLink<Simd> SiteLink;
typedef iImplGaugeField<Simd> SiteField;
@ -120,8 +118,7 @@ public:
LinkField Pmu(P.Grid());
Pmu = Zero();
for (int mu = 0; mu < Nd; mu++)
{
for (int mu = 0; mu < Nd; mu++) {
Group::GaussianFundamentalLieAlgebraMatrix(pRNG, Pmu);
RealD scale = ::sqrt(HMC_MOMENTUM_DENOMINATOR) ;
Pmu = Pmu*scale;
@ -129,11 +126,10 @@ public:
}
}
static inline Field projectForce(Field &P)
{
Field ret(P.Grid());
Group::taProj(P, ret);
return ret;
static inline Field projectForce(Field &P) {
Field ret(P.Grid());
Group::taProj(P, ret);
return ret;
}
static inline void update_field(Field& P, Field& U, double ep){
@ -164,23 +160,19 @@ public:
return Hsum.real();
}
static inline void Project(Field &U)
{
Group::ProjectOnSpecialGroup(U);
static inline void Project(Field &U) {
Group::ProjectOnSpecialGroup(U);
}
static inline void HotConfiguration(GridParallelRNG &pRNG, Field &U)
{
Group::HotConfiguration(pRNG, U);
static inline void HotConfiguration(GridParallelRNG &pRNG, Field &U) {
Group::HotConfiguration(pRNG, U);
}
static inline void TepidConfiguration(GridParallelRNG &pRNG, Field &U)
{
static inline void TepidConfiguration(GridParallelRNG &pRNG, Field &U) {
Group::TepidConfiguration(pRNG, U);
}
static inline void ColdConfiguration(GridParallelRNG &pRNG, Field &U)
{
static inline void ColdConfiguration(GridParallelRNG &pRNG, Field &U) {
Group::ColdConfiguration(pRNG, U);
}

4
Grid/qcd/hmc/GenericHMCrunner.h
View File

@ -228,11 +228,11 @@ using GenericHMCRunnerHirep =
// sp2n
template <template <typename, typename, typename> class Integrator>
using GenericSp2nHMCRunner = HMCWrapperTemplate<SpPeriodicGimplR, Integrator>;
using GenericSpHMCRunner = HMCWrapperTemplate<SpPeriodicGimplR, Integrator>;
template <class RepresentationsPolicy,
template <typename, typename, typename> class Integrator>
using GenericSp2nHMCRunnerHirep =
using GenericSpHMCRunnerHirep =
HMCWrapperTemplate<SpPeriodicGimplR, Integrator, RepresentationsPolicy>;

116
Grid/qcd/utils/GaugeGroup.h
View File

@ -29,8 +29,8 @@ See the full license in the file "LICENSE" in the top level distribution
directory
*************************************************************************************/
/* END LEGAL */
#ifndef QCD_UTIL_SUN_H
#define QCD_UTIL_SUN_H
#ifndef QCD_UTIL_GAUGEGROUP_H
#define QCD_UTIL_GAUGEGROUP_H
// Important detail: nvcc requires all template parameters to have names.
// This is the only reason why the second template parameter has a name.
@ -145,6 +145,12 @@ class GaugeGroup {
typedef Lattice<vSU2MatrixF> LatticeSU2MatrixF;
typedef Lattice<vSU2MatrixD> LatticeSU2MatrixD;
// Private implementation details are specified in the following files:
// Grid/qcd/utils/SUn.impl
// Grid/qcd/utils/SUn.impl
// The public part of the interface follows below and refers to these
// private member functions.
#include "Grid/qcd/utils/SUn.impl"
#include "Grid/qcd/utils/Sp2n.impl"
@ -161,7 +167,7 @@ class GaugeGroup {
static void testGenerators(void) { testGenerators(group_name()); }
static void printGenerators(void) {
for (int gen = 0; gen < AdjointDimension; gen++) {
for (int gen = 0; gen < AlgebraDimension; gen++) {
Matrix ta;
generator(gen, ta);
std::cout << GridLogMessage << "Nc = " << ncolour << " t_" << gen
@ -187,12 +193,11 @@ class GaugeGroup {
LatticeMatrixType lie(grid);
LatticeMatrixType la(grid);
ComplexD ci(0.0, scale);
// ComplexD cone(1.0, 0.0);
MatrixType ta;
lie = Zero();
for (int a = 0; a < AdjointDimension; a++) {
for (int a = 0; a < AlgebraDimension; a++) {
random(pRNG, ca);
ca = (ca + conjugate(ca)) * 0.5;
@ -217,7 +222,7 @@ class GaugeGroup {
Matrix ta;
out = Zero();
for (int a = 0; a < AdjointDimension; a++) {
for (int a = 0; a < AlgebraDimension; a++) {
gaussian(pRNG, ca);
generator(a, ta);
la = toComplex(ca) * ta;
@ -235,7 +240,7 @@ class GaugeGroup {
Matrix ta;
out = Zero();
for (int a = 0; a < AdjointDimension; a++) {
for (int a = 0; a < AlgebraDimension; a++) {
generator(a, ta);
la = peekColour(h, a) * timesI(ta) * scale;
out += la;
@ -250,7 +255,7 @@ class GaugeGroup {
h_out = Zero();
Matrix Ta;
for (int a = 0; a < AdjointDimension; a++) {
for (int a = 0; a < AlgebraDimension; a++) {
generator(a, Ta);
pokeColour(h_out, -2.0 * (trace(timesI(Ta) * in)) * scale, a);
}
@ -322,7 +327,7 @@ class GaugeGroup {
}
}
template <int N> // Projects on the general groups U(N), Sp(2N)xZ2 i.e. determinant il allowed a complex phase.
template <int N> // Projects on the general groups U(N), Sp(2N)xZ2 i.e. determinant is allowed a complex phase.
static void ProjectOnGeneralGroup(Lattice<iVector<iScalar<iMatrix<vComplexD, N> >, Nd> > &U) {
for (int mu = 0; mu < Nd; mu++) {
auto Umu = PeekIndex<LorentzIndex>(U, mu);
@ -350,7 +355,7 @@ class GaugeGroup {
return ProjectOnGeneralGroup(arg, group_name());
}
template <int N> // Projects on SU(N), Sp(2N)
template <int N> // Projects on SU(N), Sp(2N), with unit determinant, by first projecting on general group and then enforcing unit determinant
static void ProjectOnSpecialGroup(Lattice<iScalar<iScalar<iMatrix<vComplexD, N> > > > &Umu) {
Umu = ProjectOnGeneralGroup(Umu);
auto det = Determinant(Umu);
@ -375,6 +380,24 @@ class GaugeGroup {
}
};
template <int ncolour>
using SU = GaugeGroup<ncolour, GroupName::SU>;
template <int ncolour>
using Sp = GaugeGroup<ncolour, GroupName::Sp>;
typedef SU<2> SU2;
typedef SU<3> SU3;
typedef SU<4> SU4;
typedef SU<5> SU5;
typedef SU<Nc> FundamentalMatrices;
typedef Sp<2> Sp2;
typedef Sp<4> Sp4;
typedef Sp<6> Sp6;
typedef Sp<8> Sp8;
template <int N>
LatticeComplexD Determinant(
@ -403,32 +426,10 @@ LatticeComplexD Determinant(
}
template <int N>
static void ProjectSUn(
Lattice<iScalar<iScalar<iMatrix<vComplexD, N> > > > &Umu) {
Umu = ProjectOnGroup(Umu);
auto det = Determinant(Umu);
det = conjugate(det);
for (int i = 0; i < N; i++) {
auto element = PeekIndex<ColourIndex>(Umu, N - 1, i);
element = element * det;
PokeIndex<ColourIndex>(Umu, element, Nc - 1, i);
}
}
using ProjectSUn = typename GaugeGroup<N,GroupName::SU>::ProjectOnSpecialGroup;
template <int N>
static void ProjectSUn(
Lattice<iVector<iScalar<iMatrix<vComplexD, N> >, Nd> > &U) {
GridBase *grid = U.Grid();
// Reunitarise
for (int mu = 0; mu < Nd; mu++) {
auto Umu = PeekIndex<LorentzIndex>(U, mu);
Umu = ProjectOnGroup(Umu);
ProjectSUn(Umu);
PokeIndex<LorentzIndex>(U, Umu, mu);
}
}
using ProjectSpn = typename GaugeGroup<N,GroupName::Sp>::ProjectOnSpecialGroup;
// Explicit specialisation for SU(3).
static void ProjectSU3(
@ -463,50 +464,5 @@ static void ProjectSU3(
}
}
template <int ncolour>
using SU = GaugeGroup<ncolour, GroupName::SU>;
template <int ncolour>
using Sp = GaugeGroup<ncolour, GroupName::Sp>;
typedef SU<2> SU2;
typedef SU<3> SU3;
typedef SU<4> SU4;
typedef SU<5> SU5;
typedef SU<Nc> FundamentalMatrices;
template <int N>
static void ProjectSp2n(
Lattice<iScalar<iScalar<iMatrix<vComplexD, N> > > > &Umu) {
Umu = ProjectOnSpGroup(Umu);
auto det = Determinant(Umu); // ok ?
det = conjugate(det);
for (int i = 0; i < N; i++) {
auto element = PeekIndex<ColourIndex>(Umu, N - 1, i);
element = element * det;
PokeIndex<ColourIndex>(Umu, element, Nc - 1, i);
}
}
template <int N>
static void ProjectSp2n(
Lattice<iVector<iScalar<iMatrix<vComplexD, N> >, Nd> > &U) {
GridBase *grid = U.Grid();
for (int mu = 0; mu < Nd; mu++) {
auto Umu = PeekIndex<LorentzIndex>(U, mu);
Umu = ProjectOnSpGroup(Umu);
ProjectSp2n(Umu);
PokeIndex<LorentzIndex>(U, Umu, mu);
}
}
typedef Sp<2> Sp2;
typedef Sp<4> Sp4;
typedef Sp<6> Sp6;
typedef Sp<8> Sp8;
NAMESPACE_END(Grid);
#endif

2
tests/sp2n/Test_hmc_Sp_WF_2_Fund_3_2AS.cc
View File

@ -7,7 +7,7 @@ int main(int argc, char **argv) {
Grid_init(&argc, &argv);
typedef GenericSp2nHMCRunnerHirep<TheRepresentations, MinimumNorm2> HMCWrapper;
typedef GenericSpHMCRunnerHirep<TheRepresentations, MinimumNorm2> HMCWrapper;
typedef SpWilsonTwoIndexAntiSymmetricImplR TwoIndexFermionImplPolicy;
typedef SpWilsonTwoIndexAntiSymmetricFermionD TwoIndexFermionAction;

2
tests/sp2n/Test_hmc_Sp_Wilson2ASFermionGauge.cc
View File

@ -9,7 +9,7 @@ int main(int argc, char **argv) {
Grid_init(&argc, &argv);
typedef GenericSp2nHMCRunnerHirep<TheRepresentations, MinimumNorm2>
typedef GenericSpHMCRunnerHirep<TheRepresentations, MinimumNorm2>
HMCWrapper;
typedef SpWilsonTwoIndexAntiSymmetricImplR FermionImplPolicy;
typedef SpWilsonTwoIndexAntiSymmetricFermionD FermionAction;

2
tests/sp2n/Test_hmc_Sp_WilsonFundFermionGauge.cc
View File

@ -7,7 +7,7 @@ int main(int argc, char **argv) {
Grid_init(&argc, &argv);
typedef GenericSp2nHMCRunnerHirep<TheRepresentations, MinimumNorm2> HMCWrapper; // ok
typedef GenericSpHMCRunnerHirep<TheRepresentations, MinimumNorm2> HMCWrapper; // ok
typedef SpWilsonImplR FermionImplPolicy; // ok
typedef SpWilsonFermionD FermionAction; // ok
typedef typename FermionAction::FermionField FermionField; // ok?

2
tests/sp2n/Test_hmc_Sp_pureGaugeWilson.cc
View File

@ -37,7 +37,7 @@ int main(int argc, char **argv)
Grid_init(&argc, &argv);
GridLogLayout();
typedef GenericSp2nHMCRunner<MinimumNorm2> HMCWrapper;
typedef GenericSpHMCRunner<MinimumNorm2> HMCWrapper;
HMCWrapper TheHMC;