1
0
mirror of https://github.com/paboyle/Grid.git synced 2025-04-04 19:25:56 +01:00

Improvement to sequential conserved current insertion tests

This commit is contained in:
Lanny91 2017-05-12 16:30:43 +01:00
parent c2010f21ab
commit 34332fe393
3 changed files with 131 additions and 101 deletions

View File

@ -41,7 +41,6 @@ BEGIN_HADRONS_NAMESPACE
* options:
- q: point source propagator, 5D if available (string)
- q4d: 4D point source propagator, duplicate of q if q is 4D (string)
- qSeq: result of sequential insertion of conserved current using q (string)
- action: action used for computation of q (string)
- origin: string giving point source origin of q (string)
@ -60,7 +59,6 @@ class TestSeqConservedPar: Serializable
public:
GRID_SERIALIZABLE_CLASS_MEMBERS(TestSeqConservedPar,
std::string, q,
std::string, q4d,
std::string, qSeq,
std::string, action,
std::string, origin,
@ -103,8 +101,7 @@ TTestSeqConserved<FImpl>::TTestSeqConserved(const std::string name)
template <typename FImpl>
std::vector<std::string> TTestSeqConserved<FImpl>::getInput(void)
{
std::vector<std::string> in = {par().q, par().q4d,
par().qSeq, par().action};
std::vector<std::string> in = {par().q, par().qSeq, par().action};
return in;
}
@ -121,7 +118,11 @@ std::vector<std::string> TTestSeqConserved<FImpl>::getOutput(void)
template <typename FImpl>
void TTestSeqConserved<FImpl>::setup(void)
{
auto Ls = env().getObjectLs(par().q);
if (Ls != env().getObjectLs(par().action))
{
HADRON_ERROR("Ls mismatch between quark action and propagator");
}
}
// execution ///////////////////////////////////////////////////////////////////
@ -130,33 +131,43 @@ void TTestSeqConserved<FImpl>::execute(void)
{
PropagatorField tmp(env().getGrid());
PropagatorField &q = *env().template getObject<PropagatorField>(par().q);
PropagatorField &q4d = *env().template getObject<PropagatorField>(par().q4d);
PropagatorField &qSeq = *env().template getObject<PropagatorField>(par().qSeq);
FMat &act = *(env().template getObject<FMat>(par().action));
Gamma g5(Gamma::Algebra::Gamma5);
Gamma::Algebra gA = (par().curr == Current::Axial) ?
Gamma::Algebra::Gamma5 :
Gamma::Algebra::Identity;
Gamma g(gA);
SitePropagator qSite;
LatticeComplex c(env().getGrid());
Complex seq_res, check_res;
std::vector<TComplex> check_buf;
Complex test_S, test_V, check_S, check_V;
std::vector<SitePropagator> check_buf;
// Check sequential insertion of current gives same result as conserved
// current sink upon contraction. Assume q uses a point source.
std::vector<int> siteCoord;
siteCoord = strToVec<int>(par().origin);
peekSite(qSite, q, siteCoord);
seq_res = trace(g5*qSite);
peekSite(qSite, qSeq, siteCoord);
test_S = trace(qSite*g);
test_V = trace(qSite*g*Gamma::gmu[par().mu]);
act.ContractConservedCurrent(q, q, tmp, par().curr, par().mu);
c = trace(tmp);
sliceSum(c, check_buf, Tp);
check_res = TensorRemove(check_buf[par().t_J]);
sliceSum(tmp, check_buf, Tp);
check_S = TensorRemove(trace(check_buf[par().t_J]*g));
check_V = TensorRemove(trace(check_buf[par().t_J]*g*Gamma::gmu[par().mu]));
LOG(Message) << "Test S = " << abs(test_S) << std::endl;
LOG(Message) << "Test V = " << abs(test_V) << std::endl;
LOG(Message) << "Check S = " << abs(check_S) << std::endl;
LOG(Message) << "Check V = " << abs(check_V) << std::endl;
// Check difference = 0
check_res -= seq_res;
check_S -= test_S;
check_V -= test_V;
LOG(Message) << "Consistency check for sequential conserved "
<< par().curr << " current insertion = " << abs(check_res)
<< std::endl;
<< par().curr << " current insertion: " << std::endl;
LOG(Message) << "Check S = " << abs(check_S) << std::endl;
LOG(Message) << "Check V = " << abs(check_V) << std::endl;
}
END_MODULE_NAMESPACE

View File

@ -538,39 +538,30 @@ inline void makeWITest(Application &application, std::string &modName,
}
/*******************************************************************************
* Name: makeSeqTest
* Purpose: Create module to test sequential insertion of conserved current
* and add to application module.
* Name: makeSeqCurrComparison
* Purpose: Create module to compare sequential insertion of conserved current
* against sink contraction and add to application module.
* Parameters: application - main application that stores modules.
* modName - name of module to create.
* propName - 4D quark propagator.
* seqProp - 4D quark propagator with sequential insertion of
* propName - quark propagator (point source), 5D if available.
* seqName - 4D quark propagator with sequential insertion of
* conserved current.
* actionName - action used to compute quark propagators.
* origin - origin of point source propagator.
* t_J - time at which sequential current is inserted.
* mu - Lorentz index of sequential current.
* curr - type of conserved current inserted.
* Ls - length of 5th dimension (default = 1).
* Returns: None.
******************************************************************************/
inline void makeSeqTest(Application &application, std::string &modName,
std::string &propName, std::string &seqName,
std::string &actionName, std::string &origin,
unsigned int t_J, unsigned int mu, Current curr,
unsigned int Ls = 1)
inline void makeSeqCurrComparison(Application &application, std::string &modName,
std::string &propName, std::string &seqName,
std::string &actionName, std::string &origin,
unsigned int t_J, unsigned int mu, Current curr)
{
if (!(Environment::getInstance().hasModule(modName)))
{
MUtilities::TestSeqConserved::Par seqPar;
if (Ls > 1)
{
seqPar.q = LABEL_5D(propName);
}
else
{
seqPar.q = propName;
}
seqPar.q4d = propName;
seqPar.q = propName;
seqPar.qSeq = seqName;
seqPar.action = actionName;
seqPar.origin = origin;

View File

@ -30,84 +30,112 @@
using namespace Grid;
using namespace Hadrons;
inline void setupSeqCurrTests(Application &application, std::string modStem,
std::string &pointProp, std::string &seqStem,
std::string &actionName, std::string &solverName,
std::string &origin, Current curr,
unsigned int t_J, unsigned int mu,
unsigned int Ls = 1)
{
std::string modName = ADD_INDEX(modStem, mu);
std::string seqProp = ADD_INDEX(seqStem, mu);
std::string seqSrc = seqProp + "_src";
// 5D actions require 5D propagator as input for conserved current
// insertions.
std::string propIn;
if (Ls > 1)
{
propIn = LABEL_5D(pointProp);
}
else
{
propIn = pointProp;
}
makeConservedSequentialSource(application, seqSrc, propIn,
actionName, t_J, curr, mu);
makePropagator(application, seqProp, seqSrc, solverName);
makeSeqCurrComparison(application, modName, propIn, seqProp,
actionName, origin, t_J, mu, curr);
}
inline void setupWardIdentityTests(Application &application,
std::string &actionName,
double mass,
unsigned int Ls = 1,
bool perform_axial_tests = false)
{
// solver
std::string solverName = actionName + "_CG";
makeRBPrecCGSolver(application, solverName, actionName);
unsigned int nt = GridDefaultLatt()[Tp];
unsigned int t_J = nt/2;
/***************************************************************************
* Conserved current sink contractions: use a single point propagator for
* the Ward Identity test.
**************************************************************************/
std::string pointProp = actionName + "_q_0";
std::string origin = "0 0 0 0";
std::string modName = actionName + " Ward Identity Test";
MAKE_POINT_PROP(origin, pointProp, solverName);
makeWITest(application, modName, pointProp, actionName, mass, Ls);
/***************************************************************************
* Conserved current tests with sequential insertion of vector/axial
* current. If above Ward Identity passes, sufficient to test sequential
* insertion of conserved current agrees with contracted version.
**************************************************************************/
// Compare sequential insertion to contraction. Should be enough to perform
// for time and one space component.
std::string seqStem = ADD_INDEX(pointProp + "seq_V", t_J);
std::string modStem = actionName + " Vector Sequential Test mu";
setupSeqCurrTests(application, modStem, pointProp, seqStem, actionName,
solverName, origin, Current::Vector, t_J, Tp, Ls);
setupSeqCurrTests(application, modStem, pointProp, seqStem, actionName,
solverName, origin, Current::Vector, t_J, Xp, Ls);
// Perform axial tests only if partially-conserved axial current exists for
// the action.
if (perform_axial_tests)
{
seqStem = ADD_INDEX(pointProp + "seq_A", t_J);
modStem = actionName + " Axial Sequential Test mu";
setupSeqCurrTests(application, modStem, pointProp, seqStem, actionName,
solverName, origin, Current::Axial, t_J, Tp, Ls);
setupSeqCurrTests(application, modStem, pointProp, seqStem, actionName,
solverName, origin, Current::Axial, t_J, Xp, Ls);
}
}
int main(int argc, char *argv[])
{
// initialization //////////////////////////////////////////////////////////
Grid_init(&argc, &argv);
HadronsLogError.Active(GridLogError.isActive());
HadronsLogWarning.Active(GridLogWarning.isActive());
HadronsLogMessage.Active(GridLogMessage.isActive());
HadronsLogIterative.Active(GridLogIterative.isActive());
HadronsLogDebug.Active(GridLogDebug.isActive());
LOG(Message) << "Grid initialized" << std::endl;
HADRONS_DEFAULT_INIT;
// run setup ///////////////////////////////////////////////////////////////
Application application;
unsigned int nt = GridDefaultLatt()[Tp];
double mass = 0.04;
double M5 = 1.8;
unsigned int Ls = 12;
// global parameters
Application::GlobalPar globalPar;
globalPar.trajCounter.start = 1500;
globalPar.trajCounter.end = 1520;
globalPar.trajCounter.step = 20;
globalPar.seed = "1 2 3 4";
globalPar.genetic.maxGen = 1000;
globalPar.genetic.maxCstGen = 200;
globalPar.genetic.popSize = 20;
globalPar.genetic.mutationRate = .1;
application.setPar(globalPar);
HADRONS_DEFAULT_GLOBALS(application);
// gauge field
application.createModule<MGauge::Unit>("gauge");
std::string gaugeField = "gauge";
application.createModule<MGauge::Unit>(gaugeField);
// action
// Setup each action and the conserved current tests relevant to it.
std::string actionName = "DWF";
MAction::DWF::Par actionPar;
actionPar.gauge = "gauge";
actionPar.Ls = Ls;
actionPar.M5 = 1.8;
actionPar.mass = mass;
application.createModule<MAction::DWF>(actionName, actionPar);
makeDWFAction(application, actionName, gaugeField, mass, M5, Ls);
setupWardIdentityTests(application, actionName, mass, Ls, true);
// solver
std::string solverName = "CG";
MSolver::RBPrecCG::Par solverPar;
solverPar.action = actionName;
solverPar.residual = 1.0e-8;
application.createModule<MSolver::RBPrecCG>(solverName,
solverPar);
// Conserved current sink contractions: use a single point propagator.
std::string pointProp = "q_0";
std::string pos = "0 0 0 0";
std::string modName = "Ward Identity Test";
MAKE_POINT_PROP(pos, pointProp, solverName);
makeWITest(application, modName, pointProp, actionName, mass, Ls);
// Conserved current contractions with sequential insertion of vector/axial
// current.
std::string mom = ZERO_MOM;
unsigned int t_J = nt/2;
std::string seqPropA = ADD_INDEX(pointProp + "_seq_A", t_J);
std::string seqPropV = ADD_INDEX(pointProp + "_seq_V", t_J);
std::string seqSrcA = seqPropA + "_src";
std::string seqSrcV = seqPropV + "_src";
std::string point5d = LABEL_5D(pointProp);
makeConservedSequentialSource(application, seqSrcA, point5d,
actionName, t_J, Current::Axial, Tp, mom);
makePropagator(application, seqPropA, seqSrcA, solverName);
makeConservedSequentialSource(application, seqSrcV, point5d,
actionName, t_J, Current::Vector, Tp, mom);
makePropagator(application, seqPropV, seqSrcV, solverName);
std::string modNameA = "Axial Sequential Test";
std::string modNameV = "Vector Sequential Test";
makeSeqTest(application, modNameA, pointProp, seqPropA,
actionName, pos, t_J, Tp, Current::Axial, Ls);
makeSeqTest(application, modNameV, pointProp, seqPropV,
actionName, pos, t_J, Tp, Current::Vector, Ls);
actionName = "Wilson";
makeWilsonAction(application, actionName, gaugeField, mass);
setupWardIdentityTests(application, actionName, mass);
// execution
application.saveParameterFile("ConservedCurrentTest.xml");