2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
5 * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
6 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
7 * and including many others, as listed in the AUTHORS file in the
8 * top-level source directory and at http://www.gromacs.org.
10 * GROMACS is free software; you can redistribute it and/or
11 * modify it under the terms of the GNU Lesser General Public License
12 * as published by the Free Software Foundation; either version 2.1
13 * of the License, or (at your option) any later version.
15 * GROMACS is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 * Lesser General Public License for more details.
20 * You should have received a copy of the GNU Lesser General Public
21 * License along with GROMACS; if not, see
22 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
23 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
25 * If you want to redistribute modifications to GROMACS, please
26 * consider that scientific software is very special. Version
27 * control is crucial - bugs must be traceable. We will be happy to
28 * consider code for inclusion in the official distribution, but
29 * derived work must not be called official GROMACS. Details are found
30 * in the README & COPYING files - if they are missing, get the
31 * official version at http://www.gromacs.org.
33 * To help us fund GROMACS development, we humbly ask that you cite
34 * the research papers on the package. Check out http://www.gromacs.org.
39 * Tests for the mdrun -rerun functionality
41 * \author Mark Abraham <mark.j.abraham@gmail.com>
42 * \ingroup module_mdrun_integration_tests
48 #include "gromacs/topology/ifunc.h"
49 #include "gromacs/trajectory/energyframe.h"
50 #include "gromacs/utility/stringutil.h"
52 #include "simulatorcomparison.h"
60 /*! \brief Test fixture base for mdrun -rerun
62 * This test ensures mdrun can run a simulation, writing a trajectory
63 * and matching energies, and reproduce the same energies from a rerun
64 * to within a tight tolerance. It says nothing about whether a rerun
65 * can reproduce energies from a trajectory generated with older code,
66 * since that is not a useful property. Whether mdrun produced correct
67 * energies then and now needs different kinds of testing, but if
68 * true, this test ensures the rerun has the expected property.
70 * The limitations of mdrun and its output means that reproducing the
71 * same energies is currently only meaningful for integration without
72 * thermostats or barostats, however the present form of the test
73 * infrastructure has in-principle support for such, if that is ever
76 * We should also not compare pressure, because with constraints the
77 * non-search steps need a much larger tolerance, and per Redmine 1868
78 * we should stop computing pressure in reruns anyway.
80 * Similarly, per 1868, in the present implementation the kinetic
81 * energy quantities are not generally reproducible, either.
83 * The choices for tolerance are arbitrary but sufficient. Rerun does
84 * pair search every frame, so it cannot in general exactly reproduce
85 * quantities from a normal run, because the accumulation order
86 * differs. (Nor does it reproduce pair-search frames exactly,
88 class MdrunRerunTest :
89 public MdrunTestFixture,
90 public ::testing::WithParamInterface<std::tuple<std::string, std::string>>
93 //! Trajectory components to compare
94 static const TrajectoryFrameMatchSettings trajectoryMatchSettings;
97 // Compare box, positions and forces, but not velocities
98 // (velocities are ignored in reruns)
99 const TrajectoryFrameMatchSettings MdrunRerunTest::trajectoryMatchSettings = {
103 ComparisonConditions::MustCompare,
104 ComparisonConditions::NoComparison,
105 ComparisonConditions::MustCompare
108 TEST_P(MdrunRerunTest, WithinTolerances)
110 auto params = GetParam();
111 auto simulationName = std::get<0>(params);
112 auto integrator = std::get<1>(params);
114 formatString("Comparing normal and rerun of simulation '%s' "
115 "with integrator '%s'",
116 simulationName.c_str(), integrator.c_str()));
118 auto mdpFieldValues =
119 prepareMdpFieldValues(simulationName.c_str(), integrator.c_str(), "no", "no");
121 // bd is much less reproducible in a rerun than the other integrators
122 const int toleranceScaleFactor = (integrator == "bd") ? 2 : 1;
123 EnergyTermsToCompare energyTermsToCompare{ {
124 { interaction_function[F_EPOT].longname,
125 relativeToleranceAsPrecisionDependentUlp(10.0, 24 * toleranceScaleFactor,
126 40 * toleranceScaleFactor) },
129 // Specify how trajectory frame matching must work
130 TrajectoryComparison trajectoryComparison{ trajectoryMatchSettings,
131 TrajectoryComparison::s_defaultTrajectoryTolerances };
133 int numWarningsToTolerate = 0;
134 executeRerunTest(&fileManager_, &runner_, simulationName, numWarningsToTolerate, mdpFieldValues,
135 energyTermsToCompare, trajectoryComparison);
138 // TODO The time for OpenCL kernel compilation means these tests time
139 // out. Once that compilation is cached for the whole process, these
140 // tests can run in such configurations.
141 #if GMX_GPU != GMX_GPU_OPENCL
142 INSTANTIATE_TEST_CASE_P(
143 NormalMdrunIsReproduced,
145 ::testing::Combine(::testing::Values("argon12", "tip3p5", "alanine_vsite_vacuo"),
146 ::testing::Values("md", "md-vv", "bd", "sd")));
148 INSTANTIATE_TEST_CASE_P(
149 DISABLED_NormalMdrunIsReproduced,
151 ::testing::Combine(::testing::Values("argon12", "tip3p5", "alanine_vsite_vacuo"),
152 ::testing::Values("md", "md-vv", "bd", "sd")));
155 class MdrunRerunFreeEnergyTest :
156 public MdrunTestFixture,
157 public ::testing::WithParamInterface<std::tuple<std::string, std::string, int>>
161 TEST_P(MdrunRerunFreeEnergyTest, WithinTolerances)
163 auto params = GetParam();
164 auto simulationName = std::get<0>(params);
165 auto integrator = std::get<1>(params);
166 auto initLambdaState = std::get<2>(params);
168 formatString("Comparing normal and rerun of simulation '%s' "
169 "with integrator '%s' for initial lambda state %d",
170 simulationName.c_str(), integrator.c_str(), initLambdaState));
172 auto mdpFieldValues =
173 prepareMdpFieldValues(simulationName.c_str(), integrator.c_str(), "no", "no");
174 mdpFieldValues["other"] += formatString("\ninit-lambda-state = %d", initLambdaState);
176 EnergyTermsToCompare energyTermsToCompare{
177 { { interaction_function[F_EPOT].longname, relativeToleranceAsPrecisionDependentUlp(10.0, 24, 32) },
178 { interaction_function[F_DVDL_COUL].longname, relativeToleranceAsPrecisionDependentUlp(1.0, 8, 8) },
179 { interaction_function[F_DVDL_VDW].longname, relativeToleranceAsPrecisionDependentUlp(1.0, 8, 8) },
180 { interaction_function[F_DVDL_BONDED].longname,
181 relativeToleranceAsPrecisionDependentUlp(1.0, 8, 8) },
182 { interaction_function[F_DVDL_RESTRAINT].longname,
183 relativeToleranceAsPrecisionDependentUlp(1.0, 8, 8) } }
186 // Specify how trajectory frame matching must work
187 TrajectoryComparison trajectoryComparison{ MdrunRerunTest::trajectoryMatchSettings,
188 TrajectoryComparison::s_defaultTrajectoryTolerances };
190 // The md integrator triggers a warning for nearly decoupled
191 // states, which we need to suppress. TODO sometimes?
192 int numWarningsToTolerate = (integrator == "md") ? 1 : 0;
193 executeRerunTest(&fileManager_, &runner_, simulationName, numWarningsToTolerate, mdpFieldValues,
194 energyTermsToCompare, trajectoryComparison);
197 // TODO The time for OpenCL kernel compilation means these tests time
198 // out. Once that compilation is cached for the whole process, these
199 // tests can run in such configurations.
200 #if GMX_GPU != GMX_GPU_OPENCL
201 INSTANTIATE_TEST_CASE_P(MdrunIsReproduced,
202 MdrunRerunFreeEnergyTest,
203 ::testing::Combine(::testing::Values("nonanol_vacuo"),
204 ::testing::Values("md", "md-vv", "sd"),
205 ::testing::Range(0, 11)));
207 INSTANTIATE_TEST_CASE_P(DISABLED_MdrunIsReproduced,
208 MdrunRerunFreeEnergyTest,
209 ::testing::Combine(::testing::Values("nonanol_vacuo"),
210 ::testing::Values("md", "md-vv", "sd"),
211 ::testing::Range(0, 11)));