src/programs/mdrun/tests/replicaexchange.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   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.
   9  *
  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.
  14  *
  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.
  19  *
  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.
  24  *
  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.
  32  *
  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.
  35  */
  36
  37 /*! \internal \file
  38  * \brief
  39  * Tests for the mdrun replica-exchange functionality
  40  *
  41  * \author Mark Abraham <mark.j.abraham@gmail.com>
  42  * \ingroup module_mdrun_integration_tests
  43  */
  44 #include "gmxpre.h"
  45
  46 #include "config.h"
  47
  48 #include <regex>
  49
  50 #include <gtest/gtest.h>
  51
  52 #include "gromacs/mdtypes/md_enums.h"
  53 #include "gromacs/topology/ifunc.h"
  54 #include "gromacs/utility/basenetwork.h"
  55 #include "gromacs/utility/filestream.h"
  56 #include "gromacs/utility/path.h"
  57 #include "gromacs/utility/stringutil.h"
  58
  59 #include "testutils/refdata.h"
  60 #include "testutils/testfilemanager.h"
  61
  62 #include "energycomparison.h"
  63 #include "multisimtest.h"
  64 #include "trajectorycomparison.h"
  65
  66 namespace gmx
  67 {
  68 namespace test
  69 {
  70
  71 //! Convenience typedef
  72 typedef MultiSimTest ReplicaExchangeEnsembleTest;
  73
  74 TEST_P(ReplicaExchangeEnsembleTest, ExitsNormally)
  75 {
  76     mdrunCaller_->addOption("-replex", 1);
  77     runExitsNormallyTest();
  78 }
  79
  80 /* Note, not all preprocessor implementations nest macro expansions
  81    the same way / at all, if we would try to duplicate less code. */
  82
  83 #if GMX_LIB_MPI
  84 INSTANTIATE_TEST_SUITE_P(
  85         WithDifferentControlVariables,
  86         ReplicaExchangeEnsembleTest,
  87         ::testing::Combine(::testing::Values(NumRanksPerSimulation(1), NumRanksPerSimulation(2)),
  88                            ::testing::Values(IntegrationAlgorithm::MD),
  89                            ::testing::Values(TemperatureCoupling::VRescale),
  90                            ::testing::Values(PressureCoupling::No, PressureCoupling::Berendsen)));
  91 #else
  92 INSTANTIATE_TEST_SUITE_P(
  93         DISABLED_WithDifferentControlVariables,
  94         ReplicaExchangeEnsembleTest,
  95         ::testing::Combine(::testing::Values(NumRanksPerSimulation(1), NumRanksPerSimulation(2)),
  96                            ::testing::Values(IntegrationAlgorithm::MD),
  97                            ::testing::Values(TemperatureCoupling::VRescale),
  98                            ::testing::Values(PressureCoupling::No, PressureCoupling::Berendsen)));
  99 #endif
 100
 101 //! Convenience typedef
 102 typedef MultiSimTest ReplicaExchangeTerminationTest;
 103
 104 TEST_P(ReplicaExchangeTerminationTest, WritesCheckpointAfterMaxhTerminationAndThenRestarts)
 105 {
 106     mdrunCaller_->addOption("-replex", 1);
 107     runMaxhTest();
 108 }
 109
 110 INSTANTIATE_TEST_SUITE_P(InNvt,
 111                          ReplicaExchangeTerminationTest,
 112                          ::testing::Combine(::testing::Values(NumRanksPerSimulation(1)),
 113                                             ::testing::Values(IntegrationAlgorithm::MD),
 114                                             ::testing::Values(TemperatureCoupling::VRescale),
 115                                             ::testing::Values(PressureCoupling::No)));
 116
 117 /*! \brief Return replica exchange related output from logfile
 118  *
 119  * All replica exchange related output in log files start with 'Repl',
 120  * making extraction easy. This function also removes the printing of
 121  * energy differences, as the log files are compared exactly, and
 122  * energy differences will slightly vary between runs.
 123  *
 124  * \param logFileName  Name of log file
 125  * \return  Replica exchange related output in log file
 126  */
 127 static std::string getReplicaExchangeOutputFromLogFile(const std::string& logFileName)
 128 {
 129     TextInputFile logFile(logFileName);
 130     std::string   replExOutput;
 131     std::string   line;
 132     while (logFile.readLine(&line))
 133     {
 134         // All replica exchange output lines starts with "Repl"
 135         if (startsWith(line, "Repl"))
 136         {
 137             // This is an exact comparison, so we can't compare the energies which
 138             // are slightly different per run. Energies are tested later.
 139             const auto pos = line.find("dE_term");
 140             if (pos != std::string::npos)
 141             {
 142                 line.replace(line.begin() + pos, line.end(), "[ not checked ]\n");
 143             }
 144             replExOutput.append(line);
 145         }
 146     }
 147     return replExOutput;
 148 }
 149
 150 //! Convenience typedef
 151 typedef MultiSimTest ReplicaExchangeRegressionTest;
 152
 153 /* Run replica exchange simulations, compare to reference data
 154  *
 155  * Reference data generated by
 156  *
 157  * GROMACS version:    2022-dev
 158  * Precision:          single and double (separate reference data)
 159  * Memory model:       64 bit
 160  * MPI library:        MPI
 161  * OpenMP support:     enabled (GMX_OPENMP_MAX_THREADS = 64)
 162  * GPU support:        disabled
 163  * SIMD instructions:  AVX2_256
 164  * FFT library:        fftw-3.3.9-sse2-avx
 165  * RDTSCP usage:       enabled
 166  * TNG support:        enabled
 167  * Hwloc support:      disabled
 168  * Tracing support:    disabled
 169  * C compiler:         /usr/local/bin/mpicc Clang 8.0.1
 170  * C++ compiler:       /usr/local/bin/mpic++ Clang 8.0.1
 171  *
 172  */
 173 TEST_P(ReplicaExchangeRegressionTest, WithinTolerances)
 174 {
 175     if (!mpiSetupValid())
 176     {
 177         // Can't test multi-sim without multiple simulations
 178         return;
 179     }
 180
 181     if (size_ != 4)
 182     {
 183         // Results are depending on number of ranks, and we can't have reference
 184         // data for all cases. Restricting the regression tests to runs with 4 ranks.
 185         // This allows testing 4 replicas with single rank, or 2 replicas with 2 ranks each.
 186         return;
 187     }
 188     const auto& tcoupl = std::get<2>(GetParam());
 189     const auto& pcoupl = std::get<3>(GetParam());
 190
 191     const int numSteps       = 16;
 192     const int exchangePeriod = 4;
 193     // grompp warns about generating velocities and using parrinello-rahman
 194     const int maxWarnings = (pcoupl == PressureCoupling::ParrinelloRahman ? 1 : 0);
 195
 196     mdrunCaller_->addOption("-replex", exchangePeriod);
 197     // Seeds need to be reproducible for regression, but can be different per simulation
 198     mdrunCaller_->addOption("-reseed", 98713 + simulationNumber_);
 199
 200     SimulationRunner runner(&fileManager_);
 201     runner.useTopGroAndNdxFromDatabase("tip3p5");
 202
 203     runGrompp(&runner, numSteps, true, maxWarnings);
 204     ASSERT_EQ(0, runner.callMdrun(*mdrunCaller_));
 205
 206 #if GMX_LIB_MPI
 207     // Make sure all simulations are finished before checking the results.
 208     MPI_Barrier(MdrunTestFixtureBase::communicator_);
 209 #endif
 210
 211     // We only test simulation results on one rank to avoid problems with reference file access.
 212     if (rank_ == 0)
 213     {
 214         // Create reference data helper object
 215         TestReferenceData refData;
 216
 217         // Specify how energy trajectory comparison must work
 218         const auto hasConservedField =
 219                 !(tcoupl == TemperatureCoupling::No && pcoupl == PressureCoupling::No);
 220         // Tolerances copied from simulator tests
 221         EnergyTermsToCompare energyTermsToCompare{ {
 222                 { interaction_function[F_EPOT].longname,
 223                   relativeToleranceAsPrecisionDependentUlp(60.0, 200, 160) },
 224                 { interaction_function[F_EKIN].longname,
 225                   relativeToleranceAsPrecisionDependentUlp(60.0, 200, 160) },
 226         } };
 227         if (hasConservedField)
 228         {
 229             energyTermsToCompare.emplace(interaction_function[F_ECONSERVED].longname,
 230                                          relativeToleranceAsPrecisionDependentUlp(50.0, 100, 80));
 231         }
 232         if (pcoupl != PressureCoupling::No)
 233         {
 234             energyTermsToCompare.emplace("Volume",
 235                                          relativeToleranceAsPrecisionDependentUlp(10.0, 200, 160));
 236         }
 237
 238         // Specify how trajectory frame matching must work.
 239         const TrajectoryFrameMatchSettings trajectoryMatchSettings{ true,
 240                                                                     true,
 241                                                                     true,
 242                                                                     ComparisonConditions::MustCompare,
 243                                                                     ComparisonConditions::MustCompare,
 244                                                                     ComparisonConditions::MustCompare,
 245                                                                     MaxNumFrames::compareAllFrames() };
 246         TrajectoryTolerances trajectoryTolerances = TrajectoryComparison::s_defaultTrajectoryTolerances;
 247         // By default, velocity tolerance is MUCH tighter than force tolerance
 248         trajectoryTolerances.velocities = trajectoryTolerances.forces;
 249         // Build the functor that will compare reference and test
 250         // trajectory frames in the chosen way.
 251         TrajectoryComparison trajectoryComparison{ trajectoryMatchSettings, trajectoryTolerances };
 252
 253         // Loop over simulations
 254         for (int simulationNumber = 0; simulationNumber < (size_ / numRanksPerSimulation_);
 255              simulationNumber++)
 256         {
 257             TestReferenceChecker simulationChecker(refData.rootChecker().checkCompound(
 258                     "Simulation", formatString("Replica %d", simulationNumber)));
 259
 260             const auto logFileName =
 261                     std::regex_replace(runner.logFileName_,
 262                                        std::regex(formatString("sim_%d", simulationNumber_)),
 263                                        formatString("sim_%d", simulationNumber));
 264             const auto energyFileName =
 265                     std::regex_replace(runner.edrFileName_,
 266                                        std::regex(formatString("sim_%d", simulationNumber_)),
 267                                        formatString("sim_%d", simulationNumber));
 268             const auto trajectoryFileName =
 269                     std::regex_replace(runner.fullPrecisionTrajectoryFileName_,
 270                                        std::regex(formatString("sim_%d", simulationNumber_)),
 271                                        formatString("sim_%d", simulationNumber));
 272
 273             // Check log replica exchange related output (contains exchange statistics)
 274             auto replicaExchangeOutputChecker =
 275                     simulationChecker.checkCompound("ReplExOutput", "Output");
 276             const auto replExOutput = getReplicaExchangeOutputFromLogFile(logFileName);
 277             replicaExchangeOutputChecker.checkTextBlock(replExOutput, "Replica Exchange Output");
 278
 279             // Check that the energies agree with the refdata within tolerance.
 280             checkEnergiesAgainstReferenceData(energyFileName, energyTermsToCompare, &simulationChecker);
 281
 282             // Check that the trajectories agree with the refdata within tolerance.
 283             checkTrajectoryAgainstReferenceData(trajectoryFileName, trajectoryComparison, &simulationChecker);
 284
 285         } // end loop over simulations
 286     }     // end testing simulations on one rank
 287
 288 #if GMX_LIB_MPI
 289     // Make sure testing is complete before returning - ranks delete temporary files on exit
 290     MPI_Barrier(MdrunTestFixtureBase::communicator_);
 291 #endif
 292 }
 293
 294 /*! \brief Helper struct printing custom test name
 295  *
 296  * Regression test results not only depend on the test parameters, but
 297  * also on the total number of ranks and the precision. Names must
 298  * reflect that to identify correct reference data.
 299  */
 300 struct PrintReplicaExchangeParametersToString
 301 {
 302     template<class ParamType>
 303     std::string operator()(const testing::TestParamInfo<ParamType>& parameter) const
 304     {
 305         auto testIdentifier =
 306                 formatString("ReplExRegression_%s_%s_%s_%dRanks_%dRanksPerSimulation_%s",
 307                              enumValueToString(std::get<1>(parameter.param)),
 308                              enumValueToString(std::get<2>(parameter.param)),
 309                              enumValueToString(std::get<3>(parameter.param)),
 310                              gmx_node_num(),
 311                              static_cast<int>(std::get<0>(parameter.param)),
 312                              GMX_DOUBLE ? "d" : "s");
 313         // Valid GTest names cannot include hyphens
 314         testIdentifier.erase(std::remove(testIdentifier.begin(), testIdentifier.end(), '-'),
 315                              testIdentifier.end());
 316         return testIdentifier;
 317     }
 318 };
 319
 320 #if GMX_LIB_MPI
 321 INSTANTIATE_TEST_SUITE_P(
 322         ReplicaExchangeIsEquivalentToReferenceLeapFrog,
 323         ReplicaExchangeRegressionTest,
 324         ::testing::Combine(::testing::Values(NumRanksPerSimulation(1), NumRanksPerSimulation(2)),
 325                            ::testing::Values(IntegrationAlgorithm::MD),
 326                            ::testing::Values(TemperatureCoupling::VRescale, TemperatureCoupling::NoseHoover),
 327                            ::testing::Values(PressureCoupling::CRescale, PressureCoupling::ParrinelloRahman)),
 328         PrintReplicaExchangeParametersToString());
 329 INSTANTIATE_TEST_SUITE_P(ReplicaExchangeIsEquivalentToReferenceVelocityVerlet,
 330                          ReplicaExchangeRegressionTest,
 331                          ::testing::Combine(::testing::Values(NumRanksPerSimulation(1),
 332                                                               NumRanksPerSimulation(2)),
 333                                             ::testing::Values(IntegrationAlgorithm::VV),
 334                                             ::testing::Values(TemperatureCoupling::NoseHoover),
 335                                             ::testing::Values(PressureCoupling::No)),
 336                          PrintReplicaExchangeParametersToString());
 337 #else
 338 INSTANTIATE_TEST_SUITE_P(
 339         DISABLED_ReplicaExchangeIsEquivalentToReferenceLeapFrog,
 340         ReplicaExchangeRegressionTest,
 341         ::testing::Combine(::testing::Values(NumRanksPerSimulation(1), NumRanksPerSimulation(2)),
 342                            ::testing::Values(IntegrationAlgorithm::MD),
 343                            ::testing::Values(TemperatureCoupling::VRescale, TemperatureCoupling::NoseHoover),
 344                            ::testing::Values(PressureCoupling::CRescale, PressureCoupling::ParrinelloRahman)),
 345         PrintReplicaExchangeParametersToString());
 346 INSTANTIATE_TEST_SUITE_P(DISABLED_ReplicaExchangeIsEquivalentToReferenceVelocityVerlet,
 347                          ReplicaExchangeRegressionTest,
 348                          ::testing::Combine(::testing::Values(NumRanksPerSimulation(1),
 349                                                               NumRanksPerSimulation(2)),
 350                                             ::testing::Values(IntegrationAlgorithm::VV),
 351                                             ::testing::Values(TemperatureCoupling::NoseHoover),
 352                                             ::testing::Values(PressureCoupling::No)),
 353                          PrintReplicaExchangeParametersToString());
 354 #endif
 355 } // namespace test
 356 } // namespace gmx