From c7411433871452b75f08052a3f1ac3fed2dce13c Mon Sep 17 00:00:00 2001 From: Mark Abraham Date: Tue, 29 Dec 2015 18:39:44 +0100 Subject: [PATCH] Updated mdrun comparison infrastructure The functionality for running grompp to produce a .tpr file useful for comparing two mdrun calculations works better as free functions than a base class. This patch is essentially a rewrite of the interface, though the core functionality is roughly the same. Used raw string literals for nicer looking embedded strings. Refs #1868 Change-Id: I4b5137b34c3da286e344809b71afd985d9c7bfab --- src/programs/mdrun/tests/CMakeLists.txt | 2 +- ...uncomparisonfixture.cpp => mdruncomparison.cpp} | 168 +++++++++++---------- src/programs/mdrun/tests/mdruncomparison.h | 110 ++++++++++++++ src/programs/mdrun/tests/mdruncomparisonfixture.h | 143 ------------------ 4 files changed, 196 insertions(+), 227 deletions(-) rename src/programs/mdrun/tests/{mdruncomparisonfixture.cpp => mdruncomparison.cpp} (54%) create mode 100644 src/programs/mdrun/tests/mdruncomparison.h delete mode 100644 src/programs/mdrun/tests/mdruncomparisonfixture.h diff --git a/src/programs/mdrun/tests/CMakeLists.txt b/src/programs/mdrun/tests/CMakeLists.txt index 7bdc45edef..778fd80411 100644 --- a/src/programs/mdrun/tests/CMakeLists.txt +++ b/src/programs/mdrun/tests/CMakeLists.txt @@ -35,7 +35,7 @@ # make an "object library" for code that we re-use for both kinds of tests add_library(mdrun_test_objlib OBJECT energyreader.cpp - mdruncomparisonfixture.cpp + mdruncomparison.cpp moduletest.cpp terminationhelper.cpp # PME tests diff --git a/src/programs/mdrun/tests/mdruncomparisonfixture.cpp b/src/programs/mdrun/tests/mdruncomparison.cpp similarity index 54% rename from src/programs/mdrun/tests/mdruncomparisonfixture.cpp rename to src/programs/mdrun/tests/mdruncomparison.cpp index 14f63886a5..165417bd26 100644 --- a/src/programs/mdrun/tests/mdruncomparisonfixture.cpp +++ b/src/programs/mdrun/tests/mdruncomparison.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2016, by the GROMACS development team, led by + * Copyright (c) 2016,2018, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -34,14 +34,14 @@ */ /*! \internal \file * \brief - * Implements classes in mdruncomparisonfixture.h. + * Implements declarations from in mdruncomparison.h. * * \author Mark Abraham * \ingroup module_mdrun_integration_tests */ #include "gmxpre.h" -#include "mdruncomparisonfixture.h" +#include "mdruncomparison.h" #include #include @@ -57,17 +57,13 @@ namespace gmx namespace test { -MdrunComparisonFixture::~MdrunComparisonFixture() -{ -} - namespace { //! Helper typedef -typedef std::map MdpFileValues; +using MdpFileValues = std::map; -//! Database of .mdp strings that supports MdrunComparisonFixture::prepareSimulation() +//! Database of .mdp strings that supports prepareDefaultMdpValues() MdpFileValues mdpFileValueDatabase_g { // Simple system with 12 argon atoms, fairly widely separated @@ -135,35 +131,36 @@ MdpFileValues mdpFileValueDatabase_g }, { "other", - "free-energy = yes\n" - "sc-alpha = 0.5\n" - "sc-r-power = 6\n" - "nstdhdl = 4\n" - "init-lambda-state = 3\n" - "fep_lambdas = 0.00 0.50 1.00 1.00 1.00\n" - "vdw_lambdas = 0.00 0.00 0.00 0.50 1.00\n" - "couple-moltype = nonanol\n" - "couple-lambda0 = vdw-q\n" - "couple-lambda1 = none\n" - "couple-intramol = yes\n" + R"(free-energy = yes + sc-alpha = 0.5 + sc-r-power = 6 + nstdhdl = 4 + init-lambda-state = 3 + fep_lambdas = 0.00 0.50 1.00 1.00 1.00 + vdw_lambdas = 0.00 0.00 0.00 0.50 1.00 + couple-moltype = nonanol + couple-lambda0 = vdw-q + couple-lambda1 = none + couple-intramol = yes)" } } } }; -} // namespace - -MdrunComparisonFixture::MdpFieldValues MdrunComparisonFixture::prepareMdpFieldValues(const char *simulationName) +/*! \brief Prepare default .mdp values + * + * Insert suitable .mdp defaults, so that \c mdpFileValueDatabase_g + * does not need to specify repetitive values. This works because + * std::map.insert() does not over-write elements that already exist. + * + * \todo ideally some of these default values should be the same as + * grompp uses, and sourced from the same place, but that code is a + * bit of a jungle until we transition to using IMdpOptions more. + * + * \throws std::bad_alloc if out of memory + * std::out_of_range if \c simulationName is not in the database */ +MdpFieldValues prepareDefaultMdpFieldValues(const char *simulationName) { - /* Insert suitable .mdp defaults, so that the database set up - * above does not need to specify redundant values. This works - * because std::map.insert() does not over-write elements that - * already exist. - * - * TODO ideally some of these default values should be the same as - * grompp uses, and sourced from the same place, but that code is - * a bit of a jungle. */ - - typedef MdpFieldValues::value_type MdpField; + using MdpField = MdpFieldValues::value_type; auto &mdpFieldValues = mdpFileValueDatabase_g.at(simulationName); mdpFieldValues.insert(MdpField("nsteps", "16")); @@ -176,10 +173,25 @@ MdrunComparisonFixture::MdpFieldValues MdrunComparisonFixture::prepareMdpFieldVa return mdpFieldValues; } -void MdrunComparisonFixture::prepareMdpFile(const MdpFieldValues &mdpFieldValues, - const char *integrator, - const char *tcoupl, - const char *pcoupl) +} // namespace + +MdpFieldValues +prepareMdpFieldValues(const char *simulationName, + const char *integrator, + const char *tcoupl, + const char *pcoupl) +{ + using MdpField = MdpFieldValues::value_type; + + auto mdpFieldValues = prepareDefaultMdpFieldValues(simulationName); + mdpFieldValues.insert(MdpField("integrator", integrator)); + mdpFieldValues.insert(MdpField("tcoupl", tcoupl)); + mdpFieldValues.insert(MdpField("pcoupl", pcoupl)); + return mdpFieldValues; +} + +std::string +prepareMdpFileContents(const MdpFieldValues &mdpFieldValues) { /* Set up an .mdp file that permits a highly reproducible * simulation. The format string needs to be configured with @@ -200,52 +212,42 @@ void MdrunComparisonFixture::prepareMdpFile(const MdpFieldValues &mdpFieldValues * energies were not computed with those from rerun on the same * coordinates. */ - runner_.useStringAsMdpFile(formatString("rcoulomb = 0.7\n" - "rvdw = 0.7\n" - "rlist = -1\n" - "bd-fric = 1000\n" - "verlet-buffer-tolerance = 0.000001\n" - "nsteps = %s\n" - "nstenergy = 4\n" - "nstlist = 8\n" - "nstxout = 4\n" - "nstvout = 4\n" - "nstfout = 4\n" - "integrator = %s\n" - "ld-seed = 234262\n" - "tcoupl = %s\n" - "ref-t = %s\n" - "tau-t = 1\n" - "tc-grps = System\n" - "pcoupl = %s\n" - "pcoupltype = isotropic\n" - "ref-p = 1\n" - "tau-p = %s\n" - "compressibility = %s\n" - "constraints = %s\n" - "constraint-algorithm = lincs\n" - "lincs-order = 2\n" - "lincs-iter = 5\n" - "%s", - mdpFieldValues.at("nsteps").c_str(), - integrator, tcoupl, - mdpFieldValues.at("ref-t").c_str(), - pcoupl, - mdpFieldValues.at("tau-p").c_str(), - mdpFieldValues.at("compressibility").c_str(), - mdpFieldValues.at("constraints").c_str(), - mdpFieldValues.at("other").c_str())); -} - -void MdrunComparisonFixture::runTest(const char *simulationName, - const char *integrator, - const char *tcoupl, - const char *pcoupl, - FloatingPointTolerance tolerance) -{ - CommandLine caller; - caller.append("grompp"); - runTest(caller, simulationName, integrator, tcoupl, pcoupl, tolerance); + return formatString(R"(rcoulomb = 0.7 + rvdw = 0.7 + rlist = -1 + bd-fric = 1000 + verlet-buffer-tolerance = 0.000001 + nsteps = %s + nstenergy = 4 + nstlist = 8 + nstxout = 4 + nstvout = 4 + nstfout = 4 + integrator = %s + ld-seed = 234262 + tcoupl = %s + ref-t = %s + tau-t = 1 + tc-grps = System + pcoupl = %s + pcoupltype = isotropic + ref-p = 1 + tau-p = %s + compressibility = %s + constraints = %s + constraint-algorithm = lincs + lincs-order = 2 + lincs-iter = 5 + %s)", + mdpFieldValues.at("nsteps").c_str(), + mdpFieldValues.at("integrator").c_str(), + mdpFieldValues.at("tcoupl").c_str(), + mdpFieldValues.at("ref-t").c_str(), + mdpFieldValues.at("pcoupl").c_str(), + mdpFieldValues.at("tau-p").c_str(), + mdpFieldValues.at("compressibility").c_str(), + mdpFieldValues.at("constraints").c_str(), + mdpFieldValues.at("other").c_str()); } } // namespace test diff --git a/src/programs/mdrun/tests/mdruncomparison.h b/src/programs/mdrun/tests/mdruncomparison.h new file mode 100644 index 0000000000..f7a36f0f6b --- /dev/null +++ b/src/programs/mdrun/tests/mdruncomparison.h @@ -0,0 +1,110 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2016,2018, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS 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 + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \libinternal + * + * \brief Functionality for testing whether calls to mdrun produce the + * same energy and force quantities when they should do so. + */ +#ifndef GMX_MDRUN_TESTS_MDRUNCOMPARISON_H +#define GMX_MDRUN_TESTS_MDRUNCOMPARISON_H + +#include +#include +#include +#include +#include + +#include + +#include "testutils/cmdlinetest.h" + +#include "moduletest.h" + +namespace gmx +{ + +namespace test +{ + +//! Helper typedef +using MdpFieldValues = std::map; + +/*! \brief Set up values for an .mdp file that permits a highly + * reproducible simulation. + * + * An internal database of several kinds of simulation useful for such + * comparisons is available, whose \c simulationName keys are + * - argon12 + * - argon5832 + * - spc5 + * - spc216 + * - alanine_vsite_vacuo + * - alanine_vsite_solvated + * - nonanol + * + * Some of these systems are pretty minimal, because having + * few atoms means few interactions, highly reproducible + * forces, and allows tests to focus on the correctness of the + * implementation of high-level mdrun features. The boxes are + * of a reasonable size so that domain decomposition is + * possible. The pressure-coupling parameters are isotropic, + * and set up so that there will not be dramatic collapse of + * volume over the handful of MD steps that will be run. A + * single temperature-coupling group is used. + * + * \param[in] simulationName The name of the simulation, which indexes the database + * \param[in] integrator The integrator to use + * \param[in] tcoupl The temperature-coupling algorithm to use + * \param[in] pcoupl The pressure-coupling algorithm to use + * \return Mdp file values + * + * \throws std::bad_alloc if out of memory + * std::out_of_range if \c simulationName is not in the database */ +MdpFieldValues +prepareMdpFieldValues(const char *simulationName, + const char *integrator, + const char *tcoupl, + const char *pcoupl); + +/*! \brief Make a string containing an .mdp file from the \c mdpFieldValues. + * + * \throws std::bad_alloc if out of memory */ +std::string +prepareMdpFileContents(const MdpFieldValues &mdpFieldValues); + +} // namespace test +} // namespace gmx + +#endif diff --git a/src/programs/mdrun/tests/mdruncomparisonfixture.h b/src/programs/mdrun/tests/mdruncomparisonfixture.h deleted file mode 100644 index 10de82489a..0000000000 --- a/src/programs/mdrun/tests/mdruncomparisonfixture.h +++ /dev/null @@ -1,143 +0,0 @@ -/* - * This file is part of the GROMACS molecular simulation package. - * - * Copyright (c) 2016, by the GROMACS development team, led by - * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, - * and including many others, as listed in the AUTHORS file in the - * top-level source directory and at http://www.gromacs.org. - * - * GROMACS is free software; you can redistribute it and/or - * modify it under the terms of the GNU Lesser General Public License - * as published by the Free Software Foundation; either version 2.1 - * of the License, or (at your option) any later version. - * - * GROMACS 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 - * Lesser General Public License for more details. - * - * You should have received a copy of the GNU Lesser General Public - * License along with GROMACS; if not, see - * http://www.gnu.org/licenses, or write to the Free Software Foundation, - * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. - * - * If you want to redistribute modifications to GROMACS, please - * consider that scientific software is very special. Version - * control is crucial - bugs must be traceable. We will be happy to - * consider code for inclusion in the official distribution, but - * derived work must not be called official GROMACS. Details are found - * in the README & COPYING files - if they are missing, get the - * official version at http://www.gromacs.org. - * - * To help us fund GROMACS development, we humbly ask that you cite - * the research papers on the package. Check out http://www.gromacs.org. - */ -/*! \libinternal - * - * \brief Functionality for testing whether calls to mdrun produce the - * same energy and force quantities when they should do so. - */ -#ifndef GMX_MDRUN_TESTS_MDRUNCOMPARISONFIXTURE_H -#define GMX_MDRUN_TESTS_MDRUNCOMPARISONFIXTURE_H - -#include -#include -#include -#include - -#include - -#include "testutils/cmdlinetest.h" - -#include "moduletest.h" - -namespace gmx -{ - -namespace test -{ - -class EnergyFrame; -class TrajectoryFrame; -class FloatingPointTolerance; - -/*! \libinternal - * \brief Declares abstract base text fixture class for integration - * tests of mdrun functionality that will compare multiple calls to - * mdrun. - * - * An internal database of several kinds of simulation useful for such - * comparisons is available. - * - * Any method in this class may throw std::bad_alloc if out of memory. - * - * \ingroup module_mdrun_integration_tests - */ -class MdrunComparisonFixture : public MdrunTestFixture -{ - public: - //! Destructor - virtual ~MdrunComparisonFixture(); - //! Helper typedef - typedef std::map MdpFieldValues; - /*! \brief Prepare .mdp values to do a simulation - * - * A database of several kinds of simulation useful for different - * kinds of tests is available. - * - argon12 - * - argon5832 - * - spc5 - * - spc216 - * - alanine_vsite_vacuo - * - alanine_vsite_solvated - * - nonanol - * - * Some of these systems are pretty minimal, because having - * few atoms means few interactions, highly reproducible - * forces, and allows tests to focus on the correctness of the - * implementation of high-level mdrun features. The boxes are - * of a reasonable size so that domain decomposition is - * possible. The pressure-coupling parameters are isotropic, - * and set up so that there will not be dramatic collapse of - * volume over the handful of MD steps that will be run. A - * single temperature-coupling group is used. - * - * This is separate from prepareMdpFile, so that derived - * classes can react to the .mdp settings, e.g. by stopping a - * run after half the steps. - * - * \throws std::bad_alloc if out of memory - * std::out_of_range if \c simulationName is not in the database */ - MdpFieldValues prepareMdpFieldValues(const char *simulationName); - /*! \brief Set up an .mdp file that permits a highly reproducible - * simulation. - * - * \throws std::bad_alloc if out of memory */ - void prepareMdpFile(const MdpFieldValues &mdpFieldValues, - const char *integrator, - const char *tcoupl, - const char *pcoupl); - /*! \brief Run mdrun two ways in a test. Subclasses must override this method. - * - * It is expected that this method calls - * prepareMdpFieldValues() and prepareMdpFile() to help set up - * a call to grompp with gromppCallerRef. Then mdrun will be - * called and perhaps energies and forces compared. */ - virtual void runTest(const CommandLine &gromppCallerRef, - const char *simulationName, - const char *integrator, - const char *tcoupl, - const char *pcoupl, - FloatingPointTolerance tolerance) = 0; - //! Convenience overload of runTest() for cases that don't need to customize the command line for grompp - virtual void runTest(const char *simulationName, - const char *integrator, - const char *tcoupl, - const char *pcoupl, - FloatingPointTolerance tolerance); -}; - -} // namespace test -} // namespace gmx - -#endif -- 2.11.4.GIT