1
0
mirror of https://github.com/paboyle/Grid.git synced 2025-04-09 21:50:45 +01:00

Hadrons: first cleaning/integration of A2A/meson fields

This commit is contained in:
Antonin Portelli 2018-08-06 12:11:52 +01:00
parent d2650e89bd
commit 3f0f92cda6
7 changed files with 89 additions and 113 deletions

View File

@ -51,6 +51,16 @@ class A2AModesSchurDiagTwo
v_high_4d.resize(size_4d, grid_4d); v_high_4d.resize(size_4d, grid_4d);
} }
int get_Nh(void) const
{
return Nh;
}
int get_Nl(void) const
{
return Nl;
}
void high_modes(Field &source_5d, Field &w_source_5d, Field &source_4d, int i) void high_modes(Field &source_5d, Field &w_source_5d, Field &source_4d, int i)
{ {
int i5d; int i5d;

View File

@ -142,7 +142,6 @@ void TimeDilutedSpinColorDiagonalNoise<FImpl>::generateNoise(GridParallelRNG &rn
LatticeCoordinate(tLat, nd - 1); LatticeCoordinate(tLat, nd - 1);
bernoulli(rng, eta); bernoulli(rng, eta);
eta = (2.*eta - shift)*(1./::sqrt(2.)); eta = (2.*eta - shift)*(1./::sqrt(2.));
std::cout << eta << std::endl;
for (unsigned int t = 0; t < nt_; ++t) for (unsigned int t = 0; t < nt_; ++t)
{ {
etaCut = where((tLat == t), eta, 0.*eta); etaCut = where((tLat == t), eta, 0.*eta);
@ -154,7 +153,6 @@ void TimeDilutedSpinColorDiagonalNoise<FImpl>::generateNoise(GridParallelRNG &rn
{ {
noise[i] = zero; noise[i] = zero;
pokeColour(noise[i], etas, c); pokeColour(noise[i], etas, c);
std::cout << noise[i] << std::endl;
i++; i++;
} }
} }

View File

@ -24,8 +24,6 @@ class A2AMesonFieldPar : Serializable
int, cacheBlock, int, cacheBlock,
int, schurBlock, int, schurBlock,
int, Nmom, int, Nmom,
int, N,
int, Nl,
std::string, A2A, std::string, A2A,
std::string, output); std::string, output);
}; };
@ -82,9 +80,7 @@ TA2AMesonField<FImpl>::TA2AMesonField(const std::string name)
template <typename FImpl> template <typename FImpl>
std::vector<std::string> TA2AMesonField<FImpl>::getInput(void) std::vector<std::string> TA2AMesonField<FImpl>::getInput(void)
{ {
std::vector<std::string> in = {par().A2A + "_class"}; std::vector<std::string> in = {par().A2A};
in.push_back(par().A2A + "_w_high_4d");
in.push_back(par().A2A + "_v_high_4d");
return in; return in;
} }
@ -102,18 +98,17 @@ std::vector<std::string> TA2AMesonField<FImpl>::getOutput(void)
template <typename FImpl> template <typename FImpl>
void TA2AMesonField<FImpl>::setup(void) void TA2AMesonField<FImpl>::setup(void)
{ {
auto &a2a = envGet(A2ABase, par().A2A + "_class"); auto &a2a = envGet(A2ABase, par().A2A);
int nt = env().getDim(Tp); int Ls = env().getObjectLs(par().A2A);
int Nl = par().Nl;
int N = par().N;
int Ls_ = env().getObjectLs(par().A2A + "_class");
// Four D fields // Four D fields
envTmp(std::vector<FermionField>, "w", 1, par().schurBlock, FermionField(env().getGrid(1))); envTmp(std::vector<FermionField>, "w", 1, par().schurBlock,
envTmp(std::vector<FermionField>, "v", 1, par().schurBlock, FermionField(env().getGrid(1))); FermionField(env().getGrid()));
envTmp(std::vector<FermionField>, "v", 1, par().schurBlock,
FermionField(env().getGrid()));
// 5D tmp // 5D tmp
envTmpLat(FermionField, "tmp_5d", Ls_); envTmpLat(FermionField, "tmp_5d", Ls);
} }
@ -305,7 +300,7 @@ void TA2AMesonField<FImpl>::execute(void)
{ {
LOG(Message) << "Computing A2A meson field" << std::endl; LOG(Message) << "Computing A2A meson field" << std::endl;
auto &a2a = envGet(A2ABase, par().A2A + "_class"); auto &a2a = envGet(A2ABase, par().A2A);
// 2+6+4+4 = 16 gammas // 2+6+4+4 = 16 gammas
// Ordering defined here // Ordering defined here
@ -335,8 +330,9 @@ void TA2AMesonField<FImpl>::execute(void)
int nx = env().getDim(Xp); int nx = env().getDim(Xp);
int ny = env().getDim(Yp); int ny = env().getDim(Yp);
int nz = env().getDim(Zp); int nz = env().getDim(Zp);
int N = par().N; int Nl = a2a.get_Nl();
int Nl = par().Nl; int N = Nl + a2a.get_Nh();
int ngamma = gammas.size(); int ngamma = gammas.size();
int schurBlock = par().schurBlock; int schurBlock = par().schurBlock;
@ -407,7 +403,6 @@ void TA2AMesonField<FImpl>::execute(void)
int N_iii = MIN(N_ii-ii,cacheBlock); int N_iii = MIN(N_ii-ii,cacheBlock);
int N_jjj = MIN(N_jj-jj,cacheBlock); int N_jjj = MIN(N_jj-jj,cacheBlock);
Eigen::Tensor<ComplexD,5> mesonFieldBlocked(nmom,ngamma,nt,N_iii,N_jjj); Eigen::Tensor<ComplexD,5> mesonFieldBlocked(nmom,ngamma,nt,N_iii,N_jjj);
t_contr-=usecond(); t_contr-=usecond();

View File

@ -5,3 +5,4 @@ using namespace Hadrons;
using namespace MNoise; using namespace MNoise;
template class Grid::Hadrons::MNoise::TTimeDilutedSpinColorDiagonal<FIMPL>; template class Grid::Hadrons::MNoise::TTimeDilutedSpinColorDiagonal<FIMPL>;
template class Grid::Hadrons::MNoise::TTimeDilutedSpinColorDiagonal<ZFIMPL>;

View File

@ -31,6 +31,7 @@ public:
}; };
MODULE_REGISTER_TMP(TimeDilutedSpinColorDiagonal, TTimeDilutedSpinColorDiagonal<FIMPL>, MNoise); MODULE_REGISTER_TMP(TimeDilutedSpinColorDiagonal, TTimeDilutedSpinColorDiagonal<FIMPL>, MNoise);
MODULE_REGISTER_TMP(ZTimeDilutedSpinColorDiagonal, TTimeDilutedSpinColorDiagonal<ZFIMPL>, MNoise);
/****************************************************************************** /******************************************************************************
* TTimeDilutedSpinColorDiagonal implementation * * TTimeDilutedSpinColorDiagonal implementation *

View File

@ -4,5 +4,5 @@ using namespace Grid;
using namespace Hadrons; using namespace Hadrons;
using namespace MSolver; using namespace MSolver;
template class Grid::Hadrons::MSolver::TA2AVectors<FIMPL, HADRONS_DEFAULT_LANCZOS_NBASIS>; template class Grid::Hadrons::MSolver::TA2AVectors<FIMPL, FermionEigenPack<FIMPL>>;
template class Grid::Hadrons::MSolver::TA2AVectors<ZFIMPL, HADRONS_DEFAULT_LANCZOS_NBASIS>; template class Grid::Hadrons::MSolver::TA2AVectors<ZFIMPL, FermionEigenPack<ZFIMPL>>;

View File

@ -7,6 +7,7 @@
#include <Grid/Hadrons/Solver.hpp> #include <Grid/Hadrons/Solver.hpp>
#include <Grid/Hadrons/EigenPack.hpp> #include <Grid/Hadrons/EigenPack.hpp>
#include <Grid/Hadrons/AllToAllVectors.hpp> #include <Grid/Hadrons/AllToAllVectors.hpp>
#include <Grid/Hadrons/DilutedNoise.hpp>
BEGIN_HADRONS_NAMESPACE BEGIN_HADRONS_NAMESPACE
@ -21,26 +22,20 @@ public:
GRID_SERIALIZABLE_CLASS_MEMBERS(A2AVectorsPar, GRID_SERIALIZABLE_CLASS_MEMBERS(A2AVectorsPar,
bool, return_5d, bool, return_5d,
int, Nl, int, Nl,
int, N, std::string, noise,
std::vector<std::string>, sources,
std::string, action, std::string, action,
std::string, eigenPack, std::string, eigenPack,
std::string, solver); std::string, solver);
}; };
template <typename FImpl, int nBasis> template <typename FImpl, typename Pack>
class TA2AVectors : public Module<A2AVectorsPar> class TA2AVectors : public Module<A2AVectorsPar>
{ {
public: public:
FERM_TYPE_ALIASES(FImpl,); FERM_TYPE_ALIASES(FImpl,);
SOLVER_TYPE_ALIASES(FImpl,); SOLVER_TYPE_ALIASES(FImpl,);
typedef A2AModesSchurDiagTwo<FermionField, FMat, Solver> A2ABase;
typedef FermionEigenPack<FImpl> EPack; public:
typedef CoarseFermionEigenPack<FImpl, nBasis> CoarseEPack;
typedef A2AModesSchurDiagTwo<typename FImpl::FermionField, FMat, Solver> A2ABase;
public:
// constructor // constructor
TA2AVectors(const std::string name); TA2AVectors(const std::string name);
// destructor // destructor
@ -53,46 +48,37 @@ class TA2AVectors : public Module<A2AVectorsPar>
virtual void setup(void); virtual void setup(void);
// execution // execution
virtual void execute(void); virtual void execute(void);
private:
unsigned int Ls_;
std::string className_;
}; };
MODULE_REGISTER_TMP(A2AVectors, ARG(TA2AVectors<FIMPL, HADRONS_DEFAULT_LANCZOS_NBASIS>), MSolver); MODULE_REGISTER_TMP(A2AVectors,
MODULE_REGISTER_TMP(ZA2AVectors, ARG(TA2AVectors<ZFIMPL, HADRONS_DEFAULT_LANCZOS_NBASIS>), MSolver); ARG(TA2AVectors<FIMPL, FermionEigenPack<FIMPL>>), MSolver);
MODULE_REGISTER_TMP(ZA2AVectors,
ARG(TA2AVectors<ZFIMPL, FermionEigenPack<ZFIMPL>>), MSolver);
/****************************************************************************** /******************************************************************************
* TA2AVectors implementation * * TA2AVectors implementation *
******************************************************************************/ ******************************************************************************/
// constructor ///////////////////////////////////////////////////////////////// // constructor /////////////////////////////////////////////////////////////////
template <typename FImpl, int nBasis> template <typename FImpl, typename Pack>
TA2AVectors<FImpl, nBasis>::TA2AVectors(const std::string name) TA2AVectors<FImpl, Pack>::TA2AVectors(const std::string name)
: Module<A2AVectorsPar>(name) : Module<A2AVectorsPar>(name)
, className_ (name + "_class")
{} {}
// dependencies/products /////////////////////////////////////////////////////// // dependencies/products ///////////////////////////////////////////////////////
template <typename FImpl, int nBasis> template <typename FImpl, typename Pack>
std::vector<std::string> TA2AVectors<FImpl, nBasis>::getInput(void) std::vector<std::string> TA2AVectors<FImpl, Pack>::getInput(void)
{ {
int Nl = par().Nl; int Nl = par().Nl;
std::string sub_string = ""; std::string sub_string = "";
if (Nl > 0) sub_string = "_subtract"; if (Nl > 0) sub_string = "_subtract";
std::vector<std::string> in = {par().solver + sub_string};
int n = par().sources.size(); std::vector<std::string> in = {par().solver + sub_string, par().noise};
for (unsigned int t = 0; t < n; t += 1)
{
in.push_back(par().sources[t]);
}
return in; return in;
} }
template <typename FImpl, int nBasis> template <typename FImpl, typename Pack>
std::vector<std::string> TA2AVectors<FImpl, nBasis>::getReference(void) std::vector<std::string> TA2AVectors<FImpl, Pack>::getReference(void)
{ {
std::vector<std::string> ref = {par().action}; std::vector<std::string> ref = {par().action};
@ -104,22 +90,21 @@ std::vector<std::string> TA2AVectors<FImpl, nBasis>::getReference(void)
return ref; return ref;
} }
template <typename FImpl, int nBasis> template <typename FImpl, typename Pack>
std::vector<std::string> TA2AVectors<FImpl, nBasis>::getOutput(void) std::vector<std::string> TA2AVectors<FImpl, Pack>::getOutput(void)
{ {
std::vector<std::string> out = {getName(), className_}; std::vector<std::string> out = {getName()};
return out; return out;
} }
// setup /////////////////////////////////////////////////////////////////////// // setup ///////////////////////////////////////////////////////////////////////
template <typename FImpl, int nBasis> template <typename FImpl, typename Pack>
void TA2AVectors<FImpl, nBasis>::setup(void) void TA2AVectors<FImpl, Pack>::setup(void)
{ {
int N = par().N;
int Nl = par().Nl; int Nl = par().Nl;
int Nh = N - Nl; bool return_5d = par().return_5d;
bool return_5d = par().return_5d; auto &noise = envGet(DilutedNoise<FImpl>, par().noise);
int Ls; int Ls;
std::string sub_string = ""; std::string sub_string = "";
@ -139,94 +124,80 @@ void TA2AVectors<FImpl, nBasis>::setup(void)
if (Nl > 0) if (Nl > 0)
{ {
// Low modes // Low modes
auto &epack = envGet(EPack, par().eigenPack); auto &epack = envGet(Pack, par().eigenPack);
LOG(Message) << "Creating a2a vectors " << getName() << LOG(Message) << "Creating a2a vectors " << getName() <<
" using eigenpack '" << par().eigenPack << "' (" " using eigenpack '" << par().eigenPack << "' ("
<< epack.evec.size() << " modes)" << << epack.evec.size() << " modes)" <<
" and " << Nh << " high modes." << std::endl; " and " << noise.size() << " high modes." << std::endl;
evec = &epack.evec; evec = &epack.evec;
eval = &epack.eval; eval = &epack.eval;
} }
else else
{ {
LOG(Message) << "Creating a2a vectors " << getName() << LOG(Message) << "Creating a2a vectors " << getName() <<
" using " << Nh << " high modes only." << std::endl; " using " << noise.size() << " high modes only." << std::endl;
} }
envCreate(A2ABase, className_, Ls, envCreate(A2ABase, getName(), Ls, evec, eval, action, solver, Nl, noise.size(),
evec, eval, return_5d);
action,
solver,
Nl, Nh,
return_5d);
} }
// execution /////////////////////////////////////////////////////////////////// // execution ///////////////////////////////////////////////////////////////////
template <typename FImpl, int nBasis> template <typename FImpl, typename Pack>
void TA2AVectors<FImpl, nBasis>::execute(void) void TA2AVectors<FImpl, Pack>::execute(void)
{ {
auto &action = envGet(FMat, par().action); auto &action = envGet(FMat, par().action);
auto &noise = envGet(DilutedNoise<FImpl>, par().noise);
int Nt = env().getDim(Tp); int Ls;
int Nc = FImpl::Dimension;
int Ls_;
int Nl = par().Nl; int Nl = par().Nl;
std::string sub_string = ""; std::string sub_string = "";
if (Nl > 0) sub_string = "_subtract"; if (Nl > 0) sub_string = "_subtract";
Ls_ = env().getObjectLs(par().solver + sub_string); Ls = env().getObjectLs(par().solver + sub_string);
auto &a2areturn = envGet(A2ABase, className_); auto &a2areturn = envGet(A2ABase, getName());
// High modes // High modes
auto sources = par().sources;
int Nsrc = par().sources.size();
envGetTmp(FermionField, ferm_src); envGetTmp(FermionField, ferm_src);
envGetTmp(FermionField, unphys_ferm); envGetTmp(FermionField, unphys_ferm);
envGetTmp(FermionField, tmp); envGetTmp(FermionField, tmp);
for (unsigned int i = 0; i < noise.size(); i++)
int N_count = 0; {
for (unsigned int s = 0; s < Ns; ++s) LOG(Message) << "A2A src for noise vector " << i << std::endl;
for (unsigned int c = 0; c < Nc; ++c) // source conversion for 4D sources
for (unsigned int T = 0; T < Nsrc; T++) if (!env().isObject5d(par().noise))
{ {
auto &prop_src = envGet(PropagatorField, sources[T]); if (Ls == 1)
LOG(Message) << "A2A src for s = " << s << " , c = " << c << ", T = " << T << std::endl;
// source conversion for 4D sources
if (!env().isObject5d(sources[T]))
{ {
if (Ls_ == 1) ferm_src = noise[i];
{ tmp = ferm_src;
PropToFerm<FImpl>(ferm_src, prop_src, s, c);
tmp = ferm_src;
}
else
{
PropToFerm<FImpl>(tmp, prop_src, s, c);
action.ImportPhysicalFermionSource(tmp, ferm_src);
action.ImportUnphysicalFermion(tmp, unphys_ferm);
}
} }
// source conversion for 5D sources
else else
{ {
if (Ls_ != env().getObjectLs(sources[T])) tmp = noise[i];
{ action.ImportPhysicalFermionSource(noise[i], ferm_src);
HADRONS_ERROR(Size, "Ls mismatch between quark action and source"); action.ImportUnphysicalFermion(noise[i], unphys_ferm);
}
else
{
PropToFerm<FImpl>(ferm_src, prop_src, s, c);
action.ExportPhysicalFermionSolution(ferm_src, tmp);
unphys_ferm = ferm_src;
}
} }
LOG(Message) << "a2areturn.high_modes Ncount = " << N_count << std::endl;
a2areturn.high_modes(ferm_src, unphys_ferm, tmp, N_count);
N_count++;
} }
// source conversion for 5D sources
else
{
if (Ls != env().getObjectLs(par().noise))
{
HADRONS_ERROR(Size, "Ls mismatch between quark action and source");
}
else
{
ferm_src = noise[i];
action.ExportPhysicalFermionSolution(ferm_src, tmp);
unphys_ferm = ferm_src;
}
}
LOG(Message) << "solveHighMode i = " << i << std::endl;
a2areturn.high_modes(ferm_src, unphys_ferm, tmp, i);
}
} }
END_MODULE_NAMESPACE END_MODULE_NAMESPACE