From cd8f0eeba183c087fbbaadf3addb3761adb42ab5 Mon Sep 17 00:00:00 2001 From: Ryan Hill Date: Fri, 28 Aug 2026 00:18:04 +0100 Subject: [PATCH 1/6] Add required RbcUkqcd presets for C0 IB-disc. diagrams --- RbcUkqcd.hpp | 240 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 240 insertions(+) diff --git a/RbcUkqcd.hpp b/RbcUkqcd.hpp index 5439b7b..a4a5f00 100644 --- a/RbcUkqcd.hpp +++ b/RbcUkqcd.hpp @@ -40,6 +40,9 @@ struct RbcUkqcd unsigned int L, T, Ls; }; + inline static constexpr hadpresets::RbcUkqcd::EnsembleParameters testPar{0.005, 0.0362, 1., 2., + 8, 8, 4}; + // C0 inline static constexpr EnsembleParameters c0UnitaryPar{0.00078, 0.0362, 1.8, 2., 48, 96, 24}; inline static constexpr unsigned int c0ZMobiusLs = 10; @@ -79,6 +82,11 @@ struct RbcUkqcd inline static constexpr DeflationParameters c1m20DeflPar{3.0e-03, 5.5, 101, 200, 220, 230}; inline static constexpr DeflationParameters c1m32DeflPar{5.0e-05, 5.5, 101, 100, 110, 120}; + // Light solvers: no deflation + static inline void addLightSolver(Application &app, const RbcUkqcd::EnsembleParameters &par, + const std::string solverName, const std::string gaugeName, + const std::string boundary, const double residual); + // Light solvers: load deflation from disk static inline void addC0LightZMobiusLCDSolver(Application &app, const std::string solverName, const std::string gaugeName, @@ -87,6 +95,12 @@ struct RbcUkqcd const std::string boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void + addC0LightMADWFLCDSolver(Grid::Hadrons::Application &app, const std::string solverName, + const std::string gaugeName, const std::string gaugeTransform, + const std::string eigenpackPath, const std::string boundary = "1 1 1 -1", + const double residual = 1.0e-8); + static inline void addM0LightLCDSolver(Application &app, const std::string solverName, const std::string gaugeName, const std::string gaugeTransform, const std::string eigenpackPath, @@ -125,6 +139,10 @@ struct RbcUkqcd const std::string gaugeName, const std::string boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addC0StrangeZMobiusSolver(Grid::Hadrons::Application &app, + const std::string solverName, + const std::string gaugeName, + const std::string boundary, const double residual); static inline void addM0StrangeSolver(Application &app, const std::string solverName, const std::string gaugeName, const std::string boundary = "1 1 1 1", @@ -147,6 +165,10 @@ struct RbcUkqcd const double residual = 1.0e-8); // Charm solvers, mass is a free parameter + static inline void addUnsmearedCharmSolver(Application &app, const std::string solverName, + const std::string gaugeName, const double mass, + const int Ls, const double M5, + const std::string boundary, const double residual); static inline void addC0CharmSolver(Application &app, const std::string solverName, const std::string gaugeName, const double mass, const std::string boundary = "1 1 1 -1", @@ -159,6 +181,46 @@ struct RbcUkqcd // Implementations ///////////////////////////////////////////////////////////////////////////////// +// Light +void RbcUkqcd::addLightSolver(Application &app, const RbcUkqcd::EnsembleParameters &par, + const std::string solverName, const std::string gaugeName, + const std::string boundary, const double residual) +{ + const std::string prefix = solverName; + + // Gauge field FP32 cast + MUtilities::GaugeSinglePrecisionCast::Par gaugeCastPar; + + gaugeCastPar.field = gaugeName; + app.createModule(prefix + "_gauge_fp32", gaugeCastPar); + + // Scaled DWF action + FP32 version + MAction::ScaledDWF::Par actionPar; + + actionPar.gauge = gaugeName; + actionPar.Ls = par.Ls; + actionPar.M5 = par.M5; + actionPar.mass = par.ml; + actionPar.scale = par.scale; + actionPar.boundary = boundary; + actionPar.twist = "0. 0. 0. 0."; + app.createModule(prefix + "_dwf", actionPar); + actionPar.gauge = prefix + "_gauge_fp32"; + app.createModule(prefix + "_dwf_fp32", actionPar); + + // Mixed-precision red-black preconditionned CG + MSolver::MixedPrecisionRBPrecCG::Par solverPar; + + solverPar.innerAction = prefix + "_dwf_fp32"; + solverPar.outerAction = prefix + "_dwf"; + solverPar.maxInnerIteration = 30000; + solverPar.maxOuterIteration = 100; + solverPar.residual = residual; + solverPar.innerGuesser = ""; + solverPar.outerGuesser = ""; + app.createModule(solverName, solverPar); +} + // Light C0 (load deflation from disk) void RbcUkqcd::addC0LightZMobiusLCDSolver(Application &app, const std::string solverName, const std::string gaugeName, @@ -236,6 +298,139 @@ void RbcUkqcd::addC0LightZMobiusLCDSolver(Application &app, const std::string so app.createModule(solverName, solverPar); } +void RbcUkqcd::addC0LightMADWFLCDSolver(Grid::Hadrons::Application &app, + const std::string solverName, const std::string gaugeName, + const std::string gaugeTransform, + const std::string eigenpackPath, const std::string boundary, + const double residual) +{ + using namespace Grid; + using namespace Hadrons; + using namespace hadpresets; + + const std::string prefix = solverName; + const bool gaugeFixed = !gaugeTransform.empty(); + + // Gauge field & transform FP32 cast + MUtilities::GaugeSinglePrecisionCast::Par gaugeCastPar; + + gaugeCastPar.field = gaugeName; + app.createModule(prefix + "_gauge_fp32", gaugeCastPar); + if (gaugeFixed) + { + MUtilities::ColourMatrixSinglePrecisionCast::Par transformCastPar; + + transformCastPar.field = gaugeTransform; + app.createModule(prefix + "_gaugeTransform_fp32", + transformCastPar); + } + + // Scaled DWF action + MAction::ScaledDWF::Par actionPar; + actionPar.gauge = gaugeName; + actionPar.Ls = RbcUkqcd::c0UnitaryPar.Ls; + actionPar.M5 = RbcUkqcd::c0UnitaryPar.M5; + actionPar.mass = RbcUkqcd::c0UnitaryPar.ml; + actionPar.scale = 2.; + actionPar.boundary = boundary; + actionPar.twist = "0. 0. 0. 0."; + app.createModule(prefix + "_dwf", actionPar); + + // Single-prec ZMobius Action + MAction::ZMobiusDWFF::Par actionFPar; + actionFPar.Ls = RbcUkqcd::c0LCDPar.Ls; + actionFPar.M5 = RbcUkqcd::c0LCDPar.M5; + actionFPar.mass = RbcUkqcd::c0LCDPar.ml; + actionFPar.b = 1.; + actionFPar.c = 0.; + actionFPar.omega = std::vector>(RbcUkqcd::c0ZMobiusOmega.begin(), + RbcUkqcd::c0ZMobiusOmega.end()); + actionFPar.boundary = boundary; + actionFPar.twist = "0. 0. 0. 0."; + actionFPar.gauge = prefix + "_gauge_fp32"; + app.createModule(prefix + "_dwf_fp32", actionFPar); + + // Compressed eigenpack + MIO::LoadCoarseFermionEigenPack250F::Par epPar; + + epPar.filestem = eigenpackPath; + epPar.multiFile = true; + epPar.redBlack = true; + epPar.sizeFine = 250; + epPar.sizeCoarse = 2000; + epPar.Ls = 10; + epPar.blockSize = "4 4 4 3 10"; + epPar.orthogonalise = gaugeFixed; + epPar.gaugeXform = gaugeFixed ? (prefix + "_gaugeTransform_fp32") : ""; + app.createModule(prefix + "_epack", epPar); + + // Inner guesser + MGuesser::CoarseDeflation250F::Par iguessPar; + + iguessPar.eigenPack = prefix + "_epack"; + iguessPar.size = 2000; + app.createModule(prefix + "_iguesser", iguessPar); + + // ZMobius-MADWF red-black preconditionned CG + MSolver::ZMADWFMixedPrecCG::Par solverPar; + solverPar.innerAction = prefix + "_dwf_fp32"; + solverPar.outerAction = prefix + "_dwf"; + solverPar.maxInnerIteration = 10000; + solverPar.maxPVIteration = 1000; + solverPar.maxOuterIteration = 100; + solverPar.innerResidual = 1e-5; + solverPar.PVResidual = 1e-5; + solverPar.outerResidual = residual; + solverPar.guesser = prefix + "_iguesser"; + app.createModule(solverName, solverPar); +} + +void RbcUkqcd::addC0StrangeZMobiusSolver(Grid::Hadrons::Application &app, + const std::string solverName, const std::string gaugeName, + const std::string boundary, const double residual) +{ + using namespace Grid; + using namespace Hadrons; + using namespace hadpresets; + + const std::string prefix = solverName; + + // Gauge field & transform FP32 cast + MUtilities::GaugeSinglePrecisionCast::Par gaugeCastPar; + + gaugeCastPar.field = gaugeName; + app.createModule(prefix + "_gauge_fp32", gaugeCastPar); + + // Scaled DWF action + FP32 version + MAction::ZMobiusDWF::Par actionPar; + + actionPar.gauge = gaugeName; + actionPar.Ls = RbcUkqcd::c0LCDPar.Ls; + actionPar.M5 = RbcUkqcd::c0LCDPar.M5; + actionPar.mass = RbcUkqcd::c0LCDPar.ms; + actionPar.b = 1.; + actionPar.c = 0.; + actionPar.omega = std::vector>(RbcUkqcd::c0ZMobiusOmega.begin(), + RbcUkqcd::c0ZMobiusOmega.end()); + actionPar.boundary = boundary; + actionPar.twist = "0. 0. 0. 0."; + app.createModule(prefix + "_dwf", actionPar); + actionPar.gauge = prefix + "_gauge_fp32"; + app.createModule(prefix + "_dwf_fp32", actionPar); + + // Batched mixed-precision red-black preconditionned CG + MSolver::ZMixedPrecisionRBPrecCG::Par solverPar; + + solverPar.innerAction = prefix + "_dwf_fp32"; + solverPar.outerAction = prefix + "_dwf"; + solverPar.maxInnerIteration = 10000; + solverPar.maxOuterIteration = 100; + solverPar.residual = residual; + solverPar.innerGuesser = ""; + solverPar.outerGuesser = ""; + app.createModule(solverName, solverPar); +} + // Light M0 (load deflation from disk) void RbcUkqcd::addM0LightLCDSolver(Application &app, const std::string solverName, const std::string gaugeName, const std::string gaugeTransform, @@ -497,6 +692,51 @@ void RbcUkqcd::addC1M32StrangeSolver(Application &app, const std::string solverN RbcUkqcd::addStrangeSolver(app, RbcUkqcd::c1m32IRLPar, solverName, gaugeName, boundary, residual); } +// Charm +void RbcUkqcd::addUnsmearedCharmSolver(Grid::Hadrons::Application &app, + const std::string solverName, const std::string gaugeName, + const double mass, const int Ls, const double M5, + const std::string boundary, const double residual) +{ + using namespace Grid; + using namespace Hadrons; + + const std::string prefix = solverName; + + // Gauge field & transform FP32 cast + MUtilities::GaugeSinglePrecisionCast::Par gaugeCastPar; + + gaugeCastPar.field = gaugeName; + app.createModule(prefix + "_gauge_fp32", gaugeCastPar); + + // Scaled DWF action + FP32 version + MAction::ScaledDWF::Par actionPar; + + actionPar.gauge = gaugeName; + actionPar.Ls = Ls; + actionPar.M5 = M5; + actionPar.mass = mass; + actionPar.scale = 2.; + + actionPar.boundary = boundary; + actionPar.twist = "0. 0. 0. 0."; + app.createModule(prefix + "_dwf", actionPar); + actionPar.gauge = prefix + "_gauge_fp32"; + app.createModule(prefix + "_dwf_fp32", actionPar); + + // Batched mixed-precision red-black preconditionned CG + MSolver::MixedPrecisionRBPrecCG::Par solverPar; + + solverPar.innerAction = prefix + "_dwf_fp32"; + solverPar.outerAction = prefix + "_dwf"; + solverPar.maxInnerIteration = 10000; + solverPar.maxOuterIteration = 100; + solverPar.residual = residual; + solverPar.innerGuesser = ""; + solverPar.outerGuesser = ""; + app.createModule(solverName, solverPar); +} + // Charm C0 void RbcUkqcd::addC0CharmSolver(Application &app, const std::string solverName, const std::string gaugeName, const double mass, From bd253247c9d402ca4689774e35f61e1d1e43e678 Mon Sep 17 00:00:00 2001 From: Ryan Hill Date: Fri, 28 Aug 2026 00:18:18 +0100 Subject: [PATCH 2/6] Add Jlqcd presets --- Jlqcd.hpp | 206 ++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 206 insertions(+) create mode 100644 Jlqcd.hpp diff --git a/Jlqcd.hpp b/Jlqcd.hpp new file mode 100644 index 0000000..ffd5bb6 --- /dev/null +++ b/Jlqcd.hpp @@ -0,0 +1,206 @@ +/* + * Copyright (C) 2026 + * + * Hadrons presets for the JLQCD scaled-domain-wall ensembles used by iblep. + */ + +#pragma once + +#include +#include + +#include +#include +#include +#include + +namespace hadpresets +{ +using namespace Grid::Hadrons; + +struct Jlqcd +{ + struct EnsembleParameters + { + double ml, ms, M5, scale; + unsigned int L, T, Ls; + }; + + inline static constexpr EnsembleParameters fUd3SaPar{0.003, 0.015, 1., 2., 64, 128, 8}; + inline static constexpr EnsembleParameters mUd3SaPar{0.0042, 0.025, 1., 2., 48, 96, 8}; + inline static constexpr EnsembleParameters mUd3SbPar{0.0042, 0.018, 1., 2., 48, 96, 8}; + inline static constexpr EnsembleParameters cUd2SaPar{0.0035, 0.040, 1., 2., 32, 64, 12}; + inline static constexpr EnsembleParameters cUd2SaLPar{0.0035, 0.040, 1., 2., 48, 96, 12}; + inline static constexpr EnsembleParameters cUd3SaPar{0.007, 0.040, 1., 2., 32, 64, 12}; + inline static constexpr EnsembleParameters cUd3SbPar{0.007, 0.030, 1., 2., 32, 64, 12}; + inline static constexpr EnsembleParameters testPar{0.010, 0.030, 1., 2., 8, 16, 8}; + + static inline const EnsembleParameters &ensemble(const std::string &name); + + static inline void addMixedPrecisionSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, + const std::string &gaugeName, const double mass, + const std::string &twist, + const std::string &boundary = "1 1 1 -1", + const double residual = 1.0e-8); + + static inline void addNoFailSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary = "1 1 1 -1", + const double residual = 1.0e-20); + + static inline void addLightSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary = "1 1 1 -1", + const double residual = 1.0e-8); + + static inline void addLightSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", + const double residual = 1.0e-8); + + static inline void addStrangeSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary = "1 1 1 -1", + const double residual = 1.0e-8); + + static inline void addStrangeSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", + const double residual = 1.0e-8); + + static inline void addHeavySolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary = "1 1 1 -1", + const double residual = 1.0e-20); +}; + +inline const Jlqcd::EnsembleParameters &Jlqcd::ensemble(const std::string &name) +{ + if (name == "fUd3Sa") + return fUd3SaPar; + if (name == "mUd3Sa") + return mUd3SaPar; + if (name == "mUd3Sb") + return mUd3SbPar; + if (name == "cUd2Sa") + return cUd2SaPar; + if (name == "cUd2SaL") + return cUd2SaLPar; + if (name == "cUd3Sa") + return cUd3SaPar; + if (name == "cUd3Sb") + return cUd3SbPar; + if (name == "test") + return testPar; + + throw std::invalid_argument("unknown JLQCD ensemble preset '" + std::string(name) + "'"); +} + +inline void Jlqcd::addMixedPrecisionSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, + const std::string &gaugeName, const double mass, + const std::string &twist, const std::string &boundary, + const double residual) +{ + const std::string prefix = solverName; + + MUtilities::GaugeSinglePrecisionCast::Par gaugeCastPar; + gaugeCastPar.field = gaugeName; + app.createModule(prefix + "_gauge_fp32", gaugeCastPar); + + MAction::ScaledDWF::Par actionPar; + actionPar.gauge = gaugeName; + actionPar.Ls = par.Ls; + actionPar.M5 = par.M5; + actionPar.mass = mass; + actionPar.scale = par.scale; + actionPar.boundary = boundary; + actionPar.twist = twist; + app.createModule(prefix + "_dwf", actionPar); + + actionPar.gauge = prefix + "_gauge_fp32"; + app.createModule(prefix + "_dwf_fp32", actionPar); + + MSolver::MixedPrecisionRBPrecCG::Par solverPar; + solverPar.innerAction = prefix + "_dwf_fp32"; + solverPar.outerAction = prefix + "_dwf"; + solverPar.maxInnerIteration = 10000; + solverPar.maxOuterIteration = 100; + solverPar.residual = residual; + solverPar.innerGuesser = ""; + solverPar.outerGuesser = ""; + app.createModule(solverName, solverPar); +} + +inline void Jlqcd::addNoFailSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + const std::string prefix = solverName; + + MAction::ScaledDWF::Par actionPar; + actionPar.gauge = gaugeName; + actionPar.Ls = par.Ls; + actionPar.M5 = par.M5; + actionPar.mass = mass; + actionPar.scale = par.scale; + actionPar.boundary = boundary; + actionPar.twist = twist; + app.createModule(prefix + "_dwf", actionPar); + + MSolver::RBPrecCGNoFail::Par solverPar; + solverPar.action = prefix + "_dwf"; + solverPar.maxIteration = 30000; + solverPar.residual = residual; + solverPar.guesser = ""; + app.createModule(solverName, solverPar); +} + +inline void Jlqcd::addLightSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addMixedPrecisionSolver(app, par, solverName, gaugeName, mass, twist, boundary, residual); +} + +inline void Jlqcd::addLightSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const std::string &twist, const std::string &boundary, + const double residual) +{ + addLightSolver(app, par, solverName, gaugeName, par.ml, twist, boundary, residual); +} + +inline void Jlqcd::addStrangeSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addMixedPrecisionSolver(app, par, solverName, gaugeName, mass, twist, boundary, residual); +} + +inline void Jlqcd::addStrangeSolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const std::string &twist, const std::string &boundary, + const double residual) +{ + addStrangeSolver(app, par, solverName, gaugeName, par.ms, twist, boundary, residual); +} + +inline void Jlqcd::addHeavySolver(Application &app, const EnsembleParameters &par, + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addNoFailSolver(app, par, solverName, gaugeName, mass, twist, boundary, residual); +} +} // namespace hadpresets From 3269f159e3f5294af5076575e2acca31383ea749 Mon Sep 17 00:00:00 2001 From: Ryan Hill Date: Thu, 24 Sep 2026 14:20:14 +0100 Subject: [PATCH 3/6] Update copyright notice --- Jlqcd.hpp | 20 +++++++++++++++++++- 1 file changed, 19 insertions(+), 1 deletion(-) diff --git a/Jlqcd.hpp b/Jlqcd.hpp index ffd5bb6..9fac006 100644 --- a/Jlqcd.hpp +++ b/Jlqcd.hpp @@ -1,7 +1,25 @@ /* * Copyright (C) 2026 * - * Hadrons presets for the JLQCD scaled-domain-wall ensembles used by iblep. + * Author: Antonin Portelli + * Author: Ryan Hill + * Elements based on workflows from Teseo San Jose Perez + * + * Hadrons 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. + * + * Hadrons 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 Hadrons. If not, see . + * + * See the full license in the file "LICENSE" in the top level distribution + * directory. */ #pragma once From d5a757128f6508ad769dfc4e4a8516fde6b02ef4 Mon Sep 17 00:00:00 2001 From: Ryan Hill Date: Tue, 29 Sep 2026 08:44:45 +0100 Subject: [PATCH 4/6] Add C0L support --- RbcUkqcd.hpp | 136 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 136 insertions(+) diff --git a/RbcUkqcd.hpp b/RbcUkqcd.hpp index a4a5f00..85f6f42 100644 --- a/RbcUkqcd.hpp +++ b/RbcUkqcd.hpp @@ -60,6 +60,10 @@ struct RbcUkqcd inline static constexpr EnsembleParameters c0LCDPar{0.00078, 0.0362, 1.8, NAN, 48, 96, c0ZMobiusLs}; + // C0L + inline static constexpr EnsembleParameters c0LUnitaryPar{0.00078, 0.0362, 1.8, 2., 64, 128, 24}; + inline static constexpr EnsembleParameters c0LLCDPar{0.00078, 0.0362, 1.8, 2., 64, 128, 24}; + // M0 inline static constexpr EnsembleParameters m0UnitaryPar{0.000678, 0.02661, 1.8, 2., 64, 128, 12}; inline static constexpr EnsembleParameters m0LCDPar{0.0006203, 0.02661, 1.8, 2., 64, 128, 12}; @@ -101,6 +105,11 @@ struct RbcUkqcd const std::string eigenpackPath, const std::string boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void + addC0LLightLCDSolver(Application &app, const std::string solverName, const std::string gaugeName, + const std::string gaugeTransform, const std::string eigenpackPath, + const std::string boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addM0LightLCDSolver(Application &app, const std::string solverName, const std::string gaugeName, const std::string gaugeTransform, const std::string eigenpackPath, @@ -143,6 +152,10 @@ struct RbcUkqcd const std::string solverName, const std::string gaugeName, const std::string boundary, const double residual); + static inline void addC0LStrangeSolver(Application &app, const std::string solverName, + const std::string gaugeName, + const std::string boundary = "1 1 1 -1", + const double residual = 1.0e-8); static inline void addM0StrangeSolver(Application &app, const std::string solverName, const std::string gaugeName, const std::string boundary = "1 1 1 1", @@ -173,6 +186,10 @@ struct RbcUkqcd const std::string gaugeName, const double mass, const std::string boundary = "1 1 1 -1", const double residual = 1.0e-18); + static inline void addC0LCharmSolver(Application &app, const std::string solverName, + const std::string gaugeName, const double mass, + const std::string boundary = "1 1 1 -1", + const double residual = 1.0e-18); static inline void addM0CharmSolver(Application &app, const std::string solverName, const std::string gaugeName, const double mass, const std::string boundary = "1 1 1 1", @@ -431,6 +448,79 @@ void RbcUkqcd::addC0StrangeZMobiusSolver(Grid::Hadrons::Application &app, app.createModule(solverName, solverPar); } +// Light C0L (load deflation from disk) +void RbcUkqcd::addC0LLightLCDSolver(Application &app, const std::string solverName, + const std::string gaugeName, const std::string gaugeTransform, + const std::string eigenpackPath, const std::string boundary, + const double residual) +{ + const std::string prefix = solverName; + const bool gaugeFixed = !gaugeTransform.empty(); + + // Gauge field FP32 cast + MUtilities::GaugeSinglePrecisionCast::Par gaugeCastPar; + + gaugeCastPar.field = gaugeName; + app.createModule(prefix + "_gauge_fp32", gaugeCastPar); + if (gaugeFixed) + { + MUtilities::ColourMatrixSinglePrecisionCast::Par transformCastPar; + + transformCastPar.field = gaugeTransform; + app.createModule(prefix + "_gaugeTransform_fp32", + transformCastPar); + } + + // Scaled DWF action + FP32 version + MAction::ScaledDWF::Par actionPar; + + actionPar.gauge = gaugeName; + actionPar.Ls = RbcUkqcd::c0LLCDPar.Ls; + actionPar.M5 = RbcUkqcd::c0LLCDPar.M5; + actionPar.mass = RbcUkqcd::c0LLCDPar.ml; + actionPar.scale = RbcUkqcd::c0LLCDPar.scale; + actionPar.boundary = boundary; + actionPar.twist = "0. 0. 0. 0."; + app.createModule(prefix + "_dwf", actionPar); + actionPar.gauge = prefix + "_gauge_fp32"; + app.createModule(prefix + "_dwf_fp32", actionPar); + + // Compressed eigenpack + MIO::LoadCoarseFermionEigenPack200F::Par epPar; + + epPar.filestem = eigenpackPath; + epPar.multiFile = true; + epPar.redBlack = true; + epPar.sizeFine = 200; + epPar.sizeCoarse = 4500; + epPar.Ls = 24; + epPar.blockSize = "4 4 4 4 24"; + epPar.orthogonalise = false; + epPar.gaugeXform = gaugeFixed ? (prefix + "_gaugeTransform_fp32") : ""; + app.createModule(prefix + "_epack", epPar); + + // Inner guesser + MGuesser::CoarseDeflation200F::Par iguessPar; + + iguessPar.eigenPack = prefix + "_epack"; + iguessPar.size = 4500; + app.createModule(prefix + "_iguesser", iguessPar); + + // Batched mixed-precision red-black preconditionned CG + MSolver::MixedPrecisionRBPrecCGBatched::Par solverPar; + + solverPar.innerAction = prefix + "_dwf_fp32"; + solverPar.outerAction = prefix + "_dwf"; + solverPar.maxInnerIteration = 250; + solverPar.maxOuterIteration = 100; + solverPar.maxPatchupIteration = 10000; + solverPar.residual = residual; + solverPar.updateResidual = true; + solverPar.innerGuesser = prefix + "_iguesser"; + solverPar.outerGuesser = ""; + app.createModule(solverName, solverPar); +} + // Light M0 (load deflation from disk) void RbcUkqcd::addM0LightLCDSolver(Application &app, const std::string solverName, const std::string gaugeName, const std::string gaugeTransform, @@ -656,6 +746,14 @@ void RbcUkqcd::addC0StrangeSolver(Application &app, const std::string solverName residual); } +void RbcUkqcd::addC0LStrangeSolver(Application &app, const std::string solverName, + const std::string gaugeName, const std::string boundary, + const double residual) +{ + RbcUkqcd::addStrangeSolver(app, RbcUkqcd::c0LUnitaryPar, solverName, gaugeName, boundary, + residual); +} + void RbcUkqcd::addM0StrangeSolver(Application &app, const std::string solverName, const std::string gaugeName, const std::string boundary, const double residual) @@ -775,6 +873,44 @@ void RbcUkqcd::addC0CharmSolver(Application &app, const std::string solverName, app.createModule(solverName, solverPar); } +// Charm C0L +void RbcUkqcd::addC0LCharmSolver(Application &app, const std::string solverName, + const std::string gaugeName, const double mass, + const std::string boundary, const double residual) +{ + const std::string prefix = solverName; + + // stout smearing + MGauge::StoutSmearing::Par smearPar; + + smearPar.gauge = gaugeName; + smearPar.steps = 3; + smearPar.rho = 0.1; + smearPar.orthogDim = ""; + app.createModule(prefix + "_gauge_3stout", smearPar); + + // Scaled DWF action + FP32 version + MAction::ScaledDWF::Par actionPar; + + actionPar.gauge = prefix + "_gauge_3stout"; + actionPar.Ls = 24; + actionPar.M5 = 1.; + actionPar.mass = mass; + actionPar.scale = 2.; + actionPar.boundary = boundary; + actionPar.twist = "0. 0. 0. 0."; + app.createModule(prefix + "_dwf", actionPar); + + // Red-black preconditionned CG + MSolver::RBPrecCG::Par solverPar; + + solverPar.action = prefix + "_dwf"; + solverPar.maxIteration = 30000; + solverPar.residual = residual; + solverPar.guesser = ""; + app.createModule(solverName, solverPar); +} + // Charm M0 void RbcUkqcd::addM0CharmSolver(Application &app, const std::string solverName, const std::string gaugeName, const double mass, From 47d5d0745c28df0a17b13c6a9260ca07242cad42 Mon Sep 17 00:00:00 2001 From: Ryan Hill Date: Tue, 29 Sep 2026 08:49:29 +0100 Subject: [PATCH 5/6] Add ensemble-specific wrappers for Jlqcd solvers --- Jlqcd.hpp | 289 ++++++++++++++++++++++++++++++++++++++++++++++-------- 1 file changed, 250 insertions(+), 39 deletions(-) diff --git a/Jlqcd.hpp b/Jlqcd.hpp index 9fac006..7b53c37 100644 --- a/Jlqcd.hpp +++ b/Jlqcd.hpp @@ -44,30 +44,34 @@ struct Jlqcd unsigned int L, T, Ls; }; - inline static constexpr EnsembleParameters fUd3SaPar{0.003, 0.015, 1., 2., 64, 128, 8}; - inline static constexpr EnsembleParameters mUd3SaPar{0.0042, 0.025, 1., 2., 48, 96, 8}; - inline static constexpr EnsembleParameters mUd3SbPar{0.0042, 0.018, 1., 2., 48, 96, 8}; - inline static constexpr EnsembleParameters cUd2SaPar{0.0035, 0.040, 1., 2., 32, 64, 12}; - inline static constexpr EnsembleParameters cUd2SaLPar{0.0035, 0.040, 1., 2., 48, 96, 12}; - inline static constexpr EnsembleParameters cUd3SaPar{0.007, 0.040, 1., 2., 32, 64, 12}; - inline static constexpr EnsembleParameters cUd3SbPar{0.007, 0.030, 1., 2., 32, 64, 12}; - inline static constexpr EnsembleParameters testPar{0.010, 0.030, 1., 2., 8, 16, 8}; + inline static constexpr EnsembleParameters fUd3SaPar {0.003, 0.015, 1., 2., 64, 128, 8}; + inline static constexpr EnsembleParameters mUd3SaPar {0.0042, 0.025, 1., 2., 48, 96, 8}; + inline static constexpr EnsembleParameters mUd3SbPar {0.0042, 0.018, 1., 2., 48, 96, 8}; + inline static constexpr EnsembleParameters cUd2SaPar {0.0035, 0.040, 1., 2., 32, 64, 12}; + inline static constexpr EnsembleParameters cUd2SaLPar{0.0035, 0.040, 1., 2., 48, 96, 12}; + inline static constexpr EnsembleParameters cUd3SaPar {0.007, 0.040, 1., 2., 32, 64, 12}; + inline static constexpr EnsembleParameters cUd3SbPar {0.007, 0.030, 1., 2., 32, 64, 12}; + inline static constexpr EnsembleParameters testPar {0.010, 0.030, 1., 2., 8, 16, 8}; - static inline const EnsembleParameters &ensemble(const std::string &name); + static inline const EnsembleParameters &ensemble(const std::string& name); - static inline void addMixedPrecisionSolver(Application &app, const EnsembleParameters &par, + static inline void addMixedPrecisionSolver(Application &app, + const EnsembleParameters &par, const std::string &solverName, - const std::string &gaugeName, const double mass, + const std::string &gaugeName, + const double mass, const std::string &twist, const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); static inline void addNoFailSolver(Application &app, const EnsembleParameters &par, - const std::string &solverName, const std::string &gaugeName, - const double mass, const std::string &twist, + const std::string &solverName, + const std::string &gaugeName, const double mass, + const std::string &twist, const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + // Light solvers (no deflation) static inline void addLightSolver(Application &app, const EnsembleParameters &par, const std::string &solverName, const std::string &gaugeName, const double mass, const std::string &twist, @@ -80,52 +84,114 @@ struct Jlqcd const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addfUd3SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addmUd3SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addmUd3SbLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd2SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd2SaLLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd3SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd3SbLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + + // Strange solvers static inline void addStrangeSolver(Application &app, const EnsembleParameters &par, - const std::string &solverName, const std::string &gaugeName, - const double mass, const std::string &twist, + const std::string &solverName, + const std::string &gaugeName, const double mass, + const std::string &twist, const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); static inline void addStrangeSolver(Application &app, const EnsembleParameters &par, - const std::string &solverName, const std::string &gaugeName, + const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addfUd3SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addmUd3SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addmUd3SbStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd2SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd2SaLStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd3SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + static inline void addcUd3SbStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-8); + + // Heavy solvers static inline void addHeavySolver(Application &app, const EnsembleParameters &par, const std::string &solverName, const std::string &gaugeName, const double mass, const std::string &twist, const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + + static inline void addfUd3SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + static inline void addmUd3SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + static inline void addmUd3SbHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + static inline void addcUd2SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + static inline void addcUd2SaLHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + static inline void addcUd3SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + static inline void addcUd3SbHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist = "0. 0. 0. 0.", + const std::string &boundary = "1 1 1 -1", const double residual = 1.0e-20); + }; -inline const Jlqcd::EnsembleParameters &Jlqcd::ensemble(const std::string &name) +inline const Jlqcd::EnsembleParameters &Jlqcd::ensemble(const std::string& name) { - if (name == "fUd3Sa") - return fUd3SaPar; - if (name == "mUd3Sa") - return mUd3SaPar; - if (name == "mUd3Sb") - return mUd3SbPar; - if (name == "cUd2Sa") - return cUd2SaPar; - if (name == "cUd2SaL") - return cUd2SaLPar; - if (name == "cUd3Sa") - return cUd3SaPar; - if (name == "cUd3Sb") - return cUd3SbPar; - if (name == "test") - return testPar; + if (name == "fUd3Sa" ) return fUd3SaPar; + if (name == "mUd3Sa" ) return mUd3SaPar; + if (name == "mUd3Sb" ) return mUd3SbPar; + if (name == "cUd2Sa" ) return cUd2SaPar; + if (name == "cUd2SaL") return cUd2SaLPar; + if (name == "cUd3Sa" ) return cUd3SaPar; + if (name == "cUd3Sb" ) return cUd3SbPar; + if (name == "test" ) return testPar; throw std::invalid_argument("unknown JLQCD ensemble preset '" + std::string(name) + "'"); } inline void Jlqcd::addMixedPrecisionSolver(Application &app, const EnsembleParameters &par, - const std::string &solverName, - const std::string &gaugeName, const double mass, - const std::string &twist, const std::string &boundary, - const double residual) + const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) { const std::string prefix = solverName; @@ -182,6 +248,7 @@ inline void Jlqcd::addNoFailSolver(Application &app, const EnsembleParameters &p app.createModule(solverName, solverPar); } +// Light solvers (no deflation) inline void Jlqcd::addLightSolver(Application &app, const EnsembleParameters &par, const std::string &solverName, const std::string &gaugeName, const double mass, const std::string &twist, @@ -198,12 +265,62 @@ inline void Jlqcd::addLightSolver(Application &app, const EnsembleParameters &pa addLightSolver(app, par, solverName, gaugeName, par.ml, twist, boundary, residual); } +inline void Jlqcd::addfUd3SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addLightSolver(app, fUd3SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addmUd3SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addLightSolver(app, mUd3SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addmUd3SbLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addLightSolver(app, mUd3SbPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd2SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addLightSolver(app, cUd2SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd2SaLLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addLightSolver(app, cUd2SaLPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd3SaLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addLightSolver(app, cUd3SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd3SbLightSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addLightSolver(app, cUd3SbPar, solverName, gaugeName, twist, boundary, residual); +} + +// Strange solvers inline void Jlqcd::addStrangeSolver(Application &app, const EnsembleParameters &par, const std::string &solverName, const std::string &gaugeName, const double mass, const std::string &twist, const std::string &boundary, const double residual) { - addMixedPrecisionSolver(app, par, solverName, gaugeName, mass, twist, boundary, residual); + addMixedPrecisionSolver(app, par, solverName, gaugeName, mass, twist, boundary, residual); } inline void Jlqcd::addStrangeSolver(Application &app, const EnsembleParameters &par, @@ -211,9 +328,59 @@ inline void Jlqcd::addStrangeSolver(Application &app, const EnsembleParameters & const std::string &twist, const std::string &boundary, const double residual) { - addStrangeSolver(app, par, solverName, gaugeName, par.ms, twist, boundary, residual); + addStrangeSolver(app, par, solverName, gaugeName, par.ms, twist, boundary, residual); } +inline void Jlqcd::addfUd3SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addStrangeSolver(app, fUd3SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addmUd3SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addStrangeSolver(app, mUd3SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addmUd3SbStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addStrangeSolver(app, mUd3SbPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd2SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addStrangeSolver(app, cUd2SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd2SaLStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addStrangeSolver(app, cUd2SaLPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd3SaStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addStrangeSolver(app, cUd3SaPar, solverName, gaugeName, twist, boundary, residual); +} + +inline void Jlqcd::addcUd3SbStrangeSolver(Application &app, const std::string &solverName, + const std::string &gaugeName, const std::string &twist, + const std::string &boundary, const double residual) +{ + addStrangeSolver(app, cUd3SbPar, solverName, gaugeName, twist, boundary, residual); +} + +// Heavy solvers inline void Jlqcd::addHeavySolver(Application &app, const EnsembleParameters &par, const std::string &solverName, const std::string &gaugeName, const double mass, const std::string &twist, @@ -221,4 +388,48 @@ inline void Jlqcd::addHeavySolver(Application &app, const EnsembleParameters &pa { addNoFailSolver(app, par, solverName, gaugeName, mass, twist, boundary, residual); } + +inline void Jlqcd::addfUd3SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addHeavySolver(app, fUd3SaPar, solverName, gaugeName, mass, twist, boundary, residual); +} +inline void Jlqcd::addmUd3SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addHeavySolver(app, mUd3SaPar, solverName, gaugeName, mass, twist, boundary, residual); +} +inline void Jlqcd::addmUd3SbHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addHeavySolver(app, mUd3SbPar, solverName, gaugeName, mass, twist, boundary, residual); +} +inline void Jlqcd::addcUd2SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addHeavySolver(app, cUd2SaPar, solverName, gaugeName, mass, twist, boundary, residual); +} +inline void Jlqcd::addcUd2SaLHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addHeavySolver(app, cUd2SaLPar, solverName, gaugeName, mass, twist, boundary, residual); +} +inline void Jlqcd::addcUd3SaHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addHeavySolver(app, cUd3SaPar, solverName, gaugeName, mass, twist, boundary, residual); +} +inline void Jlqcd::addcUd3SbHeavySolver(Application &app, const std::string &solverName, const std::string &gaugeName, + const double mass, const std::string &twist, + const std::string &boundary, const double residual) +{ + addHeavySolver(app, cUd3SbPar, solverName, gaugeName, mass, twist, boundary, residual); +} + } // namespace hadpresets From 87aba30cd755a4c6574f7f56d00c10f67a646eb6 Mon Sep 17 00:00:00 2001 From: Ryan Hill Date: Tue, 29 Sep 2026 13:48:34 +0100 Subject: [PATCH 6/6] Update C0L charm Ls --- RbcUkqcd.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/RbcUkqcd.hpp b/RbcUkqcd.hpp index 85f6f42..5b8d5e1 100644 --- a/RbcUkqcd.hpp +++ b/RbcUkqcd.hpp @@ -893,7 +893,7 @@ void RbcUkqcd::addC0LCharmSolver(Application &app, const std::string solverName, MAction::ScaledDWF::Par actionPar; actionPar.gauge = prefix + "_gauge_3stout"; - actionPar.Ls = 24; + actionPar.Ls = 12; actionPar.M5 = 1.; actionPar.mass = mass; actionPar.scale = 2.;