diff --git a/extras/Hadrons/Environment.cc b/extras/Hadrons/Environment.cc index 82b0dda1..6554122e 100644 --- a/extras/Hadrons/Environment.cc +++ b/extras/Hadrons/Environment.cc @@ -270,7 +270,7 @@ int Environment::getObjectModule(const std::string name) const unsigned int Environment::getObjectLs(const unsigned int address) const { - if (hasObject(address)) + if (hasCreatedObject(address)) { return object_[address].Ls; } diff --git a/extras/Hadrons/Modules.hpp b/extras/Hadrons/Modules.hpp index 135934fb..d1f947a2 100644 --- a/extras/Hadrons/Modules.hpp +++ b/extras/Hadrons/Modules.hpp @@ -1,31 +1,3 @@ -/************************************************************************************* - -Grid physics library, www.github.com/paboyle/Grid - -Source file: extras/Hadrons/Modules.hpp - -Copyright (C) 2015-2018 - -Author: Antonin Portelli -Author: Lanny91 - -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 #include #include @@ -43,11 +15,12 @@ See the full license in the file "LICENSE" in the top level distribution directo #include #include #include +#include #include #include #include -#include #include +#include #include #include #include diff --git a/extras/Hadrons/Modules/MSolver/LocalCoherenceLanczos.hpp b/extras/Hadrons/Modules/MSolver/LocalCoherenceLanczos.hpp new file mode 100644 index 00000000..1ce87c2b --- /dev/null +++ b/extras/Hadrons/Modules/MSolver/LocalCoherenceLanczos.hpp @@ -0,0 +1,211 @@ +#ifndef Hadrons_MSolver_LocalCoherenceLanczos_hpp_ +#define Hadrons_MSolver_LocalCoherenceLanczos_hpp_ + +#include +#include +#include +#include + +#ifndef DEFAULT_LANCZOS_NBASIS +#define DEFAULT_LANCZOS_NBASIS 60 +#endif + +BEGIN_HADRONS_NAMESPACE + +/****************************************************************************** + * LocalCoherenceLanczos * + ******************************************************************************/ +BEGIN_MODULE_NAMESPACE(MSolver) + +class LocalCoherenceLanczosPar: Serializable +{ +public: + GRID_SERIALIZABLE_CLASS_MEMBERS(LocalCoherenceLanczosPar, + std::string, action, + int, doFine, + int, doCoarse, + LanczosParams, fineParams, + LanczosParams, coarseParams, + ChebyParams, smoother, + RealD, coarseRelaxTol, + std::string, blockSize); +}; + +template +class TLocalCoherenceLanczos: public Module +{ +public: + FERM_TYPE_ALIASES(FImpl,); + typedef LocalCoherenceLanczos LCL; + typedef typename LCL::FineField FineField; + typedef typename LCL::CoarseField CoarseField; + typedef SchurDiagMooeeOperator SchurFMat; +public: + // constructor + TLocalCoherenceLanczos(const std::string name); + // destructor + virtual ~TLocalCoherenceLanczos(void) = default; + // dependency relation + virtual std::vector getInput(void); + virtual std::vector getOutput(void); + // setup + virtual void setup(void); + // execution + virtual void execute(void); +private: + void makeCoarseGrid(void); +private: + std::vector coarseDim_; + int Ls_, cLs_{1}; + std::unique_ptr coarseGrid4_{nullptr}; + std::unique_ptr coarseGrid_{nullptr}; + std::unique_ptr coarseGrid4Rb_{nullptr}; + std::unique_ptr coarseGridRb_{nullptr}; + std::string fevecName_, cevecName_, fevalName_, cevalName_; +}; + +MODULE_REGISTER_NS(LocalCoherenceLanczos, + ARG(TLocalCoherenceLanczos), + MSolver); + +/****************************************************************************** + * TLocalCoherenceLanczos implementation * + ******************************************************************************/ +// constructor ///////////////////////////////////////////////////////////////// +template +TLocalCoherenceLanczos::TLocalCoherenceLanczos(const std::string name) +: Module(name) +{ + fevecName_ = getName() + "_fineEvec"; + cevecName_ = getName() + "_coarseEvec"; + fevalName_ = getName() + "_fineEval"; + cevalName_ = getName() + "_coarseEval"; +} + +// dependencies/products /////////////////////////////////////////////////////// +template +std::vector TLocalCoherenceLanczos::getInput(void) +{ + std::vector in = {par().action}; + + return in; +} + +template +std::vector TLocalCoherenceLanczos::getOutput(void) +{ + std::vector out = {fevecName_, cevecName_, fevalName_, cevalName_}; + + return out; +} + +// setup /////////////////////////////////////////////////////////////////////// +template +void TLocalCoherenceLanczos::makeCoarseGrid(void) +{ + int nd = env().getNd(); + std::vector blockSize = strToVec(par().blockSize); + auto fineDim = env().getDim(); + + Ls_ = env().getObjectLs(par().action); + env().createGrid(Ls_); + coarseDim_.resize(nd); + for (int d = 0; d < coarseDim_.size(); d++) + { + coarseDim_[d] = fineDim[d]/blockSize[d]; + if (coarseDim_[d]*blockSize[d] != fineDim[d]) + { + HADRON_ERROR(Size, "Fine dimension " + std::to_string(d) + + " (" + std::to_string(fineDim[d]) + + ") not divisible by coarse dimension (" + + std::to_string(coarseDim_[d]) + ")"); + } + } + if (blockSize.size() > nd) + { + cLs_ = Ls_/blockSize[nd]; + if (cLs_*blockSize[nd] != Ls_) + { + HADRON_ERROR(Size, "Fine Ls (" + std::to_string(Ls_) + + ") not divisible by coarse Ls (" + + std::to_string(cLs_) + ")"); + } + } + if (Ls_ > 1) + { + coarseGrid4_.reset(SpaceTimeGrid::makeFourDimGrid( + coarseDim_, GridDefaultSimd(nd, vComplex::Nsimd()), + GridDefaultMpi())); + coarseGrid4Rb_.reset(SpaceTimeGrid::makeFourDimRedBlackGrid(coarseGrid4_.get())); + coarseGrid_.reset(SpaceTimeGrid::makeFiveDimGrid(cLs_, coarseGrid4_.get())); + coarseGridRb_.reset(SpaceTimeGrid::makeFiveDimRedBlackGrid(cLs_, coarseGrid4_.get())); + } + else + { + coarseGrid_.reset(SpaceTimeGrid::makeFourDimGrid( + coarseDim_, GridDefaultSimd(nd, vComplex::Nsimd()), + GridDefaultMpi())); + coarseGridRb_.reset(SpaceTimeGrid::makeFourDimRedBlackGrid(coarseGrid_.get())); + } +} + +template +void TLocalCoherenceLanczos::setup(void) +{ + LOG(Message) << "Setting up local coherence Lanczos eigensolver for" + << " action '" << par().action << "' (" << nBasis + << " eigenvectors)..." << std::endl; + + if (!coarseGrid_) + { + makeCoarseGrid(); + } + envCreate(std::vector, fevecName_, Ls_, par().fineParams.Nm, + env().getRbGrid(Ls_)); + envCreate(std::vector, cevecName_, Ls_, par().coarseParams.Nm, + coarseGridRb_.get()); + envCreate(std::vector, fevalName_, Ls_, par().fineParams.Nm); + envCreate(std::vector, cevalName_, Ls_, par().coarseParams.Nm); + envTmp(SchurFMat, "mat", Ls_, envGet(FMat, par().action)); + envGetTmp(SchurFMat, mat); + envTmp(LCL, "solver", Ls_, env().getRbGrid(Ls_), coarseGridRb_.get(), mat, Odd, + envGet(std::vector, fevecName_), + envGet(std::vector, cevecName_), + envGet(std::vector, fevalName_), + envGet(std::vector, cevalName_)); +} + +// execution /////////////////////////////////////////////////////////////////// +template +void TLocalCoherenceLanczos::execute(void) +{ + auto &fine = par().fineParams; + auto &coarse = par().coarseParams; + + envGetTmp(LCL, solver); + if (par().doFine) + { + LOG(Message) << "Performing fine grid IRL -- Nstop= " + << fine.Nstop << ", Nk= " << fine.Nk << ", Nm= " + << fine.Nm << std::endl; + solver.calcFine(fine.Cheby, fine.Nstop, fine.Nk, fine.Nm, fine.resid, + fine.MaxIt, fine.betastp, fine.MinRes); + LOG(Message) << "Orthogonalising" << std::endl; + solver.Orthogonalise(); + } + if (par().doCoarse) + { + LOG(Message) << "Performing coarse grid IRL -- Nstop= " + << fine.Nstop << ", Nk= " << fine.Nk << ", Nm= " + << fine.Nm << std::endl; + solver.calcCoarse(coarse.Cheby, par().smoother, par().coarseRelaxTol, + coarse.Nstop, coarse.Nk, coarse.Nm, coarse.resid, + coarse.MaxIt, coarse.betastp, coarse.MinRes); + } +} + +END_MODULE_NAMESPACE + +END_HADRONS_NAMESPACE + +#endif // Hadrons_MSolver_LocalCoherenceLanczos_hpp_