From: Christian Blau Date: Fri, 20 Sep 2019 09:50:33 +0000 (+0200) Subject: Grompp error for mismatching nst for energy calulation and densityfitting X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=c87229d67b63ca58335bf87d7da36e72286383f9;p=alexxy%2Fgromacs.git Grompp error for mismatching nst for energy calulation and densityfitting When using multiple-time-stepping with density-guided simulations, the energy evaluation frequency of the density-guided simulations must match the "nstcalcenergy" option. Change-Id: Ifcde8cada8a91540d69c9f91079bfc91c4415ba5 --- diff --git a/docs/reference-manual/special/density-guided-simulation.rst b/docs/reference-manual/special/density-guided-simulation.rst index 7362c70597..301da01260 100644 --- a/docs/reference-manual/special/density-guided-simulation.rst +++ b/docs/reference-manual/special/density-guided-simulation.rst @@ -137,8 +137,8 @@ The cost of applying forces every integration step is reduced when applying the density-guided simulation forces only every :math:`N` steps. The applied force is scaled by :math:`N` to approximate the same effective Hamiltonian as when applying the forces every step, while maintaining time-reversibility and energy -conservation. Note that for this setting, the energy output frequency should -be a multiple of :math:`N`. +conservation. Note that for this setting, the energy output frequency must be a +multiple of :math:`N`. The maximal time-step should not exceed the fastest oscillation period of any atom within the map potential divided by :math:`\pi`. This oscillation period diff --git a/src/gromacs/applied_forces/densityfitting.cpp b/src/gromacs/applied_forces/densityfitting.cpp index fb5477fa82..d067bcb727 100644 --- a/src/gromacs/applied_forces/densityfitting.cpp +++ b/src/gromacs/applied_forces/densityfitting.cpp @@ -263,6 +263,12 @@ class DensityFitting final : public IMDModule }; notifier->notifier_.subscribe(readInternalParametersFunction); + // Checking for consistency with all .mdp options + const auto checkEnergyCaluclationFrequencyFunction = [this](EnergyCalculationFrequencyErrors * energyCalculationFrequencyErrors) { + densityFittingOptions_.checkEnergyCaluclationFrequency(energyCalculationFrequencyErrors); + }; + notifier->notifier_.subscribe(checkEnergyCaluclationFrequencyFunction); + // constructing local atom sets during simulation setup const auto setLocalAtomSetFunction = [this](LocalAtomSetManager *localAtomSetManager) { this->constructLocalAtomSet(localAtomSetManager); diff --git a/src/gromacs/applied_forces/densityfittingoptions.cpp b/src/gromacs/applied_forces/densityfittingoptions.cpp index 2543cd543c..bf9f89aa12 100644 --- a/src/gromacs/applied_forces/densityfittingoptions.cpp +++ b/src/gromacs/applied_forces/densityfittingoptions.cpp @@ -51,6 +51,7 @@ #include "gromacs/utility/exceptions.h" #include "gromacs/utility/keyvaluetreebuilder.h" #include "gromacs/utility/keyvaluetreetransform.h" +#include "gromacs/utility/mdmodulenotification.h" #include "gromacs/utility/strconvert.h" #include "densityfittingamplitudelookup.h" @@ -248,6 +249,17 @@ void DensityFittingOptions::readInternalParametersFromKvt(const KeyValueTreeObje [](const KeyValueTreeValue &val) { return val.cast(); }); } +void DensityFittingOptions::checkEnergyCaluclationFrequency(EnergyCalculationFrequencyErrors * energyCalculationFrequencyErrors) const +{ + if (energyCalculationFrequencyErrors->energyCalculationIntervalInSteps() % parameters_.calculationIntervalInSteps_ != 0) + { + energyCalculationFrequencyErrors->addError("nstcalcenergy (" + + toString(energyCalculationFrequencyErrors->energyCalculationIntervalInSteps()) + + ") is not a multiple of " + DensityFittingModuleInfo::name_ + "-" + c_everyNStepsTag_ + + " (" + toString(parameters_.calculationIntervalInSteps_) + ") ."); + } +} + const std::string &DensityFittingOptions::referenceDensityFileName() const { return referenceDensityFileName_; diff --git a/src/gromacs/applied_forces/densityfittingoptions.h b/src/gromacs/applied_forces/densityfittingoptions.h index d629a7905f..034f3d3a33 100644 --- a/src/gromacs/applied_forces/densityfittingoptions.h +++ b/src/gromacs/applied_forces/densityfittingoptions.h @@ -51,6 +51,7 @@ namespace gmx { +class EnergyCalculationFrequencyErrors; class IndexGroupsAndNames; class KeyValueTreeObject; class KeyValueTreeBuilder; @@ -101,6 +102,9 @@ class DensityFittingOptions final : public IMdpOptionProvider //! Return the file name of the reference density const std::string &referenceDensityFileName() const; + //! Check if input parameters are consistent with other simulation parameters + void checkEnergyCaluclationFrequency(EnergyCalculationFrequencyErrors * energyCalculationFrequencyErrors) const; + private: const std::string c_activeTag_ = "active"; diff --git a/src/gromacs/gmxpreprocess/grompp.cpp b/src/gromacs/gmxpreprocess/grompp.cpp index a3a2c3f05f..ba57fa1f0d 100644 --- a/src/gromacs/gmxpreprocess/grompp.cpp +++ b/src/gromacs/gmxpreprocess/grompp.cpp @@ -1882,7 +1882,7 @@ int gmx_grompp(int argc, char *argv[]) { fprintf(stderr, "checking input for internal consistency...\n"); } - check_ir(mdparin, ir, opts, wi); + check_ir(mdparin, mdModules.notifier(), ir, opts, wi); if (ir->ld_seed == -1) { diff --git a/src/gromacs/gmxpreprocess/readir.cpp b/src/gromacs/gmxpreprocess/readir.cpp index 8ac8242b70..339c718580 100644 --- a/src/gromacs/gmxpreprocess/readir.cpp +++ b/src/gromacs/gmxpreprocess/readir.cpp @@ -238,8 +238,8 @@ static void process_interaction_modifier(int *eintmod) } } -void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, - warninp_t wi) +void check_ir(const char *mdparin, const gmx::MdModulesNotifier &mdModulesNotifier, + t_inputrec *ir, t_gromppopts *opts, warninp_t wi) /* Check internal consistency. * NOTE: index groups are not set here yet, don't check things * like temperature coupling group options here, but in triple_check @@ -527,6 +527,17 @@ void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, check_nst("nstcalcenergy", ir->nstcalcenergy, "nstenergy", &ir->nstenergy, wi); } + + // Inquire all MdModules, if their parameters match with the energy + // calculation frequency + gmx::EnergyCalculationFrequencyErrors energyCalculationFrequencyErrors(ir->nstcalcenergy); + mdModulesNotifier.notifier_.notify(&energyCalculationFrequencyErrors); + + // Emit all errors from the energy calculation frequency checks + for (const std::string &energyFrequencyErrorMessage : energyCalculationFrequencyErrors.errorMessages()) + { + warning_error(wi, energyFrequencyErrorMessage); + } } if (ir->nsteps == 0 && !ir->bContinuation) diff --git a/src/gromacs/gmxpreprocess/readir.h b/src/gromacs/gmxpreprocess/readir.h index 3be50585b5..f5375791fb 100644 --- a/src/gromacs/gmxpreprocess/readir.h +++ b/src/gromacs/gmxpreprocess/readir.h @@ -93,8 +93,8 @@ void init_inputrec_strings(); /*! \brief Clean up object that holds strings parsed from an .mdp file */ void done_inputrec_strings(); -void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, - warninp_t wi); +void check_ir(const char *mdparin, const gmx::MdModulesNotifier &mdModulesNotifier, + t_inputrec *ir, t_gromppopts *opts, warninp_t wi); /* Validate inputrec data. * Fatal errors will be added to nerror. */ diff --git a/src/gromacs/gmxpreprocess/tests/readir.cpp b/src/gromacs/gmxpreprocess/tests/readir.cpp index a1c83c0522..a8a391ba09 100644 --- a/src/gromacs/gmxpreprocess/tests/readir.cpp +++ b/src/gromacs/gmxpreprocess/tests/readir.cpp @@ -99,7 +99,7 @@ class GetIrTest : public ::testing::Test get_ir(inputMdpFilename.c_str(), outputMdpFilename.c_str(), &mdModules_, &ir_, &opts_, WriteMdpHeader::no, wi_); - check_ir(inputMdpFilename.c_str(), &ir_, &opts_, wi_); + check_ir(inputMdpFilename.c_str(), mdModules_.notifier(), &ir_, &opts_, wi_); // Now check bool failure = warning_errors_exist(wi_); TestReferenceData data; diff --git a/src/gromacs/utility/mdmodulenotification.h b/src/gromacs/utility/mdmodulenotification.h index 7f669aa32c..654007ecae 100644 --- a/src/gromacs/utility/mdmodulenotification.h +++ b/src/gromacs/utility/mdmodulenotification.h @@ -45,6 +45,7 @@ #define GMX_MDRUNUTILITY_MDMODULENOTIFICATION_H #include +#include #include struct t_commrec; @@ -209,11 +210,50 @@ struct MdModulesEnergyOutputToDensityFittingRequestChecker bool energyOutputToDensityFitting_ = false; }; +/*! \libinternal + * \brief Collect errors for the energy calculation frequency. + * + * Collect errors regarding energy calculation frequencies as strings that then + * may be used to issue errors. + * + * \note The mdp option "nstcalcenergy" is altered after reading the .mdp input + * and only used in certain integrators, thus this class is to be used + * only after all these operations are done. + */ +class EnergyCalculationFrequencyErrors +{ + public: + //! Construct by setting the energy calculation frequency + EnergyCalculationFrequencyErrors(int64_t energyCalculationIntervalInSteps) : + energyCalculationIntervalInSteps_(energyCalculationIntervalInSteps){} + //! Return the number of steps of an energy calculation interval + std::int64_t energyCalculationIntervalInSteps() const + { + return energyCalculationIntervalInSteps_; + } + //! Collect error messages + void addError(const std::string &errorMessage) + { + errorMessages_.push_back(errorMessage); + } + //! Return error messages + const std::vector &errorMessages() const + { + return errorMessages_; + } + private: + //! The frequency of energy calculations + const std::int64_t energyCalculationIntervalInSteps_; + //! The error messages + std::vector errorMessages_; +}; + struct MdModulesNotifier { //! Register callback function types for MdModule registerMdModuleNotification< const t_commrec &, + EnergyCalculationFrequencyErrors *, IndexGroupsAndNames, KeyValueTreeObjectBuilder, const KeyValueTreeObject &, diff --git a/src/programs/mdrun/tests/densityfittingmodule.cpp b/src/programs/mdrun/tests/densityfittingmodule.cpp index 39797876ea..adf4e4d241 100644 --- a/src/programs/mdrun/tests/densityfittingmodule.cpp +++ b/src/programs/mdrun/tests/densityfittingmodule.cpp @@ -44,6 +44,9 @@ #include +#include +#include + #include "gromacs/topology/ifunc.h" #include "gromacs/utility/stringutil.h" @@ -101,6 +104,16 @@ class DensityFittingTest : public MdrunTestFixture "density-guided-simulation-reference-density-filename = %s\n", TestFileManager::getInputFilePath("ellipsoid-density.mrc").c_str() ); + //! Mdp values for md integrator with default density fitting parameters. + const std::string mdpMdDensfitYesUnsetValues = formatString( + "integrator = md\n" + "nsteps = 2\n" + "cutoff-scheme = verlet\n" + "density-guided-simulation-active = yes\n" + "density-guided-simulation-group = FirstThreeOfTwelve\n" + "density-guided-simulation-reference-density-filename = %s\n", + TestFileManager::getInputFilePath("ellipsoid-density.mrc").c_str() ); + //! Mdp values for steepest-decent energy minimization with density fitting values set to non-defaults. const std::string mdpDensiftAllDefaultsChanged_ = formatString( "density-guided-simulation-similarity-measure = relative-entropy\n" @@ -110,7 +123,11 @@ class DensityFittingTest : public MdrunTestFixture "density-guided-simulation-gaussian-transform-spreading-range-in-multiples-of-width = 6\n" "density-guided-simulation-normalize-densities = false\n" ); - + //! Set mdp values so that energy calculation interval and density guided simulation interval mismatch. + const std::string mdpEnergyAndDensityfittingIntervalMismatch_ = formatString( + "nstcalcenergy = 7\n" + "density-guided-simulation-nst = 3\n" + ); //! The command line to call mdrun CommandLine commandLineForMdrun_; }; @@ -144,6 +161,15 @@ TEST_F(DensityFittingTest, EnergyMinimizationEnergyCorrectForRelativeEntropy) checkMdrun(expectedEnergyTermMagnitude); } +/* Test that grompp exits with error message if energy evaluation frequencies + * do not match. + */ +TEST_F(DensityFittingTest, GromppErrorWhenEnergyEvaluationFrequencyMismatch) +{ + runner_.useStringAsMdpFile(mdpMdDensfitYesUnsetValues + mdpEnergyAndDensityfittingIntervalMismatch_); + + EXPECT_DEATH_IF_SUPPORTED(runner_.callGrompp(), ".*is not a multiple of density-guided-simulation-nst.*"); +} } // namespace test } // namespace gmx