From 943fa48ce48ce1f6b5886b35d23a687047a6531c Mon Sep 17 00:00:00 2001 From: Vera Guelpers Date: Fri, 14 Dec 2018 13:45:30 +0000 Subject: [PATCH] Hadrons: Kl2 contraction using sequential propagators --- Hadrons/Modules.hpp | 1 + .../Modules/MContraction/WeakMesonDecayKl2.cc | 36 +++ .../MContraction/WeakMesonDecayKl2.hpp | 213 ++++++++++++++++++ Hadrons/Modules/MFermion/FreeProp.hpp | 5 +- Hadrons/modules.inc | 2 + 5 files changed, 256 insertions(+), 1 deletion(-) create mode 100644 Hadrons/Modules/MContraction/WeakMesonDecayKl2.cc create mode 100644 Hadrons/Modules/MContraction/WeakMesonDecayKl2.hpp diff --git a/Hadrons/Modules.hpp b/Hadrons/Modules.hpp index 1305e6b7..88caf8eb 100644 --- a/Hadrons/Modules.hpp +++ b/Hadrons/Modules.hpp @@ -9,6 +9,7 @@ #include #include #include +#include #include #include #include diff --git a/Hadrons/Modules/MContraction/WeakMesonDecayKl2.cc b/Hadrons/Modules/MContraction/WeakMesonDecayKl2.cc new file mode 100644 index 00000000..78f1a6ae --- /dev/null +++ b/Hadrons/Modules/MContraction/WeakMesonDecayKl2.cc @@ -0,0 +1,36 @@ +/************************************************************************************* + +Grid physics library, www.github.com/paboyle/Grid + +Source file: Hadrons/Modules/MContraction/WeakMesonDecayKl2.cc + +Copyright (C) 2015-2018 + +Author: Antonin Portelli +Author: Vera Guelpers + +This program is free software; you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation; either version 2 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License along +with this program; if not, write to the Free Software Foundation, Inc., +51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + +See the full license in the file "LICENSE" in the top level distribution directory +*************************************************************************************/ +/* END LEGAL */ +#include + +using namespace Grid; +using namespace Hadrons; +using namespace MContraction; + +template class Grid::Hadrons::MContraction::TWeakMesonDecayKl2; + diff --git a/Hadrons/Modules/MContraction/WeakMesonDecayKl2.hpp b/Hadrons/Modules/MContraction/WeakMesonDecayKl2.hpp new file mode 100644 index 00000000..4e77f3aa --- /dev/null +++ b/Hadrons/Modules/MContraction/WeakMesonDecayKl2.hpp @@ -0,0 +1,213 @@ +/************************************************************************************* + +Grid physics library, www.github.com/paboyle/Grid + +Source file: Hadrons/Modules/MContraction/WeakMesonDecayKl2.hpp + +Copyright (C) 2015-2018 + +Author: Antonin Portelli +Author: Vera Guelpers + + +This program is free software; you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation; either version 2 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License along +with this program; if not, write to the Free Software Foundation, Inc., +51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + +See the full license in the file "LICENSE" in the top level distribution directory +*************************************************************************************/ +/* END LEGAL */ + +#ifndef Hadrons_MContraction_WeakMesonDecayKl2_hpp_ +#define Hadrons_MContraction_WeakMesonDecayKl2_hpp_ + +#include +#include +#include +#include + +BEGIN_HADRONS_NAMESPACE + +/* +* Kl2 contraction +* ----------------------------- +* +* contraction for Kl2 decay, including the lepton +* +* trace(q1*adj(q2)*g5*gL[mu]) * (gL[mu] * lepton)_{a,b} +* +* with open spinor indices (a,b) for the lepton part +* +* q1 lepton +* /------------\ /------------ +* / \ / +* / \H_W/ +* g_5 * * * +* \ / +* \ / +* \____________/ +* q2 +* +* * options: +* - q1: input propagator 1 (string) +* - q2: input propagator 2 (string) +* - lepton: input lepton (string) +*/ + +/****************************************************************************** + * TWeakMesonDecayKl2 * + ******************************************************************************/ +BEGIN_MODULE_NAMESPACE(MContraction) + +class WeakMesonDecayKl2Par: Serializable +{ +public: + GRID_SERIALIZABLE_CLASS_MEMBERS(WeakMesonDecayKl2Par, + std::string, q1, + std::string, q2, + std::string, lepton, + std::string, output); +}; + +template +class TWeakMesonDecayKl2: public Module +{ +public: + FERM_TYPE_ALIASES(FImpl,); + class Result: Serializable + { + public: + GRID_SERIALIZABLE_CLASS_MEMBERS(Result, + int, spinidx1, + int, spinidx2, + std::vector, corr); + }; +public: + // constructor + TWeakMesonDecayKl2(const std::string name); + // destructor + virtual ~TWeakMesonDecayKl2(void) {}; + // dependencies/products + virtual std::vector getInput(void); + virtual std::vector getOutput(void); +protected: + // execution + virtual void setup(void); + // execution + virtual void execute(void); +}; + +MODULE_REGISTER_TMP(WeakMesonDecayKl2, ARG(TWeakMesonDecayKl2), MContraction); + +/****************************************************************************** + * TWeakMesonDecayKl2 implementation * + ******************************************************************************/ +// constructor ///////////////////////////////////////////////////////////////// +template +TWeakMesonDecayKl2::TWeakMesonDecayKl2(const std::string name) +: Module(name) +{} + +// dependencies/products /////////////////////////////////////////////////////// +template +std::vector TWeakMesonDecayKl2::getInput(void) +{ + std::vector input = {par().q1, par().q2, par().lepton}; + + return input; +} + +template +std::vector TWeakMesonDecayKl2::getOutput(void) +{ + std::vector output = {}; + + return output; +} + +// setup //////////////////////////////////////////////////////////////////////// +template +void TWeakMesonDecayKl2::setup(void) +{ + envTmpLat(LatticeComplex, "c"); + envTmpLat(PropagatorField, "prop_buf"); + envCreateLat(PropagatorField, getName()); + envTmpLat(LatticeComplex, "buf"); +} + +// execution /////////////////////////////////////////////////////////////////// +template +void TWeakMesonDecayKl2::execute(void) +{ + LOG(Message) << "Computing QED Kl2 contractions '" << getName() << "' using" + << " quarks '" << par().q1 << "' and '" << par().q2 << "' and" + << "lepton '" << par().lepton << "'" << std::endl; + + + auto &res = envGet(PropagatorField, getName()); res = zero; + Gamma g5(Gamma::Algebra::Gamma5); + int nt = env().getDim(Tp); + + auto &q1 = envGet(PropagatorField, par().q1); + auto &q2 = envGet(PropagatorField, par().q2); + auto &lepton = envGet(PropagatorField, par().lepton); + envGetTmp(LatticeComplex, buf); + std::vector res_summed; + envGetTmp(LatticeComplex, c); + envGetTmp(PropagatorField, prop_buf); + + std::vector result; + result.resize(Ns*Ns); + for (unsigned int i = 0; i < Ns*Ns ; ++i) + { + result[i].corr.resize(nt); + } + + + for (unsigned int mu = 0; mu < 4; ++mu) + { + c = zero; + //hadronic part: trace(q1*adj(q2)*g5*gL[mu]) + c = trace(MAKE_CW_SUBDIAG(q1, q2, GammaL(Gamma::gmu[mu]))); + prop_buf = 1.; + //multiply lepton part + res += c * prop_buf * GammaL(Gamma::gmu[mu]) * lepton; + } + + //loop over spinor index of lepton part + unsigned int i = 0; + for (unsigned int s1 = 0; s1 < Ns ; ++s1) + for (unsigned int s2 = 0; s2 < Ns ; ++s2) + { + buf = peekColour(peekSpin(res,s1,s2),0,0); + + sliceSum(buf, res_summed, Tp); + + for (unsigned int t = 0; t < nt; ++t) + { + result[i].corr[t] = TensorRemove(res_summed[t]); + } + result[i].spinidx1 = s1; + result[i].spinidx2 = s2; + + i+=1; + } + + saveResult(par().output, "weakdecay", result); +} + +END_MODULE_NAMESPACE + +END_HADRONS_NAMESPACE + +#endif // Hadrons_MContraction_WeakMesonDecayKl2_hpp_ diff --git a/Hadrons/Modules/MFermion/FreeProp.hpp b/Hadrons/Modules/MFermion/FreeProp.hpp index b5d7b1a3..6b076519 100644 --- a/Hadrons/Modules/MFermion/FreeProp.hpp +++ b/Hadrons/Modules/MFermion/FreeProp.hpp @@ -168,7 +168,10 @@ void TFreeProp::execute(void) } sol = zero; std::vector twist = strToVec(par().twist); - if(twist.size() != Nd) HADRONS_ERROR(Size, "number of twist angles does not match number of dimensions"); + if(twist.size() != Nd) + { + HADRONS_ERROR(Size, "number of twist angles does not match number of dimensions"); + } mat.FreePropagator(source,sol,mass,twist); FermToProp(prop, sol, s, c); // create 4D propagators from 5D one if necessary diff --git a/Hadrons/modules.inc b/Hadrons/modules.inc index 1e1f1f4e..0ff0f77b 100644 --- a/Hadrons/modules.inc +++ b/Hadrons/modules.inc @@ -4,6 +4,7 @@ modules_cc =\ Modules/MContraction/Meson.cc \ Modules/MContraction/WeakNeutral4ptDisc.cc \ Modules/MContraction/WeakHamiltonianNonEye.cc \ + Modules/MContraction/WeakMesonDecayKl2.cc \ Modules/MContraction/A2AAslashField.cc \ Modules/MContraction/WardIdentity.cc \ Modules/MContraction/A2AMesonField.cc \ @@ -83,6 +84,7 @@ modules_hpp =\ Modules/MContraction/Gamma3pt.hpp \ Modules/MContraction/WardIdentity.hpp \ Modules/MContraction/WeakHamiltonianEye.hpp \ + Modules/MContraction/WeakMesonDecayKl2.hpp \ Modules/MFermion/FreeProp.hpp \ Modules/MFermion/GaugeProp.hpp \ Modules/MSource/SeqGamma.hpp \