From 2ce7f2b4d8578a26a259ef3a8634110a98300101 Mon Sep 17 00:00:00 2001 From: Felix Erben Date: Tue, 8 Oct 2019 13:19:47 +0100 Subject: [PATCH] suggested changes for 1st pull request implemented --- Grid/qcd/utils/BaryonUtils.h | 44 ++++++++++--------------- Hadrons/Modules/MContraction/Baryon.hpp | 18 ++++++++-- tests/hadrons/Test_hadrons_spectrum.cc | 6 ++++ 3 files changed, 39 insertions(+), 29 deletions(-) diff --git a/Grid/qcd/utils/BaryonUtils.h b/Grid/qcd/utils/BaryonUtils.h index 5d462990..b68055e0 100644 --- a/Grid/qcd/utils/BaryonUtils.h +++ b/Grid/qcd/utils/BaryonUtils.h @@ -31,8 +31,6 @@ NAMESPACE_BEGIN(Grid); -#undef DELTA_F_EQ_2 - template class BaryonUtils { @@ -43,6 +41,9 @@ public: typedef typename FImpl::SitePropagator pobj; typedef typename ComplexField::vector_object vobj; + + static constexpr int epsilon[6][3] = {{0,1,2},{1,2,0},{2,0,1},{0,2,1},{2,1,0},{1,0,2}}; + static constexpr Complex epsilon_sgn[6]= {1,1,1,-1,-1,-1}; private: template @@ -52,7 +53,7 @@ public: const Gamma GammaA, const Gamma GammaB, const int parity, - const std::vector &wick_contractions, + const int * wick_contractions, robj &result); public: static void ContractBaryons(const PropagatorField &q1_src, @@ -76,6 +77,11 @@ public: robj &result); }; +template +constexpr int BaryonUtils::epsilon[6][3]; +template +constexpr Complex BaryonUtils::epsilon_sgn[6]; + template template void BaryonUtils::baryon_site(const mobj &D1, @@ -84,16 +90,12 @@ void BaryonUtils::baryon_site(const mobj &D1, const Gamma GammaA, const Gamma GammaB, const int parity, - const std::vector &wick_contraction, + const int * wick_contraction, robj &result) { Gamma g4(Gamma::Algebra::GammaT); //needed for parity P_\pm = 0.5*(1 \pm \gamma_4) - - static const int epsilon[6][3] = {{0,1,2},{1,2,0},{2,0,1},{0,2,1},{2,1,0},{1,0,2}}; - static const Complex epsilon_sgn[6]= {1,1,1,-1,-1,-1}; - auto gD1a = GammaA * GammaA * D1; auto gD1b = GammaA * g4 * GammaA * D1; auto pD1 = 0.5* (gD1a + (double)parity * gD1b); @@ -176,7 +178,7 @@ void BaryonUtils::ContractBaryons(const PropagatorField &q1_src, const int parity, ComplexField &baryon_corr) { - std::cout << "quarks_snk " << quarks_snk[0] << quarks_snk[1] << quarks_snk[2] << std::endl; + std::cout << "Contraction <" << quarks_snk[0] << quarks_snk[1] << quarks_snk[2] << "|" << quarks_src[0] << quarks_src[1] << quarks_src[2] << ">" << std::endl; std::cout << "GammaA " << (GammaA.g) << std::endl; std::cout << "GammaB " << (GammaB.g) << std::endl; @@ -184,17 +186,10 @@ void BaryonUtils::ContractBaryons(const PropagatorField &q1_src, GridBase *grid = q1_src.Grid(); - Gamma g4(Gamma::Algebra::GammaT); //needed for parity P_\pm = 0.5*(1 \pm \gamma_4) - - static const int epsilon[6][3] = {{0,1,2},{1,2,0},{2,0,1},{0,2,1},{2,1,0},{1,0,2}}; - static const int epsilon_sgn[6]= {1,1,1,-1,-1,-1}; - std::vector wick_contraction = {0,0,0,0,0,0}; - + int wick_contraction[6]; for (int ie=0; ie < 6 ; ie++) - if (quarks_src[0] == quarks_snk[epsilon[ie][0]] && quarks_src[1] == quarks_snk[epsilon[ie][1]] && quarks_src[2] == quarks_snk[epsilon[ie][2]]) - wick_contraction[ie]=1; + wick_contraction[ie] = (quarks_src[0] == quarks_snk[epsilon[ie][0]] && quarks_src[1] == quarks_snk[epsilon[ie][1]] && quarks_src[2] == quarks_snk[epsilon[ie][2]]) ? 1 : 0; -// typedef typename ComplexField::vector_object vobj; auto vbaryon_corr= baryon_corr.View(); auto v1 = q1_src.View(); auto v2 = q2_src.View(); @@ -225,18 +220,15 @@ void BaryonUtils::ContractBaryons_Sliced(const mobj &D1, const int parity, robj &result) { + std::cout << "Contraction <" << quarks_snk[0] << quarks_snk[1] << quarks_snk[2] << "|" << quarks_src[0] << quarks_src[1] << quarks_src[2] << ">" << std::endl; + std::cout << "GammaA " << (GammaA.g) << std::endl; + std::cout << "GammaB " << (GammaB.g) << std::endl; assert(parity==1 || parity == -1 && "Parity must be +1 or -1"); - Gamma g4(Gamma::Algebra::GammaT); //needed for parity P_\pm = 0.5*(1 \pm \gamma_4) - - static const int epsilon[6][3] = {{0,1,2},{1,2,0},{2,0,1},{0,2,1},{2,1,0},{1,0,2}}; - static const int epsilon_sgn[6]= {1,1,1,-1,-1,-1}; - std::vector wick_contraction = {0,0,0,0,0,0}; - + int wick_contraction[6]; for (int ie=0; ie < 6 ; ie++) - if (quarks_src[0] == quarks_snk[epsilon[ie][0]] && quarks_src[1] == quarks_snk[epsilon[ie][1]] && quarks_src[2] == quarks_snk[epsilon[ie][2]]) - wick_contraction[ie]=1; + wick_contraction[ie] = (quarks_src[0] == quarks_snk[epsilon[ie][0]] && quarks_src[1] == quarks_snk[epsilon[ie][1]] && quarks_src[2] == quarks_snk[epsilon[ie][2]]) ? 1 : 0; result=Zero(); baryon_site(D1,D2,D3,GammaA,GammaB,parity,wick_contraction,result); diff --git a/Hadrons/Modules/MContraction/Baryon.hpp b/Hadrons/Modules/MContraction/Baryon.hpp index 7d3f61ec..097b5c71 100644 --- a/Hadrons/Modules/MContraction/Baryon.hpp +++ b/Hadrons/Modules/MContraction/Baryon.hpp @@ -71,6 +71,11 @@ public: { public: GRID_SERIALIZABLE_CLASS_MEMBERS(Result, + Gamma::Algebra, gammaA, + Gamma::Algebra, gammaB, + std::string, quarks, + std::string, prefactors, + int, parity, std::vector, corr); }; public: @@ -144,15 +149,14 @@ void TBaryon::execute(void) for (int iQ1 = 0; iQ1 < nQ; iQ1++) for (int iQ2 = 0; iQ2 < nQ; iQ2++) LOG(Message) << prefactors[iQ1]*prefactors[iQ2] << "*<" << quarks[iQ1] << "|" << quarks[iQ2] << ">" << std::endl; - LOG(Message) << " using quarks" << par().q1 << "', " << par().q2 << "', and '" + LOG(Message) << " using quarks " << par().q1 << "', " << par().q2 << "', and '" << par().q3 << "' and (Gamma^A,Gamma^B) = ( " << par().GammaA << " , " << par().GammaB << " ) and parity " << parity << " using sink " << par().sink << "." << std::endl; + envGetTmp(LatticeComplex, c); envGetTmp(LatticeComplex, c2); - Result result; int nt = env().getDim(Tp); - result.corr.resize(nt); std::vector ggA = strToVec(par().GammaA); Gamma GammaA(ggA[0]); std::vector ggB = strToVec(par().GammaB); @@ -161,6 +165,14 @@ void TBaryon::execute(void) TComplex cs; TComplex ch; + Result result; + result.corr.resize(nt); + result.gammaA = ggA[0]; + result.gammaB = ggB[0]; + result.parity = parity; + result.quarks = par().quarks; + result.prefactors = par().prefactors; + if (envHasType(SlicedPropagator1, par().q1) and envHasType(SlicedPropagator2, par().q2) and envHasType(SlicedPropagator3, par().q3)) diff --git a/tests/hadrons/Test_hadrons_spectrum.cc b/tests/hadrons/Test_hadrons_spectrum.cc index af1dccd7..6c6194c9 100644 --- a/tests/hadrons/Test_hadrons_spectrum.cc +++ b/tests/hadrons/Test_hadrons_spectrum.cc @@ -46,6 +46,7 @@ int main(int argc, char *argv[]) // run setup /////////////////////////////////////////////////////////////// Application application; std::vector flavour = {"l", "s", "c1", "c2", "c3"}; + std::vector flavour_baryon = {"l", "s", "a", "b", "c"}; //needs to be a single character std::vector mass = {.01, .04, .2 , .25 , .3 }; // global parameters @@ -134,6 +135,11 @@ int main(int argc, char *argv[]) barPar.q1 = "Qpt_" + flavour[i]; barPar.q2 = "Qpt_" + flavour[j]; barPar.q3 = "Qpt_" + flavour[k]; + barPar.GammaA = "Identity"; + barPar.GammaB = "GammaZGamma5"; //C*GammaX + barPar.quarks = flavour_baryon[i] + flavour_baryon[j] + flavour_baryon[k]; + barPar.prefactors = "1.0"; + barPar.sink = "sink"; application.createModule( "baryon_pt_" + flavour[i] + flavour[j] + flavour[k], barPar); }