2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2013,2014,2015,2016,2018 by the GROMACS development team.
5 * Copyright (c) 2019,2020,2021, 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 multi-simulation functionality
41 * \author Mark Abraham <mark.j.abraham@gmail.com>
42 * \ingroup module_mdrun_integration_tests
46 #include "multisimtest.h"
53 #include <gtest/gtest.h>
55 #include "gromacs/mdtypes/md_enums.h"
56 #include "gromacs/utility/basenetwork.h"
57 #include "gromacs/utility/gmxassert.h"
58 #include "gromacs/utility/path.h"
59 #include "gromacs/utility/real.h"
60 #include "gromacs/utility/stringutil.h"
62 #include "testutils/cmdlinetest.h"
63 #include "testutils/mpitest.h"
65 #include "moduletest.h"
66 #include "terminationhelper.h"
73 MultiSimTest::MultiSimTest() :
74 size_(gmx_node_num()),
75 rank_(gmx_node_rank()),
76 numRanksPerSimulation_(std::get<0>(GetParam())),
77 simulationNumber_(rank_ / numRanksPerSimulation_),
78 mdrunCaller_(new CommandLine)
81 // Zero or less ranks doesn't make sense
82 GMX_RELEASE_ASSERT(numRanksPerSimulation_ > 0, "Invalid number of ranks per simulation.");
84 const char* directoryNameFormat = "sim_%d";
86 // Modify the file manager to have a temporary directory unique to
87 // each simulation. No need to have a mutex on this, nobody else
88 // can access the fileManager_ yet because we only just
90 std::string originalTempDirectory = fileManager_.getOutputTempDirectory();
91 std::string newTempDirectory =
92 Path::join(originalTempDirectory, formatString(directoryNameFormat, simulationNumber_));
93 if (rank_ % numRanksPerSimulation_ == 0)
95 // Only one rank per simulation creates directory
96 Directory::create(newTempDirectory);
99 // Make sure directories got created.
100 MPI_Barrier(MdrunTestFixtureBase::communicator_);
102 fileManager_.setOutputTempDirectory(newTempDirectory);
104 mdrunCaller_->append("mdrun");
105 mdrunCaller_->addOption("-multidir");
106 for (int i = 0; i < size_ / numRanksPerSimulation_; ++i)
108 mdrunCaller_->append(Path::join(originalTempDirectory, formatString(directoryNameFormat, i)));
112 bool MultiSimTest::mpiSetupValid() const
114 // Single simulation case is not implemented in multi-sim
115 const bool haveAtLeastTwoSimulations = ((size_ / numRanksPerSimulation_) >= 2);
116 // Mdrun will throw error if simulations don't have identical number of ranks
117 const bool simulationsHaveIdenticalRankNumber = ((size_ % numRanksPerSimulation_) == 0);
119 return (haveAtLeastTwoSimulations && simulationsHaveIdenticalRankNumber);
122 void MultiSimTest::organizeMdpFile(SimulationRunner* runner,
123 IntegrationAlgorithm integrator,
124 TemperatureCoupling tcoupl,
125 PressureCoupling pcoupl,
127 bool doRegression) const
129 GMX_RELEASE_ASSERT(mpiSetupValid(), "Creating the mdp file without valid MPI setup is useless.");
130 const real baseTemperature = 298;
131 const real basePressure = 1;
132 std::string mdpFileContents = formatString(
138 "nstcalcenergy = 1\n"
142 // pressure coupling (if active)
145 "compressibility = 4.5e-5\n"
146 // velocity generation
150 // v-rescale and c-rescale use ld-seed
152 // Two systems are used: spc2 non-interacting also at cutoff 1nm
153 // tip3p5 has box length of 1.86
156 // Trajectory output if required
161 enumValueToString(integrator),
162 enumValueToString(tcoupl),
163 enumValueToString(pcoupl),
165 baseTemperature + 0.0001 * rank_,
166 basePressure * std::pow(1.01, rank_),
167 /* Set things up so that the initial KE decreases with
168 increasing replica number, so that the (identical)
169 starting PE decreases on the first step more for the
170 replicas with higher number, which will tend to force
171 replica exchange to occur. */
172 std::max(baseTemperature - 10 * rank_, real(0)),
173 // If we do regression, we need reproducible velocity
174 // generation, which can be different per simulation
175 (doRegression ? 671324 + simulationNumber_ : -1),
176 // If we do regression, we need reproducible temperature and
177 // pressure coupling, which can be different per simulation
178 (doRegression ? 51203 + simulationNumber_ : -1),
179 // If we do regression, write one intermediate point
180 (doRegression ? int(numSteps / 2) : 0),
181 (doRegression ? int(numSteps / 2) : 0),
182 (doRegression ? int(numSteps / 2) : 0),
183 // If we do regression, print energies every step so
184 // we're sure to catch the replica exchange steps
185 (doRegression ? 1 : 1000));
186 runner->useStringAsMdpFile(mdpFileContents);
189 void MultiSimTest::runGrompp(SimulationRunner* runner, int numSteps, bool doRegression, int maxWarnings) const
191 // Call grompp once per simulation
192 if (rank_ % numRanksPerSimulation_ == 0)
194 const auto& simulator = std::get<1>(GetParam());
195 const auto& tcoupl = std::get<2>(GetParam());
196 const auto& pcoupl = std::get<3>(GetParam());
197 organizeMdpFile(runner, simulator, tcoupl, pcoupl, numSteps, doRegression);
199 caller.addOption("-maxwarn", maxWarnings);
200 EXPECT_EQ(0, runner->callGromppOnThisRank(caller));
204 // Make sure simulation masters have written the .tpr file before other ranks try to read it.
205 MPI_Barrier(MdrunTestFixtureBase::communicator_);
209 void MultiSimTest::runExitsNormallyTest()
211 if (!mpiSetupValid())
213 // Can't test multi-sim without multiple simulations
217 SimulationRunner runner(&fileManager_);
218 runner.useTopGroAndNdxFromDatabase("spc2");
222 ASSERT_EQ(0, runner.callMdrun(*mdrunCaller_));
225 void MultiSimTest::runMaxhTest()
227 if (!mpiSetupValid())
229 // Can't test multi-sim without multiple simulations
233 SimulationRunner runner(&fileManager_);
234 runner.useTopGroAndNdxFromDatabase("spc2");
236 TerminationHelper helper(&fileManager_, mdrunCaller_.get(), &runner);
237 // Make sure -maxh has a chance to propagate
239 runGrompp(&runner, numSteps);
241 helper.runFirstMdrun(runner.cptFileName_);
242 helper.runSecondMdrun();