2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2019, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
38 * Tests to compare two simulators which are expected to be identical
40 * \author Mark Abraham <mark.j.abraham@gmail.com>
41 * \author Pascal Merz <pascal.merz@me.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"
61 /*! \brief Test fixture base for two equivalent simulators
63 * This test ensures that two simulator code paths (called via different mdp
64 * options and/or environment variables) yield identical coordinate, velocity,
65 * box, force and energy trajectories, up to some (arbitrary) precision.
67 * These tests are useful to check that re-implementations of existing simulators
68 * are correct, and that different code paths expected to yield identical results
71 using SimulatorComparisonTestParams =
72 std::tuple<std::tuple<std::string, std::string, std::string, std::string>, std::string>;
73 class SimulatorComparisonTest :
74 public MdrunTestFixture,
75 public ::testing::WithParamInterface<SimulatorComparisonTestParams>
79 TEST_P(SimulatorComparisonTest, WithinTolerances)
81 auto params = GetParam();
82 auto mdpParams = std::get<0>(params);
83 auto simulationName = std::get<0>(mdpParams);
84 auto integrator = std::get<1>(mdpParams);
85 auto tcoupling = std::get<2>(mdpParams);
86 auto pcoupling = std::get<3>(mdpParams);
87 auto envVariable = std::get<1>(params);
89 if (integrator == "md-vv" && pcoupling == "Parrinello-Rahman")
91 // do_md calls this MTTK, requires Nose-Hoover, and
92 // does not work with constraints or anisotropically
96 SCOPED_TRACE(formatString(
97 "Comparing two simulations of '%s' "
98 "with integrator '%s' and '%s' temperature coupling, "
99 "switching environment variable '%s'",
100 simulationName.c_str(), integrator.c_str(), tcoupling.c_str(), envVariable.c_str()));
102 auto mdpFieldValues = prepareMdpFieldValues(simulationName.c_str(), integrator.c_str(),
103 tcoupling.c_str(), pcoupling.c_str());
105 EnergyTermsToCompare energyTermsToCompare{ {
106 { interaction_function[F_EPOT].longname, relativeToleranceAsPrecisionDependentUlp(10.0, 100, 40) },
107 { interaction_function[F_EKIN].longname, relativeToleranceAsPrecisionDependentUlp(60.0, 100, 40) },
108 { interaction_function[F_PRES].longname,
109 relativeToleranceAsPrecisionDependentFloatingPoint(10.0, 0.01, 0.001) },
112 if (simulationName == "argon12")
114 // Without constraints, we can be more strict
115 energyTermsToCompare = { {
116 { interaction_function[F_EPOT].longname,
117 relativeToleranceAsPrecisionDependentUlp(10.0, 24, 40) },
118 { interaction_function[F_EKIN].longname,
119 relativeToleranceAsPrecisionDependentUlp(10.0, 24, 40) },
120 { interaction_function[F_PRES].longname,
121 relativeToleranceAsPrecisionDependentFloatingPoint(10.0, 0.001, 0.0001) },
125 // Specify how trajectory frame matching must work.
126 TrajectoryFrameMatchSettings trajectoryMatchSettings{ true,
129 ComparisonConditions::MustCompare,
130 ComparisonConditions::MustCompare,
131 ComparisonConditions::MustCompare };
132 TrajectoryTolerances trajectoryTolerances = TrajectoryComparison::s_defaultTrajectoryTolerances;
133 if (simulationName != "argon12")
135 trajectoryTolerances.velocities = trajectoryTolerances.coordinates;
138 // Build the functor that will compare reference and test
139 // trajectory frames in the chosen way.
140 TrajectoryComparison trajectoryComparison{ trajectoryMatchSettings, trajectoryTolerances };
142 int numWarningsToTolerate = 0;
143 executeSimulatorComparisonTest(envVariable, &fileManager_, &runner_, simulationName,
144 numWarningsToTolerate, mdpFieldValues, energyTermsToCompare,
145 trajectoryComparison);
148 // TODO: The time for OpenCL kernel compilation means these tests time
149 // out. Once that compilation is cached for the whole process, these
150 // tests can run in such configurations.
151 // These tests are very sensitive, so we only run them in double precision.
152 // As we change call ordering, they might actually become too strict to be useful.
153 #if GMX_GPU != GMX_GPU_OPENCL && GMX_DOUBLE
154 INSTANTIATE_TEST_CASE_P(SimulatorsAreEquivalentDefaultModular,
155 SimulatorComparisonTest,
156 ::testing::Combine(::testing::Combine(::testing::Values("argon12", "tip3p5"),
157 ::testing::Values("md-vv"),
158 ::testing::Values("no", "v-rescale"),
159 ::testing::Values("no")),
160 ::testing::Values("GMX_DISABLE_MODULAR_SIMULATOR")));
161 INSTANTIATE_TEST_CASE_P(
162 SimulatorsAreEquivalentDefaultLegacy,
163 SimulatorComparisonTest,
164 ::testing::Combine(::testing::Combine(::testing::Values("argon12", "tip3p5"),
165 ::testing::Values("md"),
166 ::testing::Values("no", "v-rescale"),
167 ::testing::Values("no", "Parrinello-Rahman")),
168 ::testing::Values("GMX_USE_MODULAR_SIMULATOR")));
170 INSTANTIATE_TEST_CASE_P(DISABLED_SimulatorsAreEquivalentDefaultModular,
171 SimulatorComparisonTest,
172 ::testing::Combine(::testing::Combine(::testing::Values("argon12", "tip3p5"),
173 ::testing::Values("md-vv"),
174 ::testing::Values("no", "v-rescale"),
175 ::testing::Values("no")),
176 ::testing::Values("GMX_DISABLE_MODULAR_SIMULATOR")));
177 INSTANTIATE_TEST_CASE_P(
178 DISABLED_SimulatorsAreEquivalentDefaultLegacy,
179 SimulatorComparisonTest,
180 ::testing::Combine(::testing::Combine(::testing::Values("argon12", "tip3p5"),
181 ::testing::Values("md"),
182 ::testing::Values("no", "v-rescale"),
183 ::testing::Values("no", "Parrinello-Rahman")),
184 ::testing::Values("GMX_USE_MODULAR_SIMULATOR")));