src/gromacs/mdrun/runner.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
   5  * Copyright (c) 2001-2004, The GROMACS development team.
   6  * Copyright (c) 2011-2019,2020, by the GROMACS development team, led by
   7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   8  * and including many others, as listed in the AUTHORS file in the
   9  * top-level source directory and at http://www.gromacs.org.
  10  *
  11  * GROMACS is free software; you can redistribute it and/or
  12  * modify it under the terms of the GNU Lesser General Public License
  13  * as published by the Free Software Foundation; either version 2.1
  14  * of the License, or (at your option) any later version.
  15  *
  16  * GROMACS is distributed in the hope that it will be useful,
  17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19  * Lesser General Public License for more details.
  20  *
  21  * You should have received a copy of the GNU Lesser General Public
  22  * License along with GROMACS; if not, see
  23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  25  *
  26  * If you want to redistribute modifications to GROMACS, please
  27  * consider that scientific software is very special. Version
  28  * control is crucial - bugs must be traceable. We will be happy to
  29  * consider code for inclusion in the official distribution, but
  30  * derived work must not be called official GROMACS. Details are found
  31  * in the README & COPYING files - if they are missing, get the
  32  * official version at http://www.gromacs.org.
  33  *
  34  * To help us fund GROMACS development, we humbly ask that you cite
  35  * the research papers on the package. Check out http://www.gromacs.org.
  36  */
  37 /*! \internal \file
  38  *
  39  * \brief Implements the MD runner routine calling all integrators.
  40  *
  41  * \author David van der Spoel <david.vanderspoel@icm.uu.se>
  42  * \ingroup module_mdrun
  43  */
  44 #include "gmxpre.h"
  45
  46 #include "runner.h"
  47
  48 #include "config.h"
  49
  50 #include <cassert>
  51 #include <cinttypes>
  52 #include <csignal>
  53 #include <cstdlib>
  54 #include <cstring>
  55
  56 #include <algorithm>
  57 #include <memory>
  58
  59 #include "gromacs/commandline/filenm.h"
  60 #include "gromacs/domdec/builder.h"
  61 #include "gromacs/domdec/domdec.h"
  62 #include "gromacs/domdec/domdec_struct.h"
  63 #include "gromacs/domdec/gpuhaloexchange.h"
  64 #include "gromacs/domdec/localatomsetmanager.h"
  65 #include "gromacs/domdec/partition.h"
  66 #include "gromacs/ewald/ewald_utils.h"
  67 #include "gromacs/ewald/pme_gpu_program.h"
  68 #include "gromacs/ewald/pme_only.h"
  69 #include "gromacs/ewald/pme_pp_comm_gpu.h"
  70 #include "gromacs/fileio/checkpoint.h"
  71 #include "gromacs/fileio/gmxfio.h"
  72 #include "gromacs/fileio/oenv.h"
  73 #include "gromacs/fileio/tpxio.h"
  74 #include "gromacs/gmxlib/network.h"
  75 #include "gromacs/gmxlib/nrnb.h"
  76 #include "gromacs/gpu_utils/device_context.h"
  77 #include "gromacs/gpu_utils/device_stream_manager.h"
  78 #include "gromacs/gpu_utils/gpu_utils.h"
  79 #include "gromacs/hardware/cpuinfo.h"
  80 #include "gromacs/hardware/detecthardware.h"
  81 #include "gromacs/hardware/printhardware.h"
  82 #include "gromacs/imd/imd.h"
  83 #include "gromacs/listed_forces/disre.h"
  84 #include "gromacs/listed_forces/gpubonded.h"
  85 #include "gromacs/listed_forces/orires.h"
  86 #include "gromacs/math/functions.h"
  87 #include "gromacs/math/utilities.h"
  88 #include "gromacs/math/vec.h"
  89 #include "gromacs/mdlib/boxdeformation.h"
  90 #include "gromacs/mdlib/broadcaststructs.h"
  91 #include "gromacs/mdlib/calc_verletbuf.h"
  92 #include "gromacs/mdlib/dispersioncorrection.h"
  93 #include "gromacs/mdlib/enerdata_utils.h"
  94 #include "gromacs/mdlib/force.h"
  95 #include "gromacs/mdlib/forcerec.h"
  96 #include "gromacs/mdlib/gmx_omp_nthreads.h"
  97 #include "gromacs/mdlib/makeconstraints.h"
  98 #include "gromacs/mdlib/md_support.h"
  99 #include "gromacs/mdlib/mdatoms.h"
 100 #include "gromacs/mdlib/membed.h"
 101 #include "gromacs/mdlib/sighandler.h"
 102 #include "gromacs/mdlib/stophandler.h"
 103 #include "gromacs/mdlib/tgroup.h"
 104 #include "gromacs/mdlib/updategroups.h"
 105 #include "gromacs/mdlib/vsite.h"
 106 #include "gromacs/mdrun/mdmodules.h"
 107 #include "gromacs/mdrun/simulationcontext.h"
 108 #include "gromacs/mdrunutility/handlerestart.h"
 109 #include "gromacs/mdrunutility/logging.h"
 110 #include "gromacs/mdrunutility/multisim.h"
 111 #include "gromacs/mdrunutility/printtime.h"
 112 #include "gromacs/mdrunutility/threadaffinity.h"
 113 #include "gromacs/mdtypes/commrec.h"
 114 #include "gromacs/mdtypes/enerdata.h"
 115 #include "gromacs/mdtypes/fcdata.h"
 116 #include "gromacs/mdtypes/forcerec.h"
 117 #include "gromacs/mdtypes/group.h"
 118 #include "gromacs/mdtypes/inputrec.h"
 119 #include "gromacs/mdtypes/interaction_const.h"
 120 #include "gromacs/mdtypes/md_enums.h"
 121 #include "gromacs/mdtypes/mdatom.h"
 122 #include "gromacs/mdtypes/mdrunoptions.h"
 123 #include "gromacs/mdtypes/observableshistory.h"
 124 #include "gromacs/mdtypes/simulation_workload.h"
 125 #include "gromacs/mdtypes/state.h"
 126 #include "gromacs/mdtypes/state_propagator_data_gpu.h"
 127 #include "gromacs/nbnxm/gpu_data_mgmt.h"
 128 #include "gromacs/nbnxm/nbnxm.h"
 129 #include "gromacs/nbnxm/pairlist_tuning.h"
 130 #include "gromacs/pbcutil/pbc.h"
 131 #include "gromacs/pulling/output.h"
 132 #include "gromacs/pulling/pull.h"
 133 #include "gromacs/pulling/pull_rotation.h"
 134 #include "gromacs/restraint/manager.h"
 135 #include "gromacs/restraint/restraintmdmodule.h"
 136 #include "gromacs/restraint/restraintpotential.h"
 137 #include "gromacs/swap/swapcoords.h"
 138 #include "gromacs/taskassignment/decidegpuusage.h"
 139 #include "gromacs/taskassignment/decidesimulationworkload.h"
 140 #include "gromacs/taskassignment/resourcedivision.h"
 141 #include "gromacs/taskassignment/taskassignment.h"
 142 #include "gromacs/taskassignment/usergpuids.h"
 143 #include "gromacs/timing/gpu_timing.h"
 144 #include "gromacs/timing/wallcycle.h"
 145 #include "gromacs/timing/wallcyclereporting.h"
 146 #include "gromacs/topology/mtop_util.h"
 147 #include "gromacs/trajectory/trajectoryframe.h"
 148 #include "gromacs/utility/basenetwork.h"
 149 #include "gromacs/utility/cstringutil.h"
 150 #include "gromacs/utility/exceptions.h"
 151 #include "gromacs/utility/fatalerror.h"
 152 #include "gromacs/utility/filestream.h"
 153 #include "gromacs/utility/gmxassert.h"
 154 #include "gromacs/utility/gmxmpi.h"
 155 #include "gromacs/utility/keyvaluetree.h"
 156 #include "gromacs/utility/logger.h"
 157 #include "gromacs/utility/loggerbuilder.h"
 158 #include "gromacs/utility/mdmodulenotification.h"
 159 #include "gromacs/utility/physicalnodecommunicator.h"
 160 #include "gromacs/utility/pleasecite.h"
 161 #include "gromacs/utility/programcontext.h"
 162 #include "gromacs/utility/smalloc.h"
 163 #include "gromacs/utility/stringutil.h"
 164
 165 #include "isimulator.h"
 166 #include "replicaexchange.h"
 167 #include "simulatorbuilder.h"
 168
 169 #if GMX_FAHCORE
 170 #    include "corewrap.h"
 171 #endif
 172
 173 namespace gmx
 174 {
 175
 176
 177 /*! \brief Manage any development feature flag variables encountered
 178  *
 179  * The use of dev features indicated by environment variables is
 180  * logged in order to ensure that runs with such features enabled can
 181  * be identified from their log and standard output. Any cross
 182  * dependencies are also checked, and if unsatisfied, a fatal error
 183  * issued.
 184  *
 185  * Note that some development features overrides are applied already here:
 186  * the GPU communication flags are set to false in non-tMPI and non-CUDA builds.
 187  *
 188  * \param[in]  mdlog                Logger object.
 189  * \param[in]  useGpuForNonbonded   True if the nonbonded task is offloaded in this run.
 190  * \param[in]  pmeRunMode           The PME run mode for this run
 191  * \returns                         The object populated with development feature flags.
 192  */
 193 static DevelopmentFeatureFlags manageDevelopmentFeatures(const gmx::MDLogger& mdlog,
 194                                                          const bool           useGpuForNonbonded,
 195                                                          const PmeRunMode     pmeRunMode)
 196 {
 197     DevelopmentFeatureFlags devFlags;
 198
 199     // Some builds of GCC 5 give false positive warnings that these
 200     // getenv results are ignored when clearly they are used.
 201 #pragma GCC diagnostic push
 202 #pragma GCC diagnostic ignored "-Wunused-result"
 203     devFlags.enableGpuBufferOps = (getenv("GMX_USE_GPU_BUFFER_OPS") != nullptr)
 204                                   && (GMX_GPU == GMX_GPU_CUDA) && useGpuForNonbonded;
 205     devFlags.forceGpuUpdateDefault = (getenv("GMX_FORCE_UPDATE_DEFAULT_GPU") != nullptr);
 206     devFlags.enableGpuHaloExchange =
 207             (getenv("GMX_GPU_DD_COMMS") != nullptr && GMX_THREAD_MPI && (GMX_GPU == GMX_GPU_CUDA));
 208     devFlags.enableGpuPmePPComm =
 209             (getenv("GMX_GPU_PME_PP_COMMS") != nullptr && GMX_THREAD_MPI && (GMX_GPU == GMX_GPU_CUDA));
 210 #pragma GCC diagnostic pop
 211
 212     if (devFlags.enableGpuBufferOps)
 213     {
 214         GMX_LOG(mdlog.warning)
 215                 .asParagraph()
 216                 .appendTextFormatted(
 217                         "This run uses the 'GPU buffer ops' feature, enabled by the "
 218                         "GMX_USE_GPU_BUFFER_OPS environment variable.");
 219     }
 220
 221     if (devFlags.forceGpuUpdateDefault)
 222     {
 223         GMX_LOG(mdlog.warning)
 224                 .asParagraph()
 225                 .appendTextFormatted(
 226                         "This run will default to '-update gpu' as requested by the "
 227                         "GMX_FORCE_UPDATE_DEFAULT_GPU environment variable. GPU update with domain "
 228                         "decomposition lacks substantial testing and should be used with caution.");
 229     }
 230
 231     if (devFlags.enableGpuHaloExchange)
 232     {
 233         if (useGpuForNonbonded)
 234         {
 235             if (!devFlags.enableGpuBufferOps)
 236             {
 237                 GMX_LOG(mdlog.warning)
 238                         .asParagraph()
 239                         .appendTextFormatted(
 240                                 "Enabling GPU buffer operations required by GMX_GPU_DD_COMMS "
 241                                 "(equivalent with GMX_USE_GPU_BUFFER_OPS=1).");
 242                 devFlags.enableGpuBufferOps = true;
 243             }
 244             GMX_LOG(mdlog.warning)
 245                     .asParagraph()
 246                     .appendTextFormatted(
 247                             "This run has requested the 'GPU halo exchange' feature, enabled by "
 248                             "the "
 249                             "GMX_GPU_DD_COMMS environment variable.");
 250         }
 251         else
 252         {
 253             GMX_LOG(mdlog.warning)
 254                     .asParagraph()
 255                     .appendTextFormatted(
 256                             "GMX_GPU_DD_COMMS environment variable detected, but the 'GPU "
 257                             "halo exchange' feature will not be enabled as nonbonded interactions "
 258                             "are not offloaded.");
 259             devFlags.enableGpuHaloExchange = false;
 260         }
 261     }
 262
 263     if (devFlags.enableGpuPmePPComm)
 264     {
 265         if (pmeRunMode == PmeRunMode::GPU)
 266         {
 267             GMX_LOG(mdlog.warning)
 268                     .asParagraph()
 269                     .appendTextFormatted(
 270                             "This run uses the 'GPU PME-PP communications' feature, enabled "
 271                             "by the GMX_GPU_PME_PP_COMMS environment variable.");
 272         }
 273         else
 274         {
 275             std::string clarification;
 276             if (pmeRunMode == PmeRunMode::Mixed)
 277             {
 278                 clarification =
 279                         "PME FFT and gather are not offloaded to the GPU (PME is running in mixed "
 280                         "mode).";
 281             }
 282             else
 283             {
 284                 clarification = "PME is not offloaded to the GPU.";
 285             }
 286             GMX_LOG(mdlog.warning)
 287                     .asParagraph()
 288                     .appendText(
 289                             "GMX_GPU_PME_PP_COMMS environment variable detected, but the "
 290                             "'GPU PME-PP communications' feature was not enabled as "
 291                             + clarification);
 292             devFlags.enableGpuPmePPComm = false;
 293         }
 294     }
 295
 296     return devFlags;
 297 }
 298
 299 /*! \brief Barrier for safe simultaneous thread access to mdrunner data
 300  *
 301  * Used to ensure that the master thread does not modify mdrunner during copy
 302  * on the spawned threads. */
 303 static void threadMpiMdrunnerAccessBarrier()
 304 {
 305 #if GMX_THREAD_MPI
 306     MPI_Barrier(MPI_COMM_WORLD);
 307 #endif
 308 }
 309
 310 Mdrunner Mdrunner::cloneOnSpawnedThread() const
 311 {
 312     auto newRunner = Mdrunner(std::make_unique<MDModules>());
 313
 314     // All runners in the same process share a restraint manager resource because it is
 315     // part of the interface to the client code, which is associated only with the
 316     // original thread. Handles to the same resources can be obtained by copy.
 317     {
 318         newRunner.restraintManager_ = std::make_unique<RestraintManager>(*restraintManager_);
 319     }
 320
 321     // Copy members of master runner.
 322     // \todo Replace with builder when Simulation context and/or runner phases are better defined.
 323     // Ref https://gitlab.com/gromacs/gromacs/-/issues/2587 and https://gitlab.com/gromacs/gromacs/-/issues/2375
 324     newRunner.hw_opt    = hw_opt;
 325     newRunner.filenames = filenames;
 326
 327     newRunner.oenv            = oenv;
 328     newRunner.mdrunOptions    = mdrunOptions;
 329     newRunner.domdecOptions   = domdecOptions;
 330     newRunner.nbpu_opt        = nbpu_opt;
 331     newRunner.pme_opt         = pme_opt;
 332     newRunner.pme_fft_opt     = pme_fft_opt;
 333     newRunner.bonded_opt      = bonded_opt;
 334     newRunner.update_opt      = update_opt;
 335     newRunner.nstlist_cmdline = nstlist_cmdline;
 336     newRunner.replExParams    = replExParams;
 337     newRunner.pforce          = pforce;
 338     // Give the spawned thread the newly created valid communicator
 339     // for the simulation.
 340     newRunner.communicator        = MPI_COMM_WORLD;
 341     newRunner.ms                  = ms;
 342     newRunner.startingBehavior    = startingBehavior;
 343     newRunner.stopHandlerBuilder_ = std::make_unique<StopHandlerBuilder>(*stopHandlerBuilder_);
 344
 345     threadMpiMdrunnerAccessBarrier();
 346
 347     return newRunner;
 348 }
 349
 350 /*! \brief The callback used for running on spawned threads.
 351  *
 352  * Obtains the pointer to the master mdrunner object from the one
 353  * argument permitted to the thread-launch API call, copies it to make
 354  * a new runner for this thread, reinitializes necessary data, and
 355  * proceeds to the simulation. */
 356 static void mdrunner_start_fn(const void* arg)
 357 {
 358     try
 359     {
 360         auto masterMdrunner = reinterpret_cast<const gmx::Mdrunner*>(arg);
 361         /* copy the arg list to make sure that it's thread-local. This
 362            doesn't copy pointed-to items, of course; fnm, cr and fplog
 363            are reset in the call below, all others should be const. */
 364         gmx::Mdrunner mdrunner = masterMdrunner->cloneOnSpawnedThread();
 365         mdrunner.mdrunner();
 366     }
 367     GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
 368 }
 369
 370
 371 void Mdrunner::spawnThreads(int numThreadsToLaunch)
 372 {
 373 #if GMX_THREAD_MPI
 374     /* now spawn new threads that start mdrunner_start_fn(), while
 375        the main thread returns. Thread affinity is handled later. */
 376     if (tMPI_Init_fn(TRUE, numThreadsToLaunch, TMPI_AFFINITY_NONE, mdrunner_start_fn,
 377                      static_cast<const void*>(this))
 378         != TMPI_SUCCESS)
 379     {
 380         GMX_THROW(gmx::InternalError("Failed to spawn thread-MPI threads"));
 381     }
 382
 383     // Give the master thread the newly created valid communicator for
 384     // the simulation.
 385     communicator = MPI_COMM_WORLD;
 386     threadMpiMdrunnerAccessBarrier();
 387 #else
 388     GMX_UNUSED_VALUE(numThreadsToLaunch);
 389     GMX_UNUSED_VALUE(mdrunner_start_fn);
 390 #endif
 391 }
 392
 393 } // namespace gmx
 394
 395 /*! \brief Initialize variables for Verlet scheme simulation */
 396 static void prepare_verlet_scheme(FILE*               fplog,
 397                                   t_commrec*          cr,
 398                                   t_inputrec*         ir,
 399                                   int                 nstlist_cmdline,
 400                                   const gmx_mtop_t*   mtop,
 401                                   const matrix        box,
 402                                   bool                makeGpuPairList,
 403                                   const gmx::CpuInfo& cpuinfo)
 404 {
 405     // We checked the cut-offs in grompp, but double-check here.
 406     // We have PME+LJcutoff kernels for rcoulomb>rvdw.
 407     if (EEL_PME_EWALD(ir->coulombtype) && ir->vdwtype == eelCUT)
 408     {
 409         GMX_RELEASE_ASSERT(ir->rcoulomb >= ir->rvdw,
 410                            "With Verlet lists and PME we should have rcoulomb>=rvdw");
 411     }
 412     else
 413     {
 414         GMX_RELEASE_ASSERT(ir->rcoulomb == ir->rvdw,
 415                            "With Verlet lists and no PME rcoulomb and rvdw should be identical");
 416     }
 417     /* For NVE simulations, we will retain the initial list buffer */
 418     if (EI_DYNAMICS(ir->eI) && ir->verletbuf_tol > 0 && !(EI_MD(ir->eI) && ir->etc == etcNO))
 419     {
 420         /* Update the Verlet buffer size for the current run setup */
 421
 422         /* Here we assume SIMD-enabled kernels are being used. But as currently
 423          * calc_verlet_buffer_size gives the same results for 4x8 and 4x4
 424          * and 4x2 gives a larger buffer than 4x4, this is ok.
 425          */
 426         ListSetupType listType =
 427                 (makeGpuPairList ? ListSetupType::Gpu : ListSetupType::CpuSimdWhenSupported);
 428         VerletbufListSetup listSetup = verletbufGetSafeListSetup(listType);
 429
 430         const real rlist_new =
 431                 calcVerletBufferSize(*mtop, det(box), *ir, ir->nstlist, ir->nstlist - 1, -1, listSetup);
 432
 433         if (rlist_new != ir->rlist)
 434         {
 435             if (fplog != nullptr)
 436             {
 437                 fprintf(fplog,
 438                         "\nChanging rlist from %g to %g for non-bonded %dx%d atom kernels\n\n",
 439                         ir->rlist, rlist_new, listSetup.cluster_size_i, listSetup.cluster_size_j);
 440             }
 441             ir->rlist = rlist_new;
 442         }
 443     }
 444
 445     if (nstlist_cmdline > 0 && (!EI_DYNAMICS(ir->eI) || ir->verletbuf_tol <= 0))
 446     {
 447         gmx_fatal(FARGS, "Can not set nstlist without %s",
 448                   !EI_DYNAMICS(ir->eI) ? "dynamics" : "verlet-buffer-tolerance");
 449     }
 450
 451     if (EI_DYNAMICS(ir->eI))
 452     {
 453         /* Set or try nstlist values */
 454         increaseNstlist(fplog, cr, ir, nstlist_cmdline, mtop, box, makeGpuPairList, cpuinfo);
 455     }
 456 }
 457
 458 /*! \brief Override the nslist value in inputrec
 459  *
 460  * with value passed on the command line (if any)
 461  */
 462 static void override_nsteps_cmdline(const gmx::MDLogger& mdlog, int64_t nsteps_cmdline, t_inputrec* ir)
 463 {
 464     assert(ir);
 465
 466     /* override with anything else than the default -2 */
 467     if (nsteps_cmdline > -2)
 468     {
 469         char sbuf_steps[STEPSTRSIZE];
 470         char sbuf_msg[STRLEN];
 471
 472         ir->nsteps = nsteps_cmdline;
 473         if (EI_DYNAMICS(ir->eI) && nsteps_cmdline != -1)
 474         {
 475             sprintf(sbuf_msg,
 476                     "Overriding nsteps with value passed on the command line: %s steps, %.3g ps",
 477                     gmx_step_str(nsteps_cmdline, sbuf_steps), fabs(nsteps_cmdline * ir->delta_t));
 478         }
 479         else
 480         {
 481             sprintf(sbuf_msg, "Overriding nsteps with value passed on the command line: %s steps",
 482                     gmx_step_str(nsteps_cmdline, sbuf_steps));
 483         }
 484
 485         GMX_LOG(mdlog.warning).asParagraph().appendText(sbuf_msg);
 486     }
 487     else if (nsteps_cmdline < -2)
 488     {
 489         gmx_fatal(FARGS, "Invalid nsteps value passed on the command line: %" PRId64, nsteps_cmdline);
 490     }
 491     /* Do nothing if nsteps_cmdline == -2 */
 492 }
 493
 494 namespace gmx
 495 {
 496
 497 /*! \brief Return whether GPU acceleration of nonbondeds is supported with the given settings.
 498  *
 499  * If not, and if a warning may be issued, logs a warning about
 500  * falling back to CPU code. With thread-MPI, only the first
 501  * call to this function should have \c issueWarning true. */
 502 static bool gpuAccelerationOfNonbondedIsUseful(const MDLogger& mdlog, const t_inputrec& ir, bool issueWarning)
 503 {
 504     bool        gpuIsUseful = true;
 505     std::string warning;
 506
 507     if (ir.opts.ngener - ir.nwall > 1)
 508     {
 509         /* The GPU code does not support more than one energy group.
 510          * If the user requested GPUs explicitly, a fatal error is given later.
 511          */
 512         gpuIsUseful = false;
 513         warning =
 514                 "Multiple energy groups is not implemented for GPUs, falling back to the CPU. "
 515                 "For better performance, run on the GPU without energy groups and then do "
 516                 "gmx mdrun -rerun option on the trajectory with an energy group .tpr file.";
 517     }
 518
 519     if (EI_TPI(ir.eI))
 520     {
 521         gpuIsUseful = false;
 522         warning     = "TPI is not implemented for GPUs.";
 523     }
 524
 525     if (!gpuIsUseful && issueWarning)
 526     {
 527         GMX_LOG(mdlog.warning).asParagraph().appendText(warning);
 528     }
 529
 530     return gpuIsUseful;
 531 }
 532
 533 //! Initializes the logger for mdrun.
 534 static gmx::LoggerOwner buildLogger(FILE* fplog, const bool isSimulationMasterRank)
 535 {
 536     gmx::LoggerBuilder builder;
 537     if (fplog != nullptr)
 538     {
 539         builder.addTargetFile(gmx::MDLogger::LogLevel::Info, fplog);
 540     }
 541     if (isSimulationMasterRank)
 542     {
 543         builder.addTargetStream(gmx::MDLogger::LogLevel::Warning, &gmx::TextOutputFile::standardError());
 544     }
 545     return builder.build();
 546 }
 547
 548 //! Make a TaskTarget from an mdrun argument string.
 549 static TaskTarget findTaskTarget(const char* optionString)
 550 {
 551     TaskTarget returnValue = TaskTarget::Auto;
 552
 553     if (strncmp(optionString, "auto", 3) == 0)
 554     {
 555         returnValue = TaskTarget::Auto;
 556     }
 557     else if (strncmp(optionString, "cpu", 3) == 0)
 558     {
 559         returnValue = TaskTarget::Cpu;
 560     }
 561     else if (strncmp(optionString, "gpu", 3) == 0)
 562     {
 563         returnValue = TaskTarget::Gpu;
 564     }
 565     else
 566     {
 567         GMX_ASSERT(false, "Option string should have been checked for sanity already");
 568     }
 569
 570     return returnValue;
 571 }
 572
 573 //! Finish run, aggregate data to print performance info.
 574 static void finish_run(FILE*                     fplog,
 575                        const gmx::MDLogger&      mdlog,
 576                        const t_commrec*          cr,
 577                        const t_inputrec*         inputrec,
 578                        t_nrnb                    nrnb[],
 579                        gmx_wallcycle_t           wcycle,
 580                        gmx_walltime_accounting_t walltime_accounting,
 581                        nonbonded_verlet_t*       nbv,
 582                        const gmx_pme_t*          pme,
 583                        gmx_bool                  bWriteStat)
 584 {
 585     double delta_t = 0;
 586     double nbfs = 0, mflop = 0;
 587     double elapsed_time, elapsed_time_over_all_ranks, elapsed_time_over_all_threads,
 588             elapsed_time_over_all_threads_over_all_ranks;
 589     /* Control whether it is valid to print a report. Only the
 590        simulation master may print, but it should not do so if the run
 591        terminated e.g. before a scheduled reset step. This is
 592        complicated by the fact that PME ranks are unaware of the
 593        reason why they were sent a pmerecvqxFINISH. To avoid
 594        communication deadlocks, we always do the communication for the
 595        report, even if we've decided not to write the report, because
 596        how long it takes to finish the run is not important when we've
 597        decided not to report on the simulation performance.
 598
 599        Further, we only report performance for dynamical integrators,
 600        because those are the only ones for which we plan to
 601        consider doing any optimizations. */
 602     bool printReport = EI_DYNAMICS(inputrec->eI) && SIMMASTER(cr);
 603
 604     if (printReport && !walltime_accounting_get_valid_finish(walltime_accounting))
 605     {
 606         GMX_LOG(mdlog.warning)
 607                 .asParagraph()
 608                 .appendText("Simulation ended prematurely, no performance report will be written.");
 609         printReport = false;
 610     }
 611
 612     t_nrnb*                 nrnb_tot;
 613     std::unique_ptr<t_nrnb> nrnbTotalStorage;
 614     if (cr->nnodes > 1)
 615     {
 616         nrnbTotalStorage = std::make_unique<t_nrnb>();
 617         nrnb_tot         = nrnbTotalStorage.get();
 618 #if GMX_MPI
 619         MPI_Allreduce(nrnb->n, nrnb_tot->n, eNRNB, MPI_DOUBLE, MPI_SUM, cr->mpi_comm_mysim);
 620 #endif
 621     }
 622     else
 623     {
 624         nrnb_tot = nrnb;
 625     }
 626
 627     elapsed_time = walltime_accounting_get_time_since_reset(walltime_accounting);
 628     elapsed_time_over_all_threads =
 629             walltime_accounting_get_time_since_reset_over_all_threads(walltime_accounting);
 630     if (cr->nnodes > 1)
 631     {
 632 #if GMX_MPI
 633         /* reduce elapsed_time over all MPI ranks in the current simulation */
 634         MPI_Allreduce(&elapsed_time, &elapsed_time_over_all_ranks, 1, MPI_DOUBLE, MPI_SUM,
 635                       cr->mpi_comm_mysim);
 636         elapsed_time_over_all_ranks /= cr->nnodes;
 637         /* Reduce elapsed_time_over_all_threads over all MPI ranks in the
 638          * current simulation. */
 639         MPI_Allreduce(&elapsed_time_over_all_threads, &elapsed_time_over_all_threads_over_all_ranks,
 640                       1, MPI_DOUBLE, MPI_SUM, cr->mpi_comm_mysim);
 641 #endif
 642     }
 643     else
 644     {
 645         elapsed_time_over_all_ranks                  = elapsed_time;
 646         elapsed_time_over_all_threads_over_all_ranks = elapsed_time_over_all_threads;
 647     }
 648
 649     if (printReport)
 650     {
 651         print_flop(fplog, nrnb_tot, &nbfs, &mflop);
 652     }
 653
 654     if (thisRankHasDuty(cr, DUTY_PP) && DOMAINDECOMP(cr))
 655     {
 656         print_dd_statistics(cr, inputrec, fplog);
 657     }
 658
 659     /* TODO Move the responsibility for any scaling by thread counts
 660      * to the code that handled the thread region, so that there's a
 661      * mechanism to keep cycle counting working during the transition
 662      * to task parallelism. */
 663     int nthreads_pp  = gmx_omp_nthreads_get(emntNonbonded);
 664     int nthreads_pme = gmx_omp_nthreads_get(emntPME);
 665     wallcycle_scale_by_num_threads(wcycle, thisRankHasDuty(cr, DUTY_PME) && !thisRankHasDuty(cr, DUTY_PP),
 666                                    nthreads_pp, nthreads_pme);
 667     auto cycle_sum(wallcycle_sum(cr, wcycle));
 668
 669     if (printReport)
 670     {
 671         auto nbnxn_gpu_timings =
 672                 (nbv != nullptr && nbv->useGpu()) ? Nbnxm::gpu_get_timings(nbv->gpu_nbv) : nullptr;
 673         gmx_wallclock_gpu_pme_t pme_gpu_timings = {};
 674
 675         if (pme_gpu_task_enabled(pme))
 676         {
 677             pme_gpu_get_timings(pme, &pme_gpu_timings);
 678         }
 679         wallcycle_print(fplog, mdlog, cr->nnodes, cr->npmenodes, nthreads_pp, nthreads_pme,
 680                         elapsed_time_over_all_ranks, wcycle, cycle_sum, nbnxn_gpu_timings,
 681                         &pme_gpu_timings);
 682
 683         if (EI_DYNAMICS(inputrec->eI))
 684         {
 685             delta_t = inputrec->delta_t;
 686         }
 687
 688         if (fplog)
 689         {
 690             print_perf(fplog, elapsed_time_over_all_threads_over_all_ranks, elapsed_time_over_all_ranks,
 691                        walltime_accounting_get_nsteps_done_since_reset(walltime_accounting),
 692                        delta_t, nbfs, mflop);
 693         }
 694         if (bWriteStat)
 695         {
 696             print_perf(stderr, elapsed_time_over_all_threads_over_all_ranks, elapsed_time_over_all_ranks,
 697                        walltime_accounting_get_nsteps_done_since_reset(walltime_accounting),
 698                        delta_t, nbfs, mflop);
 699         }
 700     }
 701 }
 702
 703 int Mdrunner::mdrunner()
 704 {
 705     matrix                    box;
 706     t_forcerec*               fr               = nullptr;
 707     t_fcdata*                 fcd              = nullptr;
 708     real                      ewaldcoeff_q     = 0;
 709     real                      ewaldcoeff_lj    = 0;
 710     int                       nChargePerturbed = -1, nTypePerturbed = 0;
 711     gmx_wallcycle_t           wcycle;
 712     gmx_walltime_accounting_t walltime_accounting = nullptr;
 713     gmx_membed_t*             membed              = nullptr;
 714     gmx_hw_info_t*            hwinfo              = nullptr;
 715
 716     /* CAUTION: threads may be started later on in this function, so
 717        cr doesn't reflect the final parallel state right now */
 718     gmx_mtop_t mtop;
 719
 720     /* TODO: inputrec should tell us whether we use an algorithm, not a file option */
 721     const bool doEssentialDynamics = opt2bSet("-ei", filenames.size(), filenames.data());
 722     const bool doMembed            = opt2bSet("-membed", filenames.size(), filenames.data());
 723     const bool doRerun             = mdrunOptions.rerun;
 724
 725     // Handle task-assignment related user options.
 726     EmulateGpuNonbonded emulateGpuNonbonded =
 727             (getenv("GMX_EMULATE_GPU") != nullptr ? EmulateGpuNonbonded::Yes : EmulateGpuNonbonded::No);
 728
 729     std::vector<int> userGpuTaskAssignment;
 730     try
 731     {
 732         userGpuTaskAssignment = parseUserTaskAssignmentString(hw_opt.userGpuTaskAssignment);
 733     }
 734     GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
 735     auto nonbondedTarget = findTaskTarget(nbpu_opt);
 736     auto pmeTarget       = findTaskTarget(pme_opt);
 737     auto pmeFftTarget    = findTaskTarget(pme_fft_opt);
 738     auto bondedTarget    = findTaskTarget(bonded_opt);
 739     auto updateTarget    = findTaskTarget(update_opt);
 740
 741     FILE* fplog = nullptr;
 742     // If we are appending, we don't write log output because we need
 743     // to check that the old log file matches what the checkpoint file
 744     // expects. Otherwise, we should start to write log output now if
 745     // there is a file ready for it.
 746     if (logFileHandle != nullptr && startingBehavior != StartingBehavior::RestartWithAppending)
 747     {
 748         fplog = gmx_fio_getfp(logFileHandle);
 749     }
 750     const bool       isSimulationMasterRank = findIsSimulationMasterRank(ms, communicator);
 751     gmx::LoggerOwner logOwner(buildLogger(fplog, isSimulationMasterRank));
 752     gmx::MDLogger    mdlog(logOwner.logger());
 753
 754     // TODO The thread-MPI master rank makes a working
 755     // PhysicalNodeCommunicator here, but it gets rebuilt by all ranks
 756     // after the threads have been launched. This works because no use
 757     // is made of that communicator until after the execution paths
 758     // have rejoined. But it is likely that we can improve the way
 759     // this is expressed, e.g. by expressly running detection only the
 760     // master rank for thread-MPI, rather than relying on the mutex
 761     // and reference count.
 762     PhysicalNodeCommunicator physicalNodeComm(communicator, gmx_physicalnode_id_hash());
 763     hwinfo = gmx_detect_hardware(mdlog, physicalNodeComm);
 764
 765     gmx_print_detected_hardware(fplog, isSimulationMasterRank && isMasterSim(ms), mdlog, hwinfo);
 766
 767     std::vector<int> gpuIdsToUse = makeGpuIdsToUse(hwinfo->gpu_info, hw_opt.gpuIdsAvailable);
 768
 769     // Print citation requests after all software/hardware printing
 770     pleaseCiteGromacs(fplog);
 771
 772     // TODO Replace this by unique_ptr once t_inputrec is C++
 773     t_inputrec               inputrecInstance;
 774     t_inputrec*              inputrec = nullptr;
 775     std::unique_ptr<t_state> globalState;
 776
 777     auto partialDeserializedTpr = std::make_unique<PartialDeserializedTprFile>();
 778
 779     if (isSimulationMasterRank)
 780     {
 781         /* Only the master rank has the global state */
 782         globalState = std::make_unique<t_state>();
 783
 784         /* Read (nearly) all data required for the simulation
 785          * and keep the partly serialized tpr contents to send to other ranks later
 786          */
 787         *partialDeserializedTpr = read_tpx_state(ftp2fn(efTPR, filenames.size(), filenames.data()),
 788                                                  &inputrecInstance, globalState.get(), &mtop);
 789         inputrec                = &inputrecInstance;
 790     }
 791
 792     /* Check and update the hardware options for internal consistency */
 793     checkAndUpdateHardwareOptions(mdlog, &hw_opt, isSimulationMasterRank, domdecOptions.numPmeRanks,
 794                                   inputrec);
 795
 796     if (GMX_THREAD_MPI && isSimulationMasterRank)
 797     {
 798         bool useGpuForNonbonded = false;
 799         bool useGpuForPme       = false;
 800         try
 801         {
 802             GMX_RELEASE_ASSERT(inputrec != nullptr, "Keep the compiler happy");
 803
 804             // If the user specified the number of ranks, then we must
 805             // respect that, but in default mode, we need to allow for
 806             // the number of GPUs to choose the number of ranks.
 807             auto canUseGpuForNonbonded = buildSupportsNonbondedOnGpu(nullptr);
 808             useGpuForNonbonded         = decideWhetherToUseGpusForNonbondedWithThreadMpi(
 809                     nonbondedTarget, gpuIdsToUse, userGpuTaskAssignment, emulateGpuNonbonded,
 810                     canUseGpuForNonbonded,
 811                     gpuAccelerationOfNonbondedIsUseful(mdlog, *inputrec, GMX_THREAD_MPI),
 812                     hw_opt.nthreads_tmpi);
 813             useGpuForPme = decideWhetherToUseGpusForPmeWithThreadMpi(
 814                     useGpuForNonbonded, pmeTarget, gpuIdsToUse, userGpuTaskAssignment, *hwinfo,
 815                     *inputrec, mtop, hw_opt.nthreads_tmpi, domdecOptions.numPmeRanks);
 816         }
 817         GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
 818
 819         /* Determine how many thread-MPI ranks to start.
 820          *
 821          * TODO Over-writing the user-supplied value here does
 822          * prevent any possible subsequent checks from working
 823          * correctly. */
 824         hw_opt.nthreads_tmpi = get_nthreads_mpi(hwinfo, &hw_opt, gpuIdsToUse, useGpuForNonbonded,
 825                                                 useGpuForPme, inputrec, &mtop, mdlog, doMembed);
 826
 827         // Now start the threads for thread MPI.
 828         spawnThreads(hw_opt.nthreads_tmpi);
 829         // The spawned threads enter mdrunner() and execution of
 830         // master and spawned threads joins at the end of this block.
 831         physicalNodeComm = PhysicalNodeCommunicator(communicator, gmx_physicalnode_id_hash());
 832     }
 833
 834     GMX_RELEASE_ASSERT(communicator == MPI_COMM_WORLD, "Must have valid world communicator");
 835     CommrecHandle crHandle = init_commrec(communicator, ms);
 836     t_commrec*    cr       = crHandle.get();
 837     GMX_RELEASE_ASSERT(cr != nullptr, "Must have valid commrec");
 838
 839     if (PAR(cr))
 840     {
 841         /* now broadcast everything to the non-master nodes/threads: */
 842         if (!isSimulationMasterRank)
 843         {
 844             inputrec = &inputrecInstance;
 845         }
 846         init_parallel(cr, inputrec, &mtop, partialDeserializedTpr.get());
 847     }
 848     GMX_RELEASE_ASSERT(inputrec != nullptr, "All ranks should have a valid inputrec now");
 849     partialDeserializedTpr.reset(nullptr);
 850
 851     // Now the number of ranks is known to all ranks, and each knows
 852     // the inputrec read by the master rank. The ranks can now all run
 853     // the task-deciding functions and will agree on the result
 854     // without needing to communicate.
 855     const bool useDomainDecomposition = (PAR(cr) && !(EI_TPI(inputrec->eI) || inputrec->eI == eiNM));
 856
 857     // Note that these variables describe only their own node.
 858     //
 859     // Note that when bonded interactions run on a GPU they always run
 860     // alongside a nonbonded task, so do not influence task assignment
 861     // even though they affect the force calculation workload.
 862     bool useGpuForNonbonded = false;
 863     bool useGpuForPme       = false;
 864     bool useGpuForBonded    = false;
 865     bool useGpuForUpdate    = false;
 866     bool gpusWereDetected   = hwinfo->ngpu_compatible_tot > 0;
 867     try
 868     {
 869         // It's possible that there are different numbers of GPUs on
 870         // different nodes, which is the user's responsibility to
 871         // handle. If unsuitable, we will notice that during task
 872         // assignment.
 873         auto canUseGpuForNonbonded = buildSupportsNonbondedOnGpu(nullptr);
 874         useGpuForNonbonded         = decideWhetherToUseGpusForNonbonded(
 875                 nonbondedTarget, userGpuTaskAssignment, emulateGpuNonbonded, canUseGpuForNonbonded,
 876                 gpuAccelerationOfNonbondedIsUseful(mdlog, *inputrec, !GMX_THREAD_MPI), gpusWereDetected);
 877         useGpuForPme = decideWhetherToUseGpusForPme(
 878                 useGpuForNonbonded, pmeTarget, userGpuTaskAssignment, *hwinfo, *inputrec, mtop,
 879                 cr->nnodes, domdecOptions.numPmeRanks, gpusWereDetected);
 880         auto canUseGpuForBonded = buildSupportsGpuBondeds(nullptr)
 881                                   && inputSupportsGpuBondeds(*inputrec, mtop, nullptr);
 882         useGpuForBonded = decideWhetherToUseGpusForBonded(
 883                 useGpuForNonbonded, useGpuForPme, bondedTarget, canUseGpuForBonded,
 884                 EVDW_PME(inputrec->vdwtype), EEL_PME_EWALD(inputrec->coulombtype),
 885                 domdecOptions.numPmeRanks, gpusWereDetected);
 886     }
 887     GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
 888
 889     const PmeRunMode pmeRunMode = determinePmeRunMode(useGpuForPme, pmeFftTarget, *inputrec);
 890
 891     // Initialize development feature flags that enabled by environment variable
 892     // and report those features that are enabled.
 893     const DevelopmentFeatureFlags devFlags =
 894             manageDevelopmentFeatures(mdlog, useGpuForNonbonded, pmeRunMode);
 895
 896     const bool useModularSimulator = checkUseModularSimulator(
 897             false, inputrec, doRerun, mtop, ms, replExParams, nullptr, doEssentialDynamics, doMembed);
 898
 899     // Build restraints.
 900     // TODO: hide restraint implementation details from Mdrunner.
 901     // There is nothing unique about restraints at this point as far as the
 902     // Mdrunner is concerned. The Mdrunner should just be getting a sequence of
 903     // factory functions from the SimulationContext on which to call mdModules_->add().
 904     // TODO: capture all restraints into a single RestraintModule, passed to the runner builder.
 905     for (auto&& restraint : restraintManager_->getRestraints())
 906     {
 907         auto module = RestraintMDModule::create(restraint, restraint->sites());
 908         mdModules_->add(std::move(module));
 909     }
 910
 911     // TODO: Error handling
 912     mdModules_->assignOptionsToModules(*inputrec->params, nullptr);
 913     // now that the MdModules know their options, they know which callbacks to sign up to
 914     mdModules_->subscribeToSimulationSetupNotifications();
 915     const auto& mdModulesNotifier = mdModules_->notifier().simulationSetupNotifications_;
 916
 917     if (inputrec->internalParameters != nullptr)
 918     {
 919         mdModulesNotifier.notify(*inputrec->internalParameters);
 920     }
 921
 922     if (fplog != nullptr)
 923     {
 924         pr_inputrec(fplog, 0, "Input Parameters", inputrec, FALSE);
 925         fprintf(fplog, "\n");
 926     }
 927
 928     if (SIMMASTER(cr))
 929     {
 930         /* In rerun, set velocities to zero if present */
 931         if (doRerun && ((globalState->flags & (1 << estV)) != 0))
 932         {
 933             // rerun does not use velocities
 934             GMX_LOG(mdlog.info)
 935                     .asParagraph()
 936                     .appendText(
 937                             "Rerun trajectory contains velocities. Rerun does only evaluate "
 938                             "potential energy and forces. The velocities will be ignored.");
 939             for (int i = 0; i < globalState->natoms; i++)
 940             {
 941                 clear_rvec(globalState->v[i]);
 942             }
 943             globalState->flags &= ~(1 << estV);
 944         }
 945
 946         /* now make sure the state is initialized and propagated */
 947         set_state_entries(globalState.get(), inputrec, useModularSimulator);
 948     }
 949
 950     /* NM and TPI parallelize over force/energy calculations, not atoms,
 951      * so we need to initialize and broadcast the global state.
 952      */
 953     if (inputrec->eI == eiNM || inputrec->eI == eiTPI)
 954     {
 955         if (!MASTER(cr))
 956         {
 957             globalState = std::make_unique<t_state>();
 958         }
 959         broadcastStateWithoutDynamics(cr, globalState.get());
 960     }
 961
 962     /* A parallel command line option consistency check that we can
 963        only do after any threads have started. */
 964     if (!PAR(cr)
 965         && (domdecOptions.numCells[XX] > 1 || domdecOptions.numCells[YY] > 1
 966             || domdecOptions.numCells[ZZ] > 1 || domdecOptions.numPmeRanks > 0))
 967     {
 968         gmx_fatal(FARGS,
 969                   "The -dd or -npme option request a parallel simulation, "
 970 #if !GMX_MPI
 971                   "but %s was compiled without threads or MPI enabled",
 972                   output_env_get_program_display_name(oenv));
 973 #elif GMX_THREAD_MPI
 974                   "but the number of MPI-threads (option -ntmpi) is not set or is 1");
 975 #else
 976                   "but %s was not started through mpirun/mpiexec or only one rank was requested "
 977                   "through mpirun/mpiexec",
 978                   output_env_get_program_display_name(oenv));
 979 #endif
 980     }
 981
 982     if (doRerun && (EI_ENERGY_MINIMIZATION(inputrec->eI) || eiNM == inputrec->eI))
 983     {
 984         gmx_fatal(FARGS,
 985                   "The .mdp file specified an energy mininization or normal mode algorithm, and "
 986                   "these are not compatible with mdrun -rerun");
 987     }
 988
 989     if (!(EEL_PME(inputrec->coulombtype) || EVDW_PME(inputrec->vdwtype)))
 990     {
 991         if (domdecOptions.numPmeRanks > 0)
 992         {
 993             gmx_fatal_collective(FARGS, cr->mpi_comm_mysim, MASTER(cr),
 994                                  "PME-only ranks are requested, but the system does not use PME "
 995                                  "for electrostatics or LJ");
 996         }
 997
 998         domdecOptions.numPmeRanks = 0;
 999     }
1000
1001     if (useGpuForNonbonded && domdecOptions.numPmeRanks < 0)
1002     {
1003         /* With NB GPUs we don't automatically use PME-only CPU ranks. PME ranks can
1004          * improve performance with many threads per GPU, since our OpenMP
1005          * scaling is bad, but it's difficult to automate the setup.
1006          */
1007         domdecOptions.numPmeRanks = 0;
1008     }
1009     if (useGpuForPme)
1010     {
1011         if (domdecOptions.numPmeRanks < 0)
1012         {
1013             domdecOptions.numPmeRanks = 0;
1014             // TODO possibly print a note that one can opt-in for a separate PME GPU rank?
1015         }
1016         else
1017         {
1018             GMX_RELEASE_ASSERT(domdecOptions.numPmeRanks <= 1,
1019                                "PME GPU decomposition is not supported");
1020         }
1021     }
1022
1023 #if GMX_FAHCORE
1024     if (MASTER(cr))
1025     {
1026         fcRegisterSteps(inputrec->nsteps, inputrec->init_step);
1027     }
1028 #endif
1029
1030     /* NMR restraints must be initialized before load_checkpoint,
1031      * since with time averaging the history is added to t_state.
1032      * For proper consistency check we therefore need to extend
1033      * t_state here.
1034      * So the PME-only nodes (if present) will also initialize
1035      * the distance restraints.
1036      */
1037     snew(fcd, 1);
1038
1039     /* This needs to be called before read_checkpoint to extend the state */
1040     init_disres(fplog, &mtop, inputrec, cr, ms, fcd, globalState.get(), replExParams.exchangeInterval > 0);
1041
1042     init_orires(fplog, &mtop, inputrec, cr, ms, globalState.get(), &(fcd->orires));
1043
1044     auto deform = prepareBoxDeformation(globalState->box, cr, *inputrec);
1045
1046     ObservablesHistory observablesHistory = {};
1047
1048     if (startingBehavior != StartingBehavior::NewSimulation)
1049     {
1050         /* Check if checkpoint file exists before doing continuation.
1051          * This way we can use identical input options for the first and subsequent runs...
1052          */
1053         if (mdrunOptions.numStepsCommandline > -2)
1054         {
1055             /* Temporarily set the number of steps to unlimited to avoid
1056              * triggering the nsteps check in load_checkpoint().
1057              * This hack will go away soon when the -nsteps option is removed.
1058              */
1059             inputrec->nsteps = -1;
1060         }
1061
1062         load_checkpoint(opt2fn_master("-cpi", filenames.size(), filenames.data(), cr),
1063                         logFileHandle, cr, domdecOptions.numCells, inputrec, globalState.get(),
1064                         &observablesHistory, mdrunOptions.reproducible, mdModules_->notifier());
1065
1066         if (startingBehavior == StartingBehavior::RestartWithAppending && logFileHandle)
1067         {
1068             // Now we can start normal logging to the truncated log file.
1069             fplog = gmx_fio_getfp(logFileHandle);
1070             prepareLogAppending(fplog);
1071             logOwner = buildLogger(fplog, MASTER(cr));
1072             mdlog    = logOwner.logger();
1073         }
1074     }
1075
1076     if (mdrunOptions.numStepsCommandline > -2)
1077     {
1078         GMX_LOG(mdlog.info)
1079                 .asParagraph()
1080                 .appendText(
1081                         "The -nsteps functionality is deprecated, and may be removed in a future "
1082                         "version. "
1083                         "Consider using gmx convert-tpr -nsteps or changing the appropriate .mdp "
1084                         "file field.");
1085     }
1086     /* override nsteps with value set on the commandline */
1087     override_nsteps_cmdline(mdlog, mdrunOptions.numStepsCommandline, inputrec);
1088
1089     if (SIMMASTER(cr))
1090     {
1091         copy_mat(globalState->box, box);
1092     }
1093
1094     if (PAR(cr))
1095     {
1096         gmx_bcast(sizeof(box), box, cr);
1097     }
1098
1099     if (inputrec->cutoff_scheme != ecutsVERLET)
1100     {
1101         gmx_fatal(FARGS,
1102                   "This group-scheme .tpr file can no longer be run by mdrun. Please update to the "
1103                   "Verlet scheme, or use an earlier version of GROMACS if necessary.");
1104     }
1105     /* Update rlist and nstlist. */
1106     prepare_verlet_scheme(fplog, cr, inputrec, nstlist_cmdline, &mtop, box,
1107                           useGpuForNonbonded || (emulateGpuNonbonded == EmulateGpuNonbonded::Yes),
1108                           *hwinfo->cpuInfo);
1109
1110     const bool prefer1DAnd1PulseDD = (devFlags.enableGpuHaloExchange && useGpuForNonbonded);
1111     // This builder is necessary while we have multi-part construction
1112     // of DD. Before DD is constructed, we use the existence of
1113     // the builder object to indicate that further construction of DD
1114     // is needed.
1115     std::unique_ptr<DomainDecompositionBuilder> ddBuilder;
1116     if (useDomainDecomposition)
1117     {
1118         ddBuilder = std::make_unique<DomainDecompositionBuilder>(
1119                 mdlog, cr, domdecOptions, mdrunOptions, prefer1DAnd1PulseDD, mtop, *inputrec, box,
1120                 positionsFromStatePointer(globalState.get()));
1121     }
1122     else
1123     {
1124         /* PME, if used, is done on all nodes with 1D decomposition */
1125         cr->npmenodes = 0;
1126         cr->duty      = (DUTY_PP | DUTY_PME);
1127
1128         if (inputrec->pbcType == PbcType::Screw)
1129         {
1130             gmx_fatal(FARGS, "pbc=screw is only implemented with domain decomposition");
1131         }
1132     }
1133
1134     // Produce the task assignment for this rank.
1135     GpuTaskAssignmentsBuilder gpuTaskAssignmentsBuilder;
1136     GpuTaskAssignments        gpuTaskAssignments = gpuTaskAssignmentsBuilder.build(
1137             gpuIdsToUse, userGpuTaskAssignment, *hwinfo, communicator, physicalNodeComm,
1138             nonbondedTarget, pmeTarget, bondedTarget, updateTarget, useGpuForNonbonded,
1139             useGpuForPme, thisRankHasDuty(cr, DUTY_PP),
1140             // TODO cr->duty & DUTY_PME should imply that a PME
1141             // algorithm is active, but currently does not.
1142             EEL_PME(inputrec->coulombtype) && thisRankHasDuty(cr, DUTY_PME));
1143
1144     // Get the device handles for the modules, nullptr when no task is assigned.
1145     int                deviceId   = -1;
1146     DeviceInformation* deviceInfo = gpuTaskAssignments.initDevice(&deviceId);
1147
1148     // timing enabling - TODO put this in gpu_utils (even though generally this is just option handling?)
1149     bool useTiming = true;
1150     if (GMX_GPU == GMX_GPU_CUDA)
1151     {
1152         /* WARNING: CUDA timings are incorrect with multiple streams.
1153          *          This is the main reason why they are disabled by default.
1154          */
1155         // TODO: Consider turning on by default when we can detect nr of streams.
1156         useTiming = (getenv("GMX_ENABLE_GPU_TIMING") != nullptr);
1157     }
1158     else if (GMX_GPU == GMX_GPU_OPENCL)
1159     {
1160         useTiming = (getenv("GMX_DISABLE_GPU_TIMING") == nullptr);
1161     }
1162
1163     // TODO Currently this is always built, yet DD partition code
1164     // checks if it is built before using it. Probably it should
1165     // become an MDModule that is made only when another module
1166     // requires it (e.g. pull, CompEl, density fitting), so that we
1167     // don't update the local atom sets unilaterally every step.
1168     LocalAtomSetManager atomSets;
1169     if (ddBuilder)
1170     {
1171         // TODO Pass the GPU streams to ddBuilder to use in buffer
1172         // transfers (e.g. halo exchange)
1173         cr->dd = ddBuilder->build(&atomSets);
1174         // The builder's job is done, so destruct it
1175         ddBuilder.reset(nullptr);
1176         // Note that local state still does not exist yet.
1177     }
1178
1179     // The GPU update is decided here because we need to know whether the constraints or
1180     // SETTLEs can span accross the domain borders (i.e. whether or not update groups are
1181     // defined). This is only known after DD is initialized, hence decision on using GPU
1182     // update is done so late.
1183     try
1184     {
1185         const bool useUpdateGroups = cr->dd ? ddUsesUpdateGroups(*cr->dd) : false;
1186
1187         useGpuForUpdate = decideWhetherToUseGpuForUpdate(
1188                 useDomainDecomposition, useUpdateGroups, pmeRunMode, domdecOptions.numPmeRanks > 0,
1189                 useGpuForNonbonded, updateTarget, gpusWereDetected, *inputrec, mtop,
1190                 doEssentialDynamics, gmx_mtop_ftype_count(mtop, F_ORIRES) > 0,
1191                 replExParams.exchangeInterval > 0, doRerun, devFlags, mdlog);
1192     }
1193     GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
1194
1195     const bool printHostName = (cr->nnodes > 1);
1196     gpuTaskAssignments.reportGpuUsage(mdlog, printHostName, useGpuForBonded, pmeRunMode, useGpuForUpdate);
1197
1198     std::unique_ptr<DeviceStreamManager> deviceStreamManager = nullptr;
1199
1200     if (deviceInfo != nullptr)
1201     {
1202         if (DOMAINDECOMP(cr) && thisRankHasDuty(cr, DUTY_PP))
1203         {
1204             dd_setup_dlb_resource_sharing(cr, deviceId);
1205         }
1206         deviceStreamManager = std::make_unique<DeviceStreamManager>(
1207                 *deviceInfo, useGpuForPme, useGpuForNonbonded, havePPDomainDecomposition(cr),
1208                 useGpuForUpdate, useTiming);
1209     }
1210
1211     // If the user chose a task assignment, give them some hints
1212     // where appropriate.
1213     if (!userGpuTaskAssignment.empty())
1214     {
1215         gpuTaskAssignments.logPerformanceHints(mdlog, ssize(gpuIdsToUse));
1216     }
1217
1218     if (PAR(cr))
1219     {
1220         /* After possible communicator splitting in make_dd_communicators.
1221          * we can set up the intra/inter node communication.
1222          */
1223         gmx_setup_nodecomm(fplog, cr);
1224     }
1225
1226 #if GMX_MPI
1227     if (isMultiSim(ms))
1228     {
1229         GMX_LOG(mdlog.warning)
1230                 .asParagraph()
1231                 .appendTextFormatted(
1232                         "This is simulation %d out of %d running as a composite GROMACS\n"
1233                         "multi-simulation job. Setup for this simulation:\n",
1234                         ms->sim, ms->nsim);
1235     }
1236     GMX_LOG(mdlog.warning)
1237             .appendTextFormatted("Using %d MPI %s\n", cr->nnodes,
1238 #    if GMX_THREAD_MPI
1239                                  cr->nnodes == 1 ? "thread" : "threads"
1240 #    else
1241                                  cr->nnodes == 1 ? "process" : "processes"
1242 #    endif
1243             );
1244     fflush(stderr);
1245 #endif
1246
1247     // If mdrun -pin auto honors any affinity setting that already
1248     // exists. If so, it is nice to provide feedback about whether
1249     // that existing affinity setting was from OpenMP or something
1250     // else, so we run this code both before and after we initialize
1251     // the OpenMP support.
1252     gmx_check_thread_affinity_set(mdlog, &hw_opt, hwinfo->nthreads_hw_avail, FALSE);
1253     /* Check and update the number of OpenMP threads requested */
1254     checkAndUpdateRequestedNumOpenmpThreads(&hw_opt, *hwinfo, cr, ms, physicalNodeComm.size_,
1255                                             pmeRunMode, mtop, *inputrec);
1256
1257     gmx_omp_nthreads_init(mdlog, cr, hwinfo->nthreads_hw_avail, physicalNodeComm.size_,
1258                           hw_opt.nthreads_omp, hw_opt.nthreads_omp_pme, !thisRankHasDuty(cr, DUTY_PP));
1259
1260     // Enable FP exception detection, but not in
1261     // Release mode and not for compilers with known buggy FP
1262     // exception support (clang with any optimization) or suspected
1263     // buggy FP exception support (gcc 7.* with optimization).
1264 #if !defined NDEBUG                                                                         \
1265         && !((defined __clang__ || (defined(__GNUC__) && !defined(__ICC) && __GNUC__ == 7)) \
1266              && defined __OPTIMIZE__)
1267     const bool bEnableFPE = true;
1268 #else
1269     const bool bEnableFPE = false;
1270 #endif
1271     // FIXME - reconcile with gmx_feenableexcept() call from CommandLineModuleManager::run()
1272     if (bEnableFPE)
1273     {
1274         gmx_feenableexcept();
1275     }
1276
1277     /* Now that we know the setup is consistent, check for efficiency */
1278     check_resource_division_efficiency(hwinfo, gpuTaskAssignments.thisRankHasAnyGpuTask(),
1279                                        mdrunOptions.ntompOptionIsSet, cr, mdlog);
1280
1281     /* getting number of PP/PME threads on this MPI / tMPI rank.
1282        PME: env variable should be read only on one node to make sure it is
1283        identical everywhere;
1284      */
1285     const int numThreadsOnThisRank = thisRankHasDuty(cr, DUTY_PP) ? gmx_omp_nthreads_get(emntNonbonded)
1286                                                                   : gmx_omp_nthreads_get(emntPME);
1287     checkHardwareOversubscription(numThreadsOnThisRank, cr->nodeid, *hwinfo->hardwareTopology,
1288                                   physicalNodeComm, mdlog);
1289
1290     // Enable Peer access between GPUs where available
1291     // Only for DD, only master PP rank needs to perform setup, and only if thread MPI plus
1292     // any of the GPU communication features are active.
1293     if (DOMAINDECOMP(cr) && MASTER(cr) && thisRankHasDuty(cr, DUTY_PP) && GMX_THREAD_MPI
1294         && (devFlags.enableGpuHaloExchange || devFlags.enableGpuPmePPComm))
1295     {
1296         setupGpuDevicePeerAccess(gpuIdsToUse, mdlog);
1297     }
1298
1299     if (hw_opt.threadAffinity != ThreadAffinity::Off)
1300     {
1301         /* Before setting affinity, check whether the affinity has changed
1302          * - which indicates that probably the OpenMP library has changed it
1303          * since we first checked).
1304          */
1305         gmx_check_thread_affinity_set(mdlog, &hw_opt, hwinfo->nthreads_hw_avail, TRUE);
1306
1307         int numThreadsOnThisNode, intraNodeThreadOffset;
1308         analyzeThreadsOnThisNode(physicalNodeComm, numThreadsOnThisRank, &numThreadsOnThisNode,
1309                                  &intraNodeThreadOffset);
1310
1311         /* Set the CPU affinity */
1312         gmx_set_thread_affinity(mdlog, cr, &hw_opt, *hwinfo->hardwareTopology, numThreadsOnThisRank,
1313                                 numThreadsOnThisNode, intraNodeThreadOffset, nullptr);
1314     }
1315
1316     if (mdrunOptions.timingOptions.resetStep > -1)
1317     {
1318         GMX_LOG(mdlog.info)
1319                 .asParagraph()
1320                 .appendText(
1321                         "The -resetstep functionality is deprecated, and may be removed in a "
1322                         "future version.");
1323     }
1324     wcycle = wallcycle_init(fplog, mdrunOptions.timingOptions.resetStep, cr);
1325
1326     if (PAR(cr))
1327     {
1328         /* Master synchronizes its value of reset_counters with all nodes
1329          * including PME only nodes */
1330         int64_t reset_counters = wcycle_get_reset_counters(wcycle);
1331         gmx_bcast_sim(sizeof(reset_counters), &reset_counters, cr);
1332         wcycle_set_reset_counters(wcycle, reset_counters);
1333     }
1334
1335     // Membrane embedding must be initialized before we call init_forcerec()
1336     if (doMembed)
1337     {
1338         if (MASTER(cr))
1339         {
1340             fprintf(stderr, "Initializing membed");
1341         }
1342         /* Note that membed cannot work in parallel because mtop is
1343          * changed here. Fix this if we ever want to make it run with
1344          * multiple ranks. */
1345         membed = init_membed(fplog, filenames.size(), filenames.data(), &mtop, inputrec,
1346                              globalState.get(), cr, &mdrunOptions.checkpointOptions.period);
1347     }
1348
1349     const bool                   thisRankHasPmeGpuTask = gpuTaskAssignments.thisRankHasPmeGpuTask();
1350     std::unique_ptr<MDAtoms>     mdAtoms;
1351     std::unique_ptr<gmx_vsite_t> vsite;
1352     std::unique_ptr<GpuBonded>   gpuBonded;
1353
1354     t_nrnb nrnb;
1355     if (thisRankHasDuty(cr, DUTY_PP))
1356     {
1357         mdModulesNotifier.notify(*cr);
1358         mdModulesNotifier.notify(&atomSets);
1359         mdModulesNotifier.notify(inputrec->pbcType);
1360         mdModulesNotifier.notify(SimulationTimeStep{ inputrec->delta_t });
1361         /* Initiate forcerecord */
1362         fr                 = new t_forcerec;
1363         fr->forceProviders = mdModules_->initForceProviders();
1364         init_forcerec(fplog, mdlog, fr, fcd, inputrec, &mtop, cr, box,
1365                       opt2fn("-table", filenames.size(), filenames.data()),
1366                       opt2fn("-tablep", filenames.size(), filenames.data()),
1367                       opt2fns("-tableb", filenames.size(), filenames.data()), pforce);
1368
1369         // Save a handle to device stream manager to use elsewhere in the code
1370         // TODO: Forcerec is not a correct place to store it.
1371         fr->deviceStreamManager = deviceStreamManager.get();
1372
1373         if (devFlags.enableGpuPmePPComm && !thisRankHasDuty(cr, DUTY_PME))
1374         {
1375             GMX_RELEASE_ASSERT(
1376                     deviceStreamManager != nullptr,
1377                     "GPU device stream manager should be valid in order to use PME-PP direct "
1378                     "communications.");
1379             GMX_RELEASE_ASSERT(
1380                     deviceStreamManager->streamIsValid(DeviceStreamType::PmePpTransfer),
1381                     "GPU PP-PME stream should be valid in order to use GPU PME-PP direct "
1382                     "communications.");
1383             fr->pmePpCommGpu = std::make_unique<gmx::PmePpCommGpu>(
1384                     cr->mpi_comm_mysim, cr->dd->pme_nodeid, deviceStreamManager->context(),
1385                     deviceStreamManager->stream(DeviceStreamType::PmePpTransfer));
1386         }
1387
1388         fr->nbv = Nbnxm::init_nb_verlet(mdlog, inputrec, fr, cr, *hwinfo, useGpuForNonbonded,
1389                                         deviceStreamManager.get(), &mtop, box, wcycle);
1390         // TODO: Move the logic below to a GPU bonded builder
1391         if (useGpuForBonded)
1392         {
1393             GMX_RELEASE_ASSERT(deviceStreamManager != nullptr,
1394                                "GPU device stream manager should be valid in order to use GPU "
1395                                "version of bonded forces.");
1396             gpuBonded = std::make_unique<GpuBonded>(
1397                     mtop.ffparams, deviceStreamManager->context(),
1398                     deviceStreamManager->bondedStream(havePPDomainDecomposition(cr)), wcycle);
1399             fr->gpuBonded = gpuBonded.get();
1400         }
1401
1402         /* Initialize the mdAtoms structure.
1403          * mdAtoms is not filled with atom data,
1404          * as this can not be done now with domain decomposition.
1405          */
1406         mdAtoms = makeMDAtoms(fplog, mtop, *inputrec, thisRankHasPmeGpuTask);
1407         if (globalState && thisRankHasPmeGpuTask)
1408         {
1409             // The pinning of coordinates in the global state object works, because we only use
1410             // PME on GPU without DD or on a separate PME rank, and because the local state pointer
1411             // points to the global state object without DD.
1412             // FIXME: MD and EM separately set up the local state - this should happen in the same
1413             // function, which should also perform the pinning.
1414             changePinningPolicy(&globalState->x, pme_get_pinning_policy());
1415         }
1416
1417         /* Initialize the virtual site communication */
1418         vsite = initVsite(mtop, cr);
1419
1420         calc_shifts(box, fr->shift_vec);
1421
1422         /* With periodic molecules the charge groups should be whole at start up
1423          * and the virtual sites should not be far from their proper positions.
1424          */
1425         if (!inputrec->bContinuation && MASTER(cr)
1426             && !(inputrec->pbcType != PbcType::No && inputrec->bPeriodicMols))
1427         {
1428             /* Make molecules whole at start of run */
1429             if (fr->pbcType != PbcType::No)
1430             {
1431                 do_pbc_first_mtop(fplog, inputrec->pbcType, box, &mtop, globalState->x.rvec_array());
1432             }
1433             if (vsite)
1434             {
1435                 /* Correct initial vsite positions are required
1436                  * for the initial distribution in the domain decomposition
1437                  * and for the initial shell prediction.
1438                  */
1439                 constructVsitesGlobal(mtop, globalState->x);
1440             }
1441         }
1442
1443         if (EEL_PME(fr->ic->eeltype) || EVDW_PME(fr->ic->vdwtype))
1444         {
1445             ewaldcoeff_q  = fr->ic->ewaldcoeff_q;
1446             ewaldcoeff_lj = fr->ic->ewaldcoeff_lj;
1447         }
1448     }
1449     else
1450     {
1451         /* This is a PME only node */
1452
1453         GMX_ASSERT(globalState == nullptr,
1454                    "We don't need the state on a PME only rank and expect it to be unitialized");
1455
1456         ewaldcoeff_q  = calc_ewaldcoeff_q(inputrec->rcoulomb, inputrec->ewald_rtol);
1457         ewaldcoeff_lj = calc_ewaldcoeff_lj(inputrec->rvdw, inputrec->ewald_rtol_lj);
1458     }
1459
1460     gmx_pme_t* sepPmeData = nullptr;
1461     // This reference hides the fact that PME data is owned by runner on PME-only ranks and by forcerec on other ranks
1462     GMX_ASSERT(thisRankHasDuty(cr, DUTY_PP) == (fr != nullptr),
1463                "Double-checking that only PME-only ranks have no forcerec");
1464     gmx_pme_t*& pmedata = fr ? fr->pmedata : sepPmeData;
1465
1466     // TODO should live in ewald module once its testing is improved
1467     //
1468     // Later, this program could contain kernels that might be later
1469     // re-used as auto-tuning progresses, or subsequent simulations
1470     // are invoked.
1471     PmeGpuProgramStorage pmeGpuProgram;
1472     if (thisRankHasPmeGpuTask)
1473     {
1474         GMX_RELEASE_ASSERT(
1475                 (deviceStreamManager != nullptr),
1476                 "GPU device stream manager should be initialized in order to use GPU for PME.");
1477         GMX_RELEASE_ASSERT((deviceInfo != nullptr),
1478                            "GPU device should be initialized in order to use GPU for PME.");
1479         pmeGpuProgram = buildPmeGpuProgram(deviceStreamManager->context());
1480     }
1481
1482     /* Initiate PME if necessary,
1483      * either on all nodes or on dedicated PME nodes only. */
1484     if (EEL_PME(inputrec->coulombtype) || EVDW_PME(inputrec->vdwtype))
1485     {
1486         if (mdAtoms && mdAtoms->mdatoms())
1487         {
1488             nChargePerturbed = mdAtoms->mdatoms()->nChargePerturbed;
1489             if (EVDW_PME(inputrec->vdwtype))
1490             {
1491                 nTypePerturbed = mdAtoms->mdatoms()->nTypePerturbed;
1492             }
1493         }
1494         if (cr->npmenodes > 0)
1495         {
1496             /* The PME only nodes need to know nChargePerturbed(FEP on Q) and nTypePerturbed(FEP on LJ)*/
1497             gmx_bcast_sim(sizeof(nChargePerturbed), &nChargePerturbed, cr);
1498             gmx_bcast_sim(sizeof(nTypePerturbed), &nTypePerturbed, cr);
1499         }
1500
1501         if (thisRankHasDuty(cr, DUTY_PME))
1502         {
1503             try
1504             {
1505                 // TODO: This should be in the builder.
1506                 GMX_RELEASE_ASSERT(!useGpuForPme || (deviceStreamManager != nullptr),
1507                                    "Device stream manager should be valid in order to use GPU "
1508                                    "version of PME.");
1509                 GMX_RELEASE_ASSERT(
1510                         !useGpuForPme || deviceStreamManager->streamIsValid(DeviceStreamType::Pme),
1511                         "GPU PME stream should be valid in order to use GPU version of PME.");
1512
1513                 const DeviceContext* deviceContext =
1514                         useGpuForPme ? &deviceStreamManager->context() : nullptr;
1515                 const DeviceStream* pmeStream =
1516                         useGpuForPme ? &deviceStreamManager->stream(DeviceStreamType::Pme) : nullptr;
1517
1518                 pmedata = gmx_pme_init(cr, getNumPmeDomains(cr->dd), inputrec, nChargePerturbed != 0,
1519                                        nTypePerturbed != 0, mdrunOptions.reproducible, ewaldcoeff_q,
1520                                        ewaldcoeff_lj, gmx_omp_nthreads_get(emntPME), pmeRunMode,
1521                                        nullptr, deviceContext, pmeStream, pmeGpuProgram.get(), mdlog);
1522             }
1523             GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
1524         }
1525     }
1526
1527
1528     if (EI_DYNAMICS(inputrec->eI))
1529     {
1530         /* Turn on signal handling on all nodes */
1531         /*
1532          * (A user signal from the PME nodes (if any)
1533          * is communicated to the PP nodes.
1534          */
1535         signal_handler_install();
1536     }
1537
1538     pull_t* pull_work = nullptr;
1539     if (thisRankHasDuty(cr, DUTY_PP))
1540     {
1541         /* Assumes uniform use of the number of OpenMP threads */
1542         walltime_accounting = walltime_accounting_init(gmx_omp_nthreads_get(emntDefault));
1543
1544         if (inputrec->bPull)
1545         {
1546             /* Initialize pull code */
1547             pull_work = init_pull(fplog, inputrec->pull, inputrec, &mtop, cr, &atomSets,
1548                                   inputrec->fepvals->init_lambda);
1549             if (inputrec->pull->bXOutAverage || inputrec->pull->bFOutAverage)
1550             {
1551                 initPullHistory(pull_work, &observablesHistory);
1552             }
1553             if (EI_DYNAMICS(inputrec->eI) && MASTER(cr))
1554             {
1555                 init_pull_output_files(pull_work, filenames.size(), filenames.data(), oenv, startingBehavior);
1556             }
1557         }
1558
1559         std::unique_ptr<EnforcedRotation> enforcedRotation;
1560         if (inputrec->bRot)
1561         {
1562             /* Initialize enforced rotation code */
1563             enforcedRotation =
1564                     init_rot(fplog, inputrec, filenames.size(), filenames.data(), cr, &atomSets,
1565                              globalState.get(), &mtop, oenv, mdrunOptions, startingBehavior);
1566         }
1567
1568         t_swap* swap = nullptr;
1569         if (inputrec->eSwapCoords != eswapNO)
1570         {
1571             /* Initialize ion swapping code */
1572             swap = init_swapcoords(fplog, inputrec,
1573                                    opt2fn_master("-swap", filenames.size(), filenames.data(), cr),
1574                                    &mtop, globalState.get(), &observablesHistory, cr, &atomSets,
1575                                    oenv, mdrunOptions, startingBehavior);
1576         }
1577
1578         /* Let makeConstraints know whether we have essential dynamics constraints. */
1579         auto constr = makeConstraints(mtop, *inputrec, pull_work, doEssentialDynamics, fplog,
1580                                       *mdAtoms->mdatoms(), cr, ms, &nrnb, wcycle, fr->bMolPBC);
1581
1582         /* Energy terms and groups */
1583         gmx_enerdata_t enerd(mtop.groups.groups[SimulationAtomGroupType::EnergyOutput].size(),
1584                              inputrec->fepvals->n_lambda);
1585
1586         /* Kinetic energy data */
1587         gmx_ekindata_t ekind;
1588         init_ekindata(fplog, &mtop, &(inputrec->opts), &ekind, inputrec->cos_accel);
1589
1590         /* Set up interactive MD (IMD) */
1591         auto imdSession =
1592                 makeImdSession(inputrec, cr, wcycle, &enerd, ms, &mtop, mdlog,
1593                                MASTER(cr) ? globalState->x.rvec_array() : nullptr, filenames.size(),
1594                                filenames.data(), oenv, mdrunOptions.imdOptions, startingBehavior);
1595
1596         if (DOMAINDECOMP(cr))
1597         {
1598             GMX_RELEASE_ASSERT(fr, "fr was NULL while cr->duty was DUTY_PP");
1599             /* This call is not included in init_domain_decomposition mainly
1600              * because fr->cginfo_mb is set later.
1601              */
1602             dd_init_bondeds(fplog, cr->dd, mtop, vsite.get(), inputrec,
1603                             domdecOptions.checkBondedInteractions, fr->cginfo_mb);
1604         }
1605
1606         // TODO This is not the right place to manage the lifetime of
1607         // this data structure, but currently it's the easiest way to
1608         // make it work.
1609         MdrunScheduleWorkload runScheduleWork;
1610         // Also populates the simulation constant workload description.
1611         runScheduleWork.simulationWork =
1612                 createSimulationWorkload(*inputrec, useGpuForNonbonded, pmeRunMode, useGpuForBonded,
1613                                          useGpuForUpdate, devFlags.enableGpuBufferOps,
1614                                          devFlags.enableGpuHaloExchange, devFlags.enableGpuPmePPComm);
1615
1616         std::unique_ptr<gmx::StatePropagatorDataGpu> stateGpu;
1617         if (gpusWereDetected
1618             && ((useGpuForPme && thisRankHasDuty(cr, DUTY_PME))
1619                 || runScheduleWork.simulationWork.useGpuBufferOps))
1620         {
1621             GpuApiCallBehavior transferKind = (inputrec->eI == eiMD && !doRerun && !useModularSimulator)
1622                                                       ? GpuApiCallBehavior::Async
1623                                                       : GpuApiCallBehavior::Sync;
1624             GMX_RELEASE_ASSERT(deviceStreamManager != nullptr,
1625                                "GPU device stream manager should be initialized to use GPU.");
1626             stateGpu = std::make_unique<gmx::StatePropagatorDataGpu>(
1627                     *deviceStreamManager, transferKind, pme_gpu_get_block_size(fr->pmedata), wcycle);
1628             fr->stateGpu = stateGpu.get();
1629         }
1630
1631         GMX_ASSERT(stopHandlerBuilder_, "Runner must provide StopHandlerBuilder to simulator.");
1632         SimulatorBuilder simulatorBuilder;
1633
1634         // build and run simulator object based on user-input
1635         auto simulator = simulatorBuilder.build(
1636                 useModularSimulator, fplog, cr, ms, mdlog, static_cast<int>(filenames.size()),
1637                 filenames.data(), oenv, mdrunOptions, startingBehavior, vsite.get(), constr.get(),
1638                 enforcedRotation ? enforcedRotation->getLegacyEnfrot() : nullptr, deform.get(),
1639                 mdModules_->outputProvider(), mdModules_->notifier(), inputrec, imdSession.get(),
1640                 pull_work, swap, &mtop, fcd, globalState.get(), &observablesHistory, mdAtoms.get(),
1641                 &nrnb, wcycle, fr, &enerd, &ekind, &runScheduleWork, replExParams, membed,
1642                 walltime_accounting, std::move(stopHandlerBuilder_), doRerun);
1643         simulator->run();
1644
1645         if (fr->pmePpCommGpu)
1646         {
1647             // destroy object since it is no longer required. (This needs to be done while the GPU context still exists.)
1648             fr->pmePpCommGpu.reset();
1649         }
1650
1651         if (inputrec->bPull)
1652         {
1653             finish_pull(pull_work);
1654         }
1655         finish_swapcoords(swap);
1656     }
1657     else
1658     {
1659         GMX_RELEASE_ASSERT(pmedata, "pmedata was NULL while cr->duty was not DUTY_PP");
1660         /* do PME only */
1661         walltime_accounting = walltime_accounting_init(gmx_omp_nthreads_get(emntPME));
1662         gmx_pmeonly(pmedata, cr, &nrnb, wcycle, walltime_accounting, inputrec, pmeRunMode,
1663                     deviceStreamManager.get());
1664     }
1665
1666     wallcycle_stop(wcycle, ewcRUN);
1667
1668     /* Finish up, write some stuff
1669      * if rerunMD, don't write last frame again
1670      */
1671     finish_run(fplog, mdlog, cr, inputrec, &nrnb, wcycle, walltime_accounting,
1672                fr ? fr->nbv.get() : nullptr, pmedata, EI_DYNAMICS(inputrec->eI) && !isMultiSim(ms));
1673
1674     // clean up cycle counter
1675     wallcycle_destroy(wcycle);
1676
1677     deviceStreamManager.reset(nullptr);
1678     // Free PME data
1679     if (pmedata)
1680     {
1681         gmx_pme_destroy(pmedata);
1682         pmedata = nullptr;
1683     }
1684
1685     // FIXME: this is only here to manually unpin mdAtoms->chargeA_ and state->x,
1686     // before we destroy the GPU context(s) in free_gpu().
1687     // Pinned buffers are associated with contexts in CUDA.
1688     // As soon as we destroy GPU contexts after mdrunner() exits, these lines should go.
1689     mdAtoms.reset(nullptr);
1690     globalState.reset(nullptr);
1691     mdModules_.reset(nullptr); // destruct force providers here as they might also use the GPU
1692     gpuBonded.reset(nullptr);
1693     /* Free pinned buffers in *fr */
1694     delete fr;
1695     fr = nullptr;
1696
1697     if (hwinfo->gpu_info.n_dev > 0)
1698     {
1699         /* stop the GPU profiler (only CUDA) */
1700         stopGpuProfiler();
1701     }
1702
1703     /* With tMPI we need to wait for all ranks to finish deallocation before
1704      * destroying the CUDA context in free_gpu() as some tMPI ranks may be sharing
1705      * GPU and context.
1706      *
1707      * This is not a concern in OpenCL where we use one context per rank which
1708      * is freed in nbnxn_gpu_free().
1709      *
1710      * Note: it is safe to not call the barrier on the ranks which do not use GPU,
1711      * but it is easier and more futureproof to call it on the whole node.
1712      *
1713      * Note that this function needs to be called even if GPUs are not used
1714      * in this run because the PME ranks have no knowledge of whether GPUs
1715      * are used or not, but all ranks need to enter the barrier below.
1716      * \todo Remove this physical node barrier after making sure
1717      * that it's not needed anymore (with a shared GPU run).
1718      */
1719     if (GMX_THREAD_MPI)
1720     {
1721         physicalNodeComm.barrier();
1722     }
1723
1724     free_gpu(deviceInfo);
1725     sfree(fcd);
1726
1727     if (doMembed)
1728     {
1729         free_membed(membed);
1730     }
1731
1732     /* Does what it says */
1733     print_date_and_time(fplog, cr->nodeid, "Finished mdrun", gmx_gettime());
1734     walltime_accounting_destroy(walltime_accounting);
1735
1736     // Ensure log file content is written
1737     if (logFileHandle)
1738     {
1739         gmx_fio_flush(logFileHandle);
1740     }
1741
1742     /* Reset FPEs (important for unit tests) by disabling them. Assumes no
1743      * exceptions were enabled before function was called. */
1744     if (bEnableFPE)
1745     {
1746         gmx_fedisableexcept();
1747     }
1748
1749     auto rc = static_cast<int>(gmx_get_stop_condition());
1750
1751 #if GMX_THREAD_MPI
1752     /* we need to join all threads. The sub-threads join when they
1753        exit this function, but the master thread needs to be told to
1754        wait for that. */
1755     if (PAR(cr) && MASTER(cr))
1756     {
1757         tMPI_Finalize();
1758     }
1759 #endif
1760     return rc;
1761 } // namespace gmx
1762
1763 Mdrunner::~Mdrunner()
1764 {
1765     // Clean up of the Manager.
1766     // This will end up getting called on every thread-MPI rank, which is unnecessary,
1767     // but okay as long as threads synchronize some time before adding or accessing
1768     // a new set of restraints.
1769     if (restraintManager_)
1770     {
1771         restraintManager_->clear();
1772         GMX_ASSERT(restraintManager_->countRestraints() == 0,
1773                    "restraints added during runner life time should be cleared at runner "
1774                    "destruction.");
1775     }
1776 };
1777
1778 void Mdrunner::addPotential(std::shared_ptr<gmx::IRestraintPotential> puller, const std::string& name)
1779 {
1780     GMX_ASSERT(restraintManager_, "Mdrunner must have a restraint manager.");
1781     // Not sure if this should be logged through the md logger or something else,
1782     // but it is helpful to have some sort of INFO level message sent somewhere.
1783     //    std::cout << "Registering restraint named " << name << std::endl;
1784
1785     // When multiple restraints are used, it may be wasteful to register them separately.
1786     // Maybe instead register an entire Restraint Manager as a force provider.
1787     restraintManager_->addToSpec(std::move(puller), name);
1788 }
1789
1790 Mdrunner::Mdrunner(std::unique_ptr<MDModules> mdModules) : mdModules_(std::move(mdModules)) {}
1791
1792 Mdrunner::Mdrunner(Mdrunner&&) noexcept = default;
1793
1794 //NOLINTNEXTLINE(performance-noexcept-move-constructor) working around GCC bug 58265
1795 Mdrunner& Mdrunner::operator=(Mdrunner&& /*handle*/) noexcept(BUGFREE_NOEXCEPT_STRING) = default;
1796
1797 class Mdrunner::BuilderImplementation
1798 {
1799 public:
1800     BuilderImplementation() = delete;
1801     BuilderImplementation(std::unique_ptr<MDModules> mdModules, compat::not_null<SimulationContext*> context);
1802     ~BuilderImplementation();
1803
1804     BuilderImplementation& setExtraMdrunOptions(const MdrunOptions& options,
1805                                                 real                forceWarningThreshold,
1806                                                 StartingBehavior    startingBehavior);
1807
1808     void addDomdec(const DomdecOptions& options);
1809
1810     void addVerletList(int nstlist);
1811
1812     void addReplicaExchange(const ReplicaExchangeParameters& params);
1813
1814     void addNonBonded(const char* nbpu_opt);
1815
1816     void addPME(const char* pme_opt_, const char* pme_fft_opt_);
1817
1818     void addBondedTaskAssignment(const char* bonded_opt);
1819
1820     void addUpdateTaskAssignment(const char* update_opt);
1821
1822     void addHardwareOptions(const gmx_hw_opt_t& hardwareOptions);
1823
1824     void addFilenames(ArrayRef<const t_filenm> filenames);
1825
1826     void addOutputEnvironment(gmx_output_env_t* outputEnvironment);
1827
1828     void addLogFile(t_fileio* logFileHandle);
1829
1830     void addStopHandlerBuilder(std::unique_ptr<StopHandlerBuilder> builder);
1831
1832     Mdrunner build();
1833
1834 private:
1835     // Default parameters copied from runner.h
1836     // \todo Clarify source(s) of default parameters.
1837
1838     const char* nbpu_opt_    = nullptr;
1839     const char* pme_opt_     = nullptr;
1840     const char* pme_fft_opt_ = nullptr;
1841     const char* bonded_opt_  = nullptr;
1842     const char* update_opt_  = nullptr;
1843
1844     MdrunOptions mdrunOptions_;
1845
1846     DomdecOptions domdecOptions_;
1847
1848     ReplicaExchangeParameters replicaExchangeParameters_;
1849
1850     //! Command-line override for the duration of a neighbor list with the Verlet scheme.
1851     int nstlist_ = 0;
1852
1853     //! Multisim communicator handle.
1854     gmx_multisim_t* multiSimulation_;
1855
1856     //! mdrun communicator
1857     MPI_Comm communicator_ = MPI_COMM_NULL;
1858
1859     //! Print a warning if any force is larger than this (in kJ/mol nm).
1860     real forceWarningThreshold_ = -1;
1861
1862     //! Whether the simulation will start afresh, or restart with/without appending.
1863     StartingBehavior startingBehavior_ = StartingBehavior::NewSimulation;
1864
1865     //! The modules that comprise the functionality of mdrun.
1866     std::unique_ptr<MDModules> mdModules_;
1867
1868     //! \brief Parallelism information.
1869     gmx_hw_opt_t hardwareOptions_;
1870
1871     //! filename options for simulation.
1872     ArrayRef<const t_filenm> filenames_;
1873
1874     /*! \brief Handle to output environment.
1875      *
1876      * \todo gmx_output_env_t needs lifetime management.
1877      */
1878     gmx_output_env_t* outputEnvironment_ = nullptr;
1879
1880     /*! \brief Non-owning handle to MD log file.
1881      *
1882      * \todo Context should own output facilities for client.
1883      * \todo Improve log file handle management.
1884      * \internal
1885      * Code managing the FILE* relies on the ability to set it to
1886      * nullptr to check whether the filehandle is valid.
1887      */
1888     t_fileio* logFileHandle_ = nullptr;
1889
1890     /*!
1891      * \brief Builder for simulation stop signal handler.
1892      */
1893     std::unique_ptr<StopHandlerBuilder> stopHandlerBuilder_ = nullptr;
1894 };
1895
1896 Mdrunner::BuilderImplementation::BuilderImplementation(std::unique_ptr<MDModules> mdModules,
1897                                                        compat::not_null<SimulationContext*> context) :
1898     mdModules_(std::move(mdModules))
1899 {
1900     communicator_    = context->communicator_;
1901     multiSimulation_ = context->multiSimulation_.get();
1902 }
1903
1904 Mdrunner::BuilderImplementation::~BuilderImplementation() = default;
1905
1906 Mdrunner::BuilderImplementation&
1907 Mdrunner::BuilderImplementation::setExtraMdrunOptions(const MdrunOptions&    options,
1908                                                       const real             forceWarningThreshold,
1909                                                       const StartingBehavior startingBehavior)
1910 {
1911     mdrunOptions_          = options;
1912     forceWarningThreshold_ = forceWarningThreshold;
1913     startingBehavior_      = startingBehavior;
1914     return *this;
1915 }
1916
1917 void Mdrunner::BuilderImplementation::addDomdec(const DomdecOptions& options)
1918 {
1919     domdecOptions_ = options;
1920 }
1921
1922 void Mdrunner::BuilderImplementation::addVerletList(int nstlist)
1923 {
1924     nstlist_ = nstlist;
1925 }
1926
1927 void Mdrunner::BuilderImplementation::addReplicaExchange(const ReplicaExchangeParameters& params)
1928 {
1929     replicaExchangeParameters_ = params;
1930 }
1931
1932 Mdrunner Mdrunner::BuilderImplementation::build()
1933 {
1934     auto newRunner = Mdrunner(std::move(mdModules_));
1935
1936     newRunner.mdrunOptions     = mdrunOptions_;
1937     newRunner.pforce           = forceWarningThreshold_;
1938     newRunner.startingBehavior = startingBehavior_;
1939     newRunner.domdecOptions    = domdecOptions_;
1940
1941     // \todo determine an invariant to check or confirm that all gmx_hw_opt_t objects are valid
1942     newRunner.hw_opt = hardwareOptions_;
1943
1944     // No invariant to check. This parameter exists to optionally override other behavior.
1945     newRunner.nstlist_cmdline = nstlist_;
1946
1947     newRunner.replExParams = replicaExchangeParameters_;
1948
1949     newRunner.filenames = filenames_;
1950
1951     newRunner.communicator = communicator_;
1952
1953     // nullptr is a valid value for the multisim handle
1954     newRunner.ms = multiSimulation_;
1955
1956     // \todo Clarify ownership and lifetime management for gmx_output_env_t
1957     // \todo Update sanity checking when output environment has clearly specified invariants.
1958     // Initialization and default values for oenv are not well specified in the current version.
1959     if (outputEnvironment_)
1960     {
1961         newRunner.oenv = outputEnvironment_;
1962     }
1963     else
1964     {
1965         GMX_THROW(gmx::APIError(
1966                 "MdrunnerBuilder::addOutputEnvironment() is required before build()"));
1967     }
1968
1969     newRunner.logFileHandle = logFileHandle_;
1970
1971     if (nbpu_opt_)
1972     {
1973         newRunner.nbpu_opt = nbpu_opt_;
1974     }
1975     else
1976     {
1977         GMX_THROW(gmx::APIError("MdrunnerBuilder::addNonBonded() is required before build()"));
1978     }
1979
1980     if (pme_opt_ && pme_fft_opt_)
1981     {
1982         newRunner.pme_opt     = pme_opt_;
1983         newRunner.pme_fft_opt = pme_fft_opt_;
1984     }
1985     else
1986     {
1987         GMX_THROW(gmx::APIError("MdrunnerBuilder::addElectrostatics() is required before build()"));
1988     }
1989
1990     if (bonded_opt_)
1991     {
1992         newRunner.bonded_opt = bonded_opt_;
1993     }
1994     else
1995     {
1996         GMX_THROW(gmx::APIError(
1997                 "MdrunnerBuilder::addBondedTaskAssignment() is required before build()"));
1998     }
1999
2000     if (update_opt_)
2001     {
2002         newRunner.update_opt = update_opt_;
2003     }
2004     else
2005     {
2006         GMX_THROW(gmx::APIError(
2007                 "MdrunnerBuilder::addUpdateTaskAssignment() is required before build()  "));
2008     }
2009
2010
2011     newRunner.restraintManager_ = std::make_unique<gmx::RestraintManager>();
2012
2013     if (stopHandlerBuilder_)
2014     {
2015         newRunner.stopHandlerBuilder_ = std::move(stopHandlerBuilder_);
2016     }
2017     else
2018     {
2019         newRunner.stopHandlerBuilder_ = std::make_unique<StopHandlerBuilder>();
2020     }
2021
2022     return newRunner;
2023 }
2024
2025 void Mdrunner::BuilderImplementation::addNonBonded(const char* nbpu_opt)
2026 {
2027     nbpu_opt_ = nbpu_opt;
2028 }
2029
2030 void Mdrunner::BuilderImplementation::addPME(const char* pme_opt, const char* pme_fft_opt)
2031 {
2032     pme_opt_     = pme_opt;
2033     pme_fft_opt_ = pme_fft_opt;
2034 }
2035
2036 void Mdrunner::BuilderImplementation::addBondedTaskAssignment(const char* bonded_opt)
2037 {
2038     bonded_opt_ = bonded_opt;
2039 }
2040
2041 void Mdrunner::BuilderImplementation::addUpdateTaskAssignment(const char* update_opt)
2042 {
2043     update_opt_ = update_opt;
2044 }
2045
2046 void Mdrunner::BuilderImplementation::addHardwareOptions(const gmx_hw_opt_t& hardwareOptions)
2047 {
2048     hardwareOptions_ = hardwareOptions;
2049 }
2050
2051 void Mdrunner::BuilderImplementation::addFilenames(ArrayRef<const t_filenm> filenames)
2052 {
2053     filenames_ = filenames;
2054 }
2055
2056 void Mdrunner::BuilderImplementation::addOutputEnvironment(gmx_output_env_t* outputEnvironment)
2057 {
2058     outputEnvironment_ = outputEnvironment;
2059 }
2060
2061 void Mdrunner::BuilderImplementation::addLogFile(t_fileio* logFileHandle)
2062 {
2063     logFileHandle_ = logFileHandle;
2064 }
2065
2066 void Mdrunner::BuilderImplementation::addStopHandlerBuilder(std::unique_ptr<StopHandlerBuilder> builder)
2067 {
2068     stopHandlerBuilder_ = std::move(builder);
2069 }
2070
2071 MdrunnerBuilder::MdrunnerBuilder(std::unique_ptr<MDModules>           mdModules,
2072                                  compat::not_null<SimulationContext*> context) :
2073     impl_{ std::make_unique<Mdrunner::BuilderImplementation>(std::move(mdModules), context) }
2074 {
2075 }
2076
2077 MdrunnerBuilder::~MdrunnerBuilder() = default;
2078
2079 MdrunnerBuilder& MdrunnerBuilder::addSimulationMethod(const MdrunOptions&    options,
2080                                                       real                   forceWarningThreshold,
2081                                                       const StartingBehavior startingBehavior)
2082 {
2083     impl_->setExtraMdrunOptions(options, forceWarningThreshold, startingBehavior);
2084     return *this;
2085 }
2086
2087 MdrunnerBuilder& MdrunnerBuilder::addDomainDecomposition(const DomdecOptions& options)
2088 {
2089     impl_->addDomdec(options);
2090     return *this;
2091 }
2092
2093 MdrunnerBuilder& MdrunnerBuilder::addNeighborList(int nstlist)
2094 {
2095     impl_->addVerletList(nstlist);
2096     return *this;
2097 }
2098
2099 MdrunnerBuilder& MdrunnerBuilder::addReplicaExchange(const ReplicaExchangeParameters& params)
2100 {
2101     impl_->addReplicaExchange(params);
2102     return *this;
2103 }
2104
2105 MdrunnerBuilder& MdrunnerBuilder::addNonBonded(const char* nbpu_opt)
2106 {
2107     impl_->addNonBonded(nbpu_opt);
2108     return *this;
2109 }
2110
2111 MdrunnerBuilder& MdrunnerBuilder::addElectrostatics(const char* pme_opt, const char* pme_fft_opt)
2112 {
2113     // The builder method may become more general in the future, but in this version,
2114     // parameters for PME electrostatics are both required and the only parameters
2115     // available.
2116     if (pme_opt && pme_fft_opt)
2117     {
2118         impl_->addPME(pme_opt, pme_fft_opt);
2119     }
2120     else
2121     {
2122         GMX_THROW(
2123                 gmx::InvalidInputError("addElectrostatics() arguments must be non-null pointers."));
2124     }
2125     return *this;
2126 }
2127
2128 MdrunnerBuilder& MdrunnerBuilder::addBondedTaskAssignment(const char* bonded_opt)
2129 {
2130     impl_->addBondedTaskAssignment(bonded_opt);
2131     return *this;
2132 }
2133
2134 MdrunnerBuilder& MdrunnerBuilder::addUpdateTaskAssignment(const char* update_opt)
2135 {
2136     impl_->addUpdateTaskAssignment(update_opt);
2137     return *this;
2138 }
2139
2140 Mdrunner MdrunnerBuilder::build()
2141 {
2142     return impl_->build();
2143 }
2144
2145 MdrunnerBuilder& MdrunnerBuilder::addHardwareOptions(const gmx_hw_opt_t& hardwareOptions)
2146 {
2147     impl_->addHardwareOptions(hardwareOptions);
2148     return *this;
2149 }
2150
2151 MdrunnerBuilder& MdrunnerBuilder::addFilenames(ArrayRef<const t_filenm> filenames)
2152 {
2153     impl_->addFilenames(filenames);
2154     return *this;
2155 }
2156
2157 MdrunnerBuilder& MdrunnerBuilder::addOutputEnvironment(gmx_output_env_t* outputEnvironment)
2158 {
2159     impl_->addOutputEnvironment(outputEnvironment);
2160     return *this;
2161 }
2162
2163 MdrunnerBuilder& MdrunnerBuilder::addLogFile(t_fileio* logFileHandle)
2164 {
2165     impl_->addLogFile(logFileHandle);
2166     return *this;
2167 }
2168
2169 MdrunnerBuilder& MdrunnerBuilder::addStopHandlerBuilder(std::unique_ptr<StopHandlerBuilder> builder)
2170 {
2171     impl_->addStopHandlerBuilder(std::move(builder));
2172     return *this;
2173 }
2174
2175 MdrunnerBuilder::MdrunnerBuilder(MdrunnerBuilder&&) noexcept = default;
2176
2177 MdrunnerBuilder& MdrunnerBuilder::operator=(MdrunnerBuilder&&) noexcept = default;
2178
2179 } // namespace gmx