src/gromacs/modularsimulator/simulatoralgorithm.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2020,2021, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 /*! \internal \file
  36  * \brief Defines the modular simulator algorithm
  37  *
  38  * \author Pascal Merz <pascal.merz@me.com>
  39  * \ingroup module_modularsimulator
  40  */
  41
  42 #include "gmxpre.h"
  43
  44 #include "simulatoralgorithm.h"
  45
  46 #include "gromacs/commandline/filenm.h"
  47 #include "gromacs/domdec/domdec.h"
  48 #include "gromacs/ewald/pme.h"
  49 #include "gromacs/ewald/pme_load_balancing.h"
  50 #include "gromacs/ewald/pme_pp.h"
  51 #include "gromacs/gmxlib/nrnb.h"
  52 #include "gromacs/listed_forces/listed_forces.h"
  53 #include "gromacs/mdlib/checkpointhandler.h"
  54 #include "gromacs/mdlib/constr.h"
  55 #include "gromacs/mdlib/energyoutput.h"
  56 #include "gromacs/mdlib/md_support.h"
  57 #include "gromacs/mdlib/mdatoms.h"
  58 #include "gromacs/mdlib/resethandler.h"
  59 #include "gromacs/mdlib/stat.h"
  60 #include "gromacs/mdrun/replicaexchange.h"
  61 #include "gromacs/mdrun/shellfc.h"
  62 #include "gromacs/mdrunutility/freeenergy.h"
  63 #include "gromacs/mdrunutility/handlerestart.h"
  64 #include "gromacs/mdrunutility/printtime.h"
  65 #include "gromacs/mdtypes/commrec.h"
  66 #include "gromacs/mdtypes/fcdata.h"
  67 #include "gromacs/mdtypes/forcerec.h"
  68 #include "gromacs/mdtypes/inputrec.h"
  69 #include "gromacs/mdtypes/mdatom.h"
  70 #include "gromacs/mdtypes/mdrunoptions.h"
  71 #include "gromacs/mdtypes/observableshistory.h"
  72 #include "gromacs/nbnxm/nbnxm.h"
  73 #include "gromacs/timing/walltime_accounting.h"
  74 #include "gromacs/topology/topology.h"
  75 #include "gromacs/utility/cstringutil.h"
  76 #include "gromacs/utility/fatalerror.h"
  77
  78 #include "checkpointhelper.h"
  79 #include "domdechelper.h"
  80 #include "energydata.h"
  81 #include "firstorderpressurecoupling.h"
  82 #include "freeenergyperturbationdata.h"
  83 #include "modularsimulator.h"
  84 #include "pmeloadbalancehelper.h"
  85 #include "propagator.h"
  86 #include "statepropagatordata.h"
  87
  88 namespace gmx
  89 {
  90 ModularSimulatorAlgorithm::ModularSimulatorAlgorithm(std::string              topologyName,
  91                                                      FILE*                    fplog,
  92                                                      t_commrec*               cr,
  93                                                      const MDLogger&          mdlog,
  94                                                      const MdrunOptions&      mdrunOptions,
  95                                                      const t_inputrec*        inputrec,
  96                                                      t_nrnb*                  nrnb,
  97                                                      gmx_wallcycle*           wcycle,
  98                                                      t_forcerec*              fr,
  99                                                      gmx_walltime_accounting* walltime_accounting) :
 100     taskIterator_(taskQueue_.end()),
 101     statePropagatorData_(nullptr),
 102     energyData_(nullptr),
 103     freeEnergyPerturbationData_(nullptr),
 104     step_(-1),
 105     runFinished_(false),
 106     topologyName_(std::move(topologyName)),
 107     fplog(fplog),
 108     cr(cr),
 109     mdlog(mdlog),
 110     mdrunOptions(mdrunOptions),
 111     inputrec(inputrec),
 112     nrnb(nrnb),
 113     wcycle(wcycle),
 114     fr(fr),
 115     walltime_accounting(walltime_accounting)
 116 {
 117     signalHelper_ = std::make_unique<SignalHelper>();
 118 }
 119
 120 void ModularSimulatorAlgorithm::setup()
 121 {
 122     simulatorSetup();
 123     for (auto& signaller : signallerList_)
 124     {
 125         signaller->setup();
 126     }
 127     if (domDecHelper_)
 128     {
 129         domDecHelper_->setup();
 130     }
 131
 132     for (auto& element : elementSetupTeardownList_)
 133     {
 134         element->elementSetup();
 135     }
 136     statePropagatorData_->setup();
 137     if (pmeLoadBalanceHelper_)
 138     {
 139         // State must have been initialized so pmeLoadBalanceHelper_ gets a valid box
 140         pmeLoadBalanceHelper_->setup();
 141     }
 142 }
 143
 144 const SimulatorRunFunction* ModularSimulatorAlgorithm::getNextTask()
 145 {
 146     if (!taskQueue_.empty())
 147     {
 148         taskIterator_++;
 149     }
 150     if (taskIterator_ == taskQueue_.end())
 151     {
 152         if (runFinished_)
 153         {
 154             return nullptr;
 155         }
 156         updateTaskQueue();
 157         taskIterator_ = taskQueue_.begin();
 158     }
 159     return &*taskIterator_;
 160 }
 161
 162 void ModularSimulatorAlgorithm::updateTaskQueue()
 163 {
 164     // For now, we'll just clean the task queue and then re-populate
 165     // TODO: If tasks are periodic around updates of the task queue,
 166     //       we should reuse it instead
 167     taskQueue_.clear();
 168     populateTaskQueue();
 169 }
 170
 171 void ModularSimulatorAlgorithm::teardown()
 172 {
 173     for (auto& element : elementSetupTeardownList_)
 174     {
 175         element->elementTeardown();
 176     }
 177     energyData_->teardown();
 178     if (pmeLoadBalanceHelper_)
 179     {
 180         pmeLoadBalanceHelper_->teardown();
 181     }
 182     simulatorTeardown();
 183 }
 184
 185 void ModularSimulatorAlgorithm::simulatorSetup()
 186 {
 187     if (!mdrunOptions.writeConfout)
 188     {
 189         // This is on by default, and the main known use case for
 190         // turning it off is for convenience in benchmarking, which is
 191         // something that should not show up in the general user
 192         // interface.
 193         GMX_LOG(mdlog.info)
 194                 .asParagraph()
 195                 .appendText(
 196                         "The -noconfout functionality is deprecated, and "
 197                         "may be removed in a future version.");
 198     }
 199
 200     if (MASTER(cr))
 201     {
 202         char        sbuf[STEPSTRSIZE], sbuf2[STEPSTRSIZE];
 203         std::string timeString;
 204         fprintf(stderr, "starting mdrun '%s'\n", topologyName_.c_str());
 205         if (inputrec->nsteps >= 0)
 206         {
 207             timeString = formatString(
 208                     "%8.1f", static_cast<double>(inputrec->init_step + inputrec->nsteps) * inputrec->delta_t);
 209         }
 210         else
 211         {
 212             timeString = "infinite";
 213         }
 214         if (inputrec->init_step > 0)
 215         {
 216             fprintf(stderr,
 217                     "%s steps, %s ps (continuing from step %s, %8.1f ps).\n",
 218                     gmx_step_str(inputrec->init_step + inputrec->nsteps, sbuf),
 219                     timeString.c_str(),
 220                     gmx_step_str(inputrec->init_step, sbuf2),
 221                     inputrec->init_step * inputrec->delta_t);
 222         }
 223         else
 224         {
 225             fprintf(stderr, "%s steps, %s ps.\n", gmx_step_str(inputrec->nsteps, sbuf), timeString.c_str());
 226         }
 227         fprintf(fplog, "\n");
 228     }
 229
 230     walltime_accounting_start_time(walltime_accounting);
 231     wallcycle_start(wcycle, WallCycleCounter::Run);
 232     print_start(fplog, cr, walltime_accounting, "mdrun");
 233
 234     step_ = inputrec->init_step;
 235 }
 236
 237 void ModularSimulatorAlgorithm::simulatorTeardown()
 238 {
 239
 240     // Stop measuring walltime
 241     walltime_accounting_end_time(walltime_accounting);
 242
 243     if (!thisRankHasDuty(cr, DUTY_PME))
 244     {
 245         /* Tell the PME only node to finish */
 246         gmx_pme_send_finish(cr);
 247     }
 248
 249     walltime_accounting_set_nsteps_done(walltime_accounting, step_ - inputrec->init_step);
 250 }
 251
 252 void ModularSimulatorAlgorithm::preStep(Step step, Time gmx_unused time, bool isNeighborSearchingStep)
 253 {
 254     if (stopHandler_->stoppingAfterCurrentStep(isNeighborSearchingStep) && step != signalHelper_->lastStep_)
 255     {
 256         /*
 257          * Stop handler wants to stop after the current step, which was
 258          * not known when building the current task queue. This happens
 259          * e.g. when a stop is signalled by OS. We therefore want to purge
 260          * the task queue now, and re-schedule this step as last step.
 261          */
 262         // clear task queue
 263         taskQueue_.clear();
 264         // rewind step
 265         step_ = step;
 266         return;
 267     }
 268
 269     resetHandler_->setSignal(walltime_accounting);
 270     // This is a hack to avoid having to rewrite StopHandler to be a NeighborSearchSignaller
 271     // and accept the step as input. Eventually, we want to do that, but currently this would
 272     // require introducing NeighborSearchSignaller in the legacy do_md or a lot of code
 273     // duplication.
 274     stophandlerIsNSStep_    = isNeighborSearchingStep;
 275     stophandlerCurrentStep_ = step;
 276     stopHandler_->setSignal();
 277
 278     wallcycle_start(wcycle, WallCycleCounter::Step);
 279 }
 280
 281 void ModularSimulatorAlgorithm::postStep(Step step, Time gmx_unused time)
 282 {
 283     // Output stuff
 284     if (MASTER(cr))
 285     {
 286         if (do_per_step(step, inputrec->nstlog))
 287         {
 288             if (fflush(fplog) != 0)
 289             {
 290                 gmx_fatal(FARGS, "Cannot flush logfile - maybe you are out of disk space?");
 291             }
 292         }
 293     }
 294     const bool do_verbose = mdrunOptions.verbose
 295                             && (step % mdrunOptions.verboseStepPrintInterval == 0
 296                                 || step == inputrec->init_step || step == signalHelper_->lastStep_);
 297     // Print the remaining wall clock time for the run
 298     if (MASTER(cr) && (do_verbose || gmx_got_usr_signal())
 299         && !(pmeLoadBalanceHelper_ && pmeLoadBalanceHelper_->pmePrinting()))
 300     {
 301         print_time(stderr, walltime_accounting, step, inputrec, cr);
 302     }
 303
 304     double cycles = wallcycle_stop(wcycle, WallCycleCounter::Step);
 305     if (DOMAINDECOMP(cr) && wcycle)
 306     {
 307         dd_cycles_add(cr->dd, static_cast<float>(cycles), ddCyclStep);
 308     }
 309
 310     resetHandler_->resetCounters(step,
 311                                  step - inputrec->init_step,
 312                                  mdlog,
 313                                  fplog,
 314                                  cr,
 315                                  fr->nbv.get(),
 316                                  nrnb,
 317                                  fr->pmedata,
 318                                  pmeLoadBalanceHelper_ ? pmeLoadBalanceHelper_->loadBalancingObject() : nullptr,
 319                                  wcycle,
 320                                  walltime_accounting);
 321 }
 322
 323 void ModularSimulatorAlgorithm::populateTaskQueue()
 324 {
 325     /*
 326      * The registerRunFunction emplaces functions to the task queue.
 327      * All elements are owned by the ModularSimulatorAlgorithm, as is the task queue.
 328      * Elements can hence register lambdas capturing their `this` pointers without expecting
 329      * life time issues, as the task queue and the elements are in the same scope.
 330      */
 331     auto registerRunFunction = [this](SimulatorRunFunction function) {
 332         taskQueue_.emplace_back(std::move(function));
 333     };
 334
 335     Time startTime = inputrec->init_t;
 336     Time timeStep  = inputrec->delta_t;
 337     Time time      = startTime + step_ * timeStep;
 338
 339     // Run an initial call to the signallers
 340     for (auto& signaller : signallerList_)
 341     {
 342         signaller->signal(step_, time);
 343     }
 344
 345     if (checkpointHelper_)
 346     {
 347         checkpointHelper_->run(step_, time);
 348     }
 349
 350     if (pmeLoadBalanceHelper_)
 351     {
 352         pmeLoadBalanceHelper_->run(step_, time);
 353     }
 354     if (domDecHelper_)
 355     {
 356         domDecHelper_->run(step_, time);
 357     }
 358
 359     do
 360     {
 361         // local variables for lambda capturing
 362         const int  step     = step_;
 363         const bool isNSStep = step == signalHelper_->nextNSStep_;
 364
 365         // register pre-step (task queue is local, so no problem with `this`)
 366         registerRunFunction([this, step, time, isNSStep]() { preStep(step, time, isNSStep); });
 367         // register elements for step
 368         for (auto& element : elementCallList_)
 369         {
 370             element->scheduleTask(step_, time, registerRunFunction);
 371         }
 372         // register post-step (task queue is local, so no problem with `this`)
 373         registerRunFunction([this, step, time]() { postStep(step, time); });
 374
 375         // prepare next step
 376         step_++;
 377         time = startTime + step_ * timeStep;
 378         for (auto& signaller : signallerList_)
 379         {
 380             signaller->signal(step_, time);
 381         }
 382     } while (step_ != signalHelper_->nextNSStep_ && step_ <= signalHelper_->lastStep_);
 383
 384     runFinished_ = (step_ > signalHelper_->lastStep_);
 385
 386     if (runFinished_)
 387     {
 388         // task queue is local, so no problem with `this`
 389         registerRunFunction([this]() { teardown(); });
 390     }
 391 }
 392
 393 ModularSimulatorAlgorithmBuilder::ModularSimulatorAlgorithmBuilder(
 394         compat::not_null<LegacySimulatorData*>    legacySimulatorData,
 395         std::unique_ptr<ReadCheckpointDataHolder> checkpointDataHolder) :
 396     legacySimulatorData_(legacySimulatorData),
 397     signals_(std::make_unique<SimulationSignals>()),
 398     elementAdditionHelper_(this),
 399     globalCommunicationHelper_(computeGlobalCommunicationPeriod(legacySimulatorData->mdlog,
 400                                                                 legacySimulatorData->inputrec,
 401                                                                 legacySimulatorData->cr),
 402                                signals_.get()),
 403     checkpointHelperBuilder_(std::move(checkpointDataHolder),
 404                              legacySimulatorData->startingBehavior,
 405                              legacySimulatorData->cr)
 406 {
 407     if (legacySimulatorData->inputrec->efep != FreeEnergyPerturbationType::No)
 408     {
 409         freeEnergyPerturbationData_ = std::make_unique<FreeEnergyPerturbationData>(
 410                 legacySimulatorData->fplog, *legacySimulatorData->inputrec, legacySimulatorData->mdAtoms);
 411         registerExistingElement(freeEnergyPerturbationData_->element());
 412     }
 413
 414     statePropagatorData_ = std::make_unique<StatePropagatorData>(
 415             legacySimulatorData->top_global.natoms,
 416             legacySimulatorData->fplog,
 417             legacySimulatorData->cr,
 418             legacySimulatorData->state_global,
 419             legacySimulatorData->fr->nbv->useGpu(),
 420             legacySimulatorData->fr->bMolPBC,
 421             legacySimulatorData->mdrunOptions.writeConfout,
 422             opt2fn("-c", legacySimulatorData->nfile, legacySimulatorData->fnm),
 423             legacySimulatorData->inputrec,
 424             legacySimulatorData->mdAtoms->mdatoms(),
 425             legacySimulatorData->top_global);
 426     registerExistingElement(statePropagatorData_->element());
 427
 428     // Multi sim is turned off
 429     const bool simulationsShareState = false;
 430
 431     energyData_ = std::make_unique<EnergyData>(statePropagatorData_.get(),
 432                                                freeEnergyPerturbationData_.get(),
 433                                                legacySimulatorData->top_global,
 434                                                legacySimulatorData->inputrec,
 435                                                legacySimulatorData->mdAtoms,
 436                                                legacySimulatorData->enerd,
 437                                                legacySimulatorData->ekind,
 438                                                legacySimulatorData->constr,
 439                                                legacySimulatorData->fplog,
 440                                                legacySimulatorData->fr->fcdata.get(),
 441                                                legacySimulatorData->mdModulesNotifiers,
 442                                                MASTER(legacySimulatorData->cr),
 443                                                legacySimulatorData->observablesHistory,
 444                                                legacySimulatorData->startingBehavior,
 445                                                simulationsShareState);
 446     registerExistingElement(energyData_->element());
 447 }
 448
 449 ModularSimulatorAlgorithm ModularSimulatorAlgorithmBuilder::build()
 450 {
 451     if (algorithmHasBeenBuilt_)
 452     {
 453         GMX_THROW(SimulationAlgorithmSetupError(
 454                 "Tried to build ModularSimulationAlgorithm more than once."));
 455     }
 456     algorithmHasBeenBuilt_ = true;
 457
 458     // Connect propagators with thermostat / barostat
 459     for (const auto& registrationFunction : pressureTemperatureControlRegistrationFunctions_)
 460     {
 461         for (const auto& connection : propagatorConnections_)
 462         {
 463             registrationFunction(connection);
 464         }
 465     }
 466
 467     ModularSimulatorAlgorithm algorithm(*(legacySimulatorData_->top_global.name),
 468                                         legacySimulatorData_->fplog,
 469                                         legacySimulatorData_->cr,
 470                                         legacySimulatorData_->mdlog,
 471                                         legacySimulatorData_->mdrunOptions,
 472                                         legacySimulatorData_->inputrec,
 473                                         legacySimulatorData_->nrnb,
 474                                         legacySimulatorData_->wcycle,
 475                                         legacySimulatorData_->fr,
 476                                         legacySimulatorData_->walltime_accounting);
 477     registerWithInfrastructureAndSignallers(algorithm.signalHelper_.get());
 478     algorithm.statePropagatorData_        = std::move(statePropagatorData_);
 479     algorithm.energyData_                 = std::move(energyData_);
 480     algorithm.freeEnergyPerturbationData_ = std::move(freeEnergyPerturbationData_);
 481     algorithm.signals_                    = std::move(signals_);
 482     algorithm.simulationData_             = std::move(simulationData_);
 483
 484     // Multi sim is turned off
 485     const bool simulationsShareState = false;
 486
 487     // Build stop handler
 488     algorithm.stopHandler_ = legacySimulatorData_->stopHandlerBuilder->getStopHandlerMD(
 489             compat::not_null<SimulationSignal*>(
 490                     &(*globalCommunicationHelper_.simulationSignals())[eglsSTOPCOND]),
 491             simulationsShareState,
 492             MASTER(legacySimulatorData_->cr),
 493             legacySimulatorData_->inputrec->nstlist,
 494             legacySimulatorData_->mdrunOptions.reproducible,
 495             globalCommunicationHelper_.nstglobalcomm(),
 496             legacySimulatorData_->mdrunOptions.maximumHoursToRun,
 497             legacySimulatorData_->inputrec->nstlist == 0,
 498             legacySimulatorData_->fplog,
 499             algorithm.stophandlerCurrentStep_,
 500             algorithm.stophandlerIsNSStep_,
 501             legacySimulatorData_->walltime_accounting);
 502
 503     // Build reset handler
 504     const bool simulationsShareResetCounters = false;
 505     algorithm.resetHandler_                  = std::make_unique<ResetHandler>(
 506             compat::make_not_null<SimulationSignal*>(
 507                     &(*globalCommunicationHelper_.simulationSignals())[eglsRESETCOUNTERS]),
 508             simulationsShareResetCounters,
 509             legacySimulatorData_->inputrec->nsteps,
 510             MASTER(legacySimulatorData_->cr),
 511             legacySimulatorData_->mdrunOptions.timingOptions.resetHalfway,
 512             legacySimulatorData_->mdrunOptions.maximumHoursToRun,
 513             legacySimulatorData_->mdlog,
 514             legacySimulatorData_->wcycle,
 515             legacySimulatorData_->walltime_accounting);
 516
 517     // Build topology holder
 518     algorithm.topologyHolder_ = topologyHolderBuilder_.build(legacySimulatorData_->top_global,
 519                                                              legacySimulatorData_->cr,
 520                                                              legacySimulatorData_->inputrec,
 521                                                              legacySimulatorData_->fr,
 522                                                              legacySimulatorData_->mdAtoms,
 523                                                              legacySimulatorData_->constr,
 524                                                              legacySimulatorData_->vsite);
 525     registerWithInfrastructureAndSignallers(algorithm.topologyHolder_.get());
 526
 527     // Build PME load balance helper
 528     if (PmeLoadBalanceHelper::doPmeLoadBalancing(legacySimulatorData_->mdrunOptions,
 529                                                  legacySimulatorData_->inputrec,
 530                                                  legacySimulatorData_->fr))
 531     {
 532         algorithm.pmeLoadBalanceHelper_ =
 533                 std::make_unique<PmeLoadBalanceHelper>(legacySimulatorData_->mdrunOptions.verbose,
 534                                                        algorithm.statePropagatorData_.get(),
 535                                                        legacySimulatorData_->fplog,
 536                                                        legacySimulatorData_->cr,
 537                                                        legacySimulatorData_->mdlog,
 538                                                        legacySimulatorData_->inputrec,
 539                                                        legacySimulatorData_->wcycle,
 540                                                        legacySimulatorData_->fr);
 541         registerWithInfrastructureAndSignallers(algorithm.pmeLoadBalanceHelper_.get());
 542     }
 543
 544     // Build trajectory element
 545     auto trajectoryElement = trajectoryElementBuilder_.build(legacySimulatorData_->fplog,
 546                                                              legacySimulatorData_->nfile,
 547                                                              legacySimulatorData_->fnm,
 548                                                              legacySimulatorData_->mdrunOptions,
 549                                                              legacySimulatorData_->cr,
 550                                                              legacySimulatorData_->outputProvider,
 551                                                              legacySimulatorData_->mdModulesNotifiers,
 552                                                              legacySimulatorData_->inputrec,
 553                                                              legacySimulatorData_->top_global,
 554                                                              legacySimulatorData_->oenv,
 555                                                              legacySimulatorData_->wcycle,
 556                                                              legacySimulatorData_->startingBehavior,
 557                                                              simulationsShareState);
 558     registerWithInfrastructureAndSignallers(trajectoryElement.get());
 559
 560     // Build domdec helper
 561     if (DOMAINDECOMP(legacySimulatorData_->cr))
 562     {
 563         algorithm.domDecHelper_ =
 564                 domDecHelperBuilder_.build(legacySimulatorData_->mdrunOptions.verbose,
 565                                            legacySimulatorData_->mdrunOptions.verboseStepPrintInterval,
 566                                            algorithm.statePropagatorData_.get(),
 567                                            algorithm.topologyHolder_.get(),
 568                                            globalCommunicationHelper_.nstglobalcomm(),
 569                                            legacySimulatorData_->fplog,
 570                                            legacySimulatorData_->cr,
 571                                            legacySimulatorData_->mdlog,
 572                                            legacySimulatorData_->constr,
 573                                            legacySimulatorData_->inputrec,
 574                                            legacySimulatorData_->mdAtoms,
 575                                            legacySimulatorData_->nrnb,
 576                                            legacySimulatorData_->wcycle,
 577                                            legacySimulatorData_->fr,
 578                                            legacySimulatorData_->vsite,
 579                                            legacySimulatorData_->imdSession,
 580                                            legacySimulatorData_->pull_work);
 581         registerWithInfrastructureAndSignallers(algorithm.domDecHelper_.get());
 582     }
 583     // Build checkpoint helper (do this last so everyone else can be a checkpoint client!)
 584     {
 585         checkpointHelperBuilder_.setCheckpointHandler(std::make_unique<CheckpointHandler>(
 586                 compat::make_not_null<SimulationSignal*>(&(*algorithm.signals_)[eglsCHKPT]),
 587                 simulationsShareState,
 588                 legacySimulatorData_->inputrec->nstlist == 0,
 589                 MASTER(legacySimulatorData_->cr),
 590                 legacySimulatorData_->mdrunOptions.writeConfout,
 591                 legacySimulatorData_->mdrunOptions.checkpointOptions.period));
 592         algorithm.checkpointHelper_ =
 593                 checkpointHelperBuilder_.build(legacySimulatorData_->inputrec->init_step,
 594                                                trajectoryElement.get(),
 595                                                legacySimulatorData_->fplog,
 596                                                legacySimulatorData_->cr,
 597                                                legacySimulatorData_->observablesHistory,
 598                                                legacySimulatorData_->walltime_accounting,
 599                                                legacySimulatorData_->state_global,
 600                                                legacySimulatorData_->mdrunOptions.writeConfout);
 601         registerWithInfrastructureAndSignallers(algorithm.checkpointHelper_.get());
 602     }
 603
 604     // Build signallers
 605     {
 606         /* Signallers need to be called in an exact order. Some signallers are clients
 607          * of other signallers, which requires the clients signallers to be called
 608          * _after_ any signaller they are registered to - otherwise, they couldn't
 609          * adapt their behavior to the information they got signalled.
 610          *
 611          * Signallers being clients of other signallers require registration.
 612          * That registration happens during construction, which in turn means that
 613          * we want to construct the signallers in the reverse order of their later
 614          * call order.
 615          *
 616          * For the above reasons, the `addSignaller` lambda defined below emplaces
 617          * added signallers at the beginning of the signaller list, which will yield
 618          * a signaller list which is inverse to the build order (and hence equal to
 619          * the intended call order).
 620          */
 621         auto addSignaller = [this, &algorithm](auto signaller) {
 622             registerWithInfrastructureAndSignallers(signaller.get());
 623             algorithm.signallerList_.emplace(algorithm.signallerList_.begin(), std::move(signaller));
 624         };
 625         const auto* inputrec   = legacySimulatorData_->inputrec;
 626         auto        virialMode = EnergySignallerVirialMode::Off;
 627         if (inputrec->epc != PressureCoupling::No)
 628         {
 629             if (EI_VV(inputrec->eI))
 630             {
 631                 virialMode = EnergySignallerVirialMode::OnStepAndNext;
 632             }
 633             else
 634             {
 635                 virialMode = EnergySignallerVirialMode::OnStep;
 636             }
 637         }
 638         addSignaller(energySignallerBuilder_.build(
 639                 inputrec->nstcalcenergy,
 640                 computeFepPeriod(*inputrec, legacySimulatorData_->replExParams),
 641                 inputrec->nstpcouple,
 642                 virialMode));
 643         addSignaller(trajectorySignallerBuilder_.build(inputrec->nstxout,
 644                                                        inputrec->nstvout,
 645                                                        inputrec->nstfout,
 646                                                        inputrec->nstxout_compressed,
 647                                                        trajectoryElement->tngBoxOut(),
 648                                                        trajectoryElement->tngLambdaOut(),
 649                                                        trajectoryElement->tngBoxOutCompressed(),
 650                                                        trajectoryElement->tngLambdaOutCompressed(),
 651                                                        inputrec->nstenergy));
 652         addSignaller(loggingSignallerBuilder_.build(
 653                 inputrec->nstlog, inputrec->init_step, legacySimulatorData_->startingBehavior));
 654         addSignaller(lastStepSignallerBuilder_.build(
 655                 inputrec->nsteps, inputrec->init_step, algorithm.stopHandler_.get()));
 656         addSignaller(neighborSearchSignallerBuilder_.build(
 657                 inputrec->nstlist, inputrec->init_step, inputrec->init_t));
 658     }
 659
 660     // Move setup / teardown list
 661     algorithm.elementSetupTeardownList_ = std::move(setupAndTeardownList_);
 662
 663     // Create element list
 664     // Checkpoint helper needs to be in the call list (as first element!) to react to last step
 665     algorithm.elementCallList_.emplace_back(algorithm.checkpointHelper_.get());
 666     // Next, update the free energy lambda vector if needed
 667     if (algorithm.freeEnergyPerturbationData_)
 668     {
 669         algorithm.elementCallList_.emplace_back(algorithm.freeEnergyPerturbationData_->element());
 670     }
 671     // Then, move the built algorithm
 672     algorithm.elementsOwnershipList_.insert(algorithm.elementsOwnershipList_.end(),
 673                                             std::make_move_iterator(elements_.begin()),
 674                                             std::make_move_iterator(elements_.end()));
 675     algorithm.elementCallList_.insert(algorithm.elementCallList_.end(),
 676                                       std::make_move_iterator(callList_.begin()),
 677                                       std::make_move_iterator(callList_.end()));
 678     // Finally, all trajectory writing is happening after the step
 679     // (relevant data was stored by elements through energy signaller)
 680     algorithm.elementsOwnershipList_.emplace_back(std::move(trajectoryElement));
 681     algorithm.elementCallList_.emplace_back(algorithm.elementsOwnershipList_.back().get());
 682     algorithm.elementSetupTeardownList_.emplace_back(algorithm.elementsOwnershipList_.back().get());
 683
 684     algorithm.setup();
 685     return algorithm;
 686 }
 687
 688 bool ModularSimulatorAlgorithmBuilder::elementExists(const ISimulatorElement* element) const
 689 {
 690     // Check whether element exists in element list
 691     if (std::any_of(elements_.begin(), elements_.end(), [element](auto& existingElement) {
 692             return element == existingElement.get();
 693         }))
 694     {
 695         return true;
 696     }
 697     // Check whether element exists in other places controlled by *this
 698     return ((statePropagatorData_ && statePropagatorData_->element() == element)
 699             || (energyData_ && energyData_->element() == element)
 700             || (freeEnergyPerturbationData_ && freeEnergyPerturbationData_->element() == element));
 701 }
 702
 703 std::optional<SignallerCallback> ModularSimulatorAlgorithm::SignalHelper::registerLastStepCallback()
 704 {
 705     return [this](Step step, Time gmx_unused time) { this->lastStep_ = step; };
 706 }
 707
 708 std::optional<SignallerCallback> ModularSimulatorAlgorithm::SignalHelper::registerNSCallback()
 709 {
 710     return [this](Step step, Time gmx_unused time) { this->nextNSStep_ = step; };
 711 }
 712
 713 GlobalCommunicationHelper::GlobalCommunicationHelper(int nstglobalcomm, SimulationSignals* simulationSignals) :
 714     nstglobalcomm_(nstglobalcomm), simulationSignals_(simulationSignals)
 715 {
 716 }
 717
 718 int GlobalCommunicationHelper::nstglobalcomm() const
 719 {
 720     return nstglobalcomm_;
 721 }
 722
 723 SimulationSignals* GlobalCommunicationHelper::simulationSignals()
 724 {
 725     return simulationSignals_;
 726 }
 727
 728 ModularSimulatorAlgorithmBuilderHelper::ModularSimulatorAlgorithmBuilderHelper(
 729         ModularSimulatorAlgorithmBuilder* builder) :
 730     builder_(builder)
 731 {
 732 }
 733
 734 bool ModularSimulatorAlgorithmBuilderHelper::elementIsStored(const ISimulatorElement* element) const
 735 {
 736     return builder_->elementExists(element);
 737 }
 738
 739 std::optional<std::any> ModularSimulatorAlgorithmBuilderHelper::builderData(const std::string& key) const
 740 {
 741     const auto iter = builder_->builderData_.find(key);
 742     if (iter == builder_->builderData_.end())
 743     {
 744         return std::nullopt;
 745     }
 746     else
 747     {
 748         return iter->second;
 749     }
 750 }
 751
 752 void ModularSimulatorAlgorithmBuilderHelper::registerTemperaturePressureControl(
 753         std::function<void(const PropagatorConnection&)> registrationFunction)
 754 {
 755     builder_->pressureTemperatureControlRegistrationFunctions_.emplace_back(std::move(registrationFunction));
 756 }
 757
 758 void ModularSimulatorAlgorithmBuilderHelper::registerPropagator(PropagatorConnection connectionData)
 759 {
 760     builder_->propagatorConnections_.emplace_back(std::move(connectionData));
 761 }
 762
 763 } // namespace gmx