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/handlerestart.h"
63 #include "gromacs/mdrunutility/printtime.h"
64 #include "gromacs/mdtypes/commrec.h"
65 #include "gromacs/mdtypes/fcdata.h"
66 #include "gromacs/mdtypes/forcerec.h"
67 #include "gromacs/mdtypes/inputrec.h"
68 #include "gromacs/mdtypes/mdatom.h"
69 #include "gromacs/mdtypes/mdrunoptions.h"
70 #include "gromacs/mdtypes/observableshistory.h"
71 #include "gromacs/nbnxm/nbnxm.h"
72 #include "gromacs/timing/walltime_accounting.h"
73 #include "gromacs/topology/topology.h"
74 #include "gromacs/utility/cstringutil.h"
75 #include "gromacs/utility/fatalerror.h"
76
77 #include "checkpointhelper.h"
78 #include "domdechelper.h"
79 #include "energydata.h"
80 #include "freeenergyperturbationdata.h"
81 #include "modularsimulator.h"
82 #include "parrinellorahmanbarostat.h"
83 #include "pmeloadbalancehelper.h"
84 #include "propagator.h"
85 #include "statepropagatordata.h"
86 #include "velocityscalingtemperaturecoupling.h"
87
88 namespace gmx
89 {
ModularSimulatorAlgorithm(std::string topologyName,FILE * fplog,t_commrec * cr,const MDLogger & mdlog,const MdrunOptions & mdrunOptions,t_inputrec * inputrec,t_nrnb * nrnb,gmx_wallcycle * wcycle,t_forcerec * fr,gmx_walltime_accounting * walltime_accounting)90 ModularSimulatorAlgorithm::ModularSimulatorAlgorithm(std::string topologyName,
91 FILE* fplog,
92 t_commrec* cr,
93 const MDLogger& mdlog,
94 const MdrunOptions& mdrunOptions,
95 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
setup()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 : elementsOwnershipList_)
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
getNextTask()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
updateTaskQueue()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
teardown()171 void ModularSimulatorAlgorithm::teardown()
172 {
173 for (auto& element : elementsOwnershipList_)
174 {
175 element->elementTeardown();
176 }
177 energyData_->teardown();
178 if (pmeLoadBalanceHelper_)
179 {
180 pmeLoadBalanceHelper_->teardown();
181 }
182 simulatorTeardown();
183 }
184
simulatorSetup()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("%8.1f", static_cast<double>(inputrec->init_step + inputrec->nsteps)
208 * inputrec->delta_t);
209 }
210 else
211 {
212 timeString = "infinite";
213 }
214 if (inputrec->init_step > 0)
215 {
216 fprintf(stderr, "%s steps, %s ps (continuing from step %s, %8.1f ps).\n",
217 gmx_step_str(inputrec->init_step + inputrec->nsteps, sbuf), timeString.c_str(),
218 gmx_step_str(inputrec->init_step, sbuf2), inputrec->init_step * inputrec->delta_t);
219 }
220 else
221 {
222 fprintf(stderr, "%s steps, %s ps.\n", gmx_step_str(inputrec->nsteps, sbuf),
223 timeString.c_str());
224 }
225 fprintf(fplog, "\n");
226 }
227
228 walltime_accounting_start_time(walltime_accounting);
229 wallcycle_start(wcycle, ewcRUN);
230 print_start(fplog, cr, walltime_accounting, "mdrun");
231
232 step_ = inputrec->init_step;
233 }
234
simulatorTeardown()235 void ModularSimulatorAlgorithm::simulatorTeardown()
236 {
237
238 // Stop measuring walltime
239 walltime_accounting_end_time(walltime_accounting);
240
241 if (!thisRankHasDuty(cr, DUTY_PME))
242 {
243 /* Tell the PME only node to finish */
244 gmx_pme_send_finish(cr);
245 }
246
247 walltime_accounting_set_nsteps_done(walltime_accounting, step_ - inputrec->init_step);
248 }
249
preStep(Step step,Time gmx_unused time,bool isNeighborSearchingStep)250 void ModularSimulatorAlgorithm::preStep(Step step, Time gmx_unused time, bool isNeighborSearchingStep)
251 {
252 if (stopHandler_->stoppingAfterCurrentStep(isNeighborSearchingStep) && step != signalHelper_->lastStep_)
253 {
254 /*
255 * Stop handler wants to stop after the current step, which was
256 * not known when building the current task queue. This happens
257 * e.g. when a stop is signalled by OS. We therefore want to purge
258 * the task queue now, and re-schedule this step as last step.
259 */
260 // clear task queue
261 taskQueue_.clear();
262 // rewind step
263 step_ = step;
264 return;
265 }
266
267 resetHandler_->setSignal(walltime_accounting);
268 // This is a hack to avoid having to rewrite StopHandler to be a NeighborSearchSignaller
269 // and accept the step as input. Eventually, we want to do that, but currently this would
270 // require introducing NeighborSearchSignaller in the legacy do_md or a lot of code
271 // duplication.
272 stophandlerIsNSStep_ = isNeighborSearchingStep;
273 stophandlerCurrentStep_ = step;
274 stopHandler_->setSignal();
275
276 wallcycle_start(wcycle, ewcSTEP);
277 }
278
postStep(Step step,Time gmx_unused time)279 void ModularSimulatorAlgorithm::postStep(Step step, Time gmx_unused time)
280 {
281 // Output stuff
282 if (MASTER(cr))
283 {
284 if (do_per_step(step, inputrec->nstlog))
285 {
286 if (fflush(fplog) != 0)
287 {
288 gmx_fatal(FARGS, "Cannot flush logfile - maybe you are out of disk space?");
289 }
290 }
291 }
292 const bool do_verbose = mdrunOptions.verbose
293 && (step % mdrunOptions.verboseStepPrintInterval == 0
294 || step == inputrec->init_step || step == signalHelper_->lastStep_);
295 // Print the remaining wall clock time for the run
296 if (MASTER(cr) && (do_verbose || gmx_got_usr_signal())
297 && !(pmeLoadBalanceHelper_ && pmeLoadBalanceHelper_->pmePrinting()))
298 {
299 print_time(stderr, walltime_accounting, step, inputrec, cr);
300 }
301
302 double cycles = wallcycle_stop(wcycle, ewcSTEP);
303 if (DOMAINDECOMP(cr) && wcycle)
304 {
305 dd_cycles_add(cr->dd, static_cast<float>(cycles), ddCyclStep);
306 }
307
308 resetHandler_->resetCounters(
309 step, step - inputrec->init_step, mdlog, fplog, cr, fr->nbv.get(), nrnb, fr->pmedata,
310 pmeLoadBalanceHelper_ ? pmeLoadBalanceHelper_->loadBalancingObject() : nullptr, wcycle,
311 walltime_accounting);
312 }
313
populateTaskQueue()314 void ModularSimulatorAlgorithm::populateTaskQueue()
315 {
316 /*
317 * The registerRunFunction emplaces functions to the task queue.
318 * All elements are owned by the ModularSimulatorAlgorithm, as is the task queue.
319 * Elements can hence register lambdas capturing their `this` pointers without expecting
320 * life time issues, as the task queue and the elements are in the same scope.
321 */
322 auto registerRunFunction = [this](SimulatorRunFunction function) {
323 taskQueue_.emplace_back(std::move(function));
324 };
325
326 Time startTime = inputrec->init_t;
327 Time timeStep = inputrec->delta_t;
328 Time time = startTime + step_ * timeStep;
329
330 // Run an initial call to the signallers
331 for (auto& signaller : signallerList_)
332 {
333 signaller->signal(step_, time);
334 }
335
336 if (checkpointHelper_)
337 {
338 checkpointHelper_->run(step_, time);
339 }
340
341 if (pmeLoadBalanceHelper_)
342 {
343 pmeLoadBalanceHelper_->run(step_, time);
344 }
345 if (domDecHelper_)
346 {
347 domDecHelper_->run(step_, time);
348 }
349
350 do
351 {
352 // local variables for lambda capturing
353 const int step = step_;
354 const bool isNSStep = step == signalHelper_->nextNSStep_;
355
356 // register pre-step (task queue is local, so no problem with `this`)
357 registerRunFunction([this, step, time, isNSStep]() { preStep(step, time, isNSStep); });
358 // register elements for step
359 for (auto& element : elementCallList_)
360 {
361 element->scheduleTask(step_, time, registerRunFunction);
362 }
363 // register post-step (task queue is local, so no problem with `this`)
364 registerRunFunction([this, step, time]() { postStep(step, time); });
365
366 // prepare next step
367 step_++;
368 time = startTime + step_ * timeStep;
369 for (auto& signaller : signallerList_)
370 {
371 signaller->signal(step_, time);
372 }
373 } while (step_ != signalHelper_->nextNSStep_ && step_ <= signalHelper_->lastStep_);
374
375 runFinished_ = (step_ > signalHelper_->lastStep_);
376
377 if (runFinished_)
378 {
379 // task queue is local, so no problem with `this`
380 registerRunFunction([this]() { teardown(); });
381 }
382 }
383
ModularSimulatorAlgorithmBuilder(compat::not_null<LegacySimulatorData * > legacySimulatorData,std::unique_ptr<ReadCheckpointDataHolder> checkpointDataHolder)384 ModularSimulatorAlgorithmBuilder::ModularSimulatorAlgorithmBuilder(
385 compat::not_null<LegacySimulatorData*> legacySimulatorData,
386 std::unique_ptr<ReadCheckpointDataHolder> checkpointDataHolder) :
387 legacySimulatorData_(legacySimulatorData),
388 signals_(std::make_unique<SimulationSignals>()),
389 elementAdditionHelper_(this),
390 globalCommunicationHelper_(computeGlobalCommunicationPeriod(legacySimulatorData->mdlog,
391 legacySimulatorData->inputrec,
392 legacySimulatorData->cr),
393 signals_.get()),
394 checkpointHelperBuilder_(std::move(checkpointDataHolder),
395 legacySimulatorData->startingBehavior,
396 legacySimulatorData->cr)
397 {
398 if (legacySimulatorData->inputrec->efep != efepNO)
399 {
400 freeEnergyPerturbationData_ = std::make_unique<FreeEnergyPerturbationData>(
401 legacySimulatorData->fplog, legacySimulatorData->inputrec, legacySimulatorData->mdAtoms);
402 }
403
404 statePropagatorData_ = std::make_unique<StatePropagatorData>(
405 legacySimulatorData->top_global->natoms, legacySimulatorData->fplog, legacySimulatorData->cr,
406 legacySimulatorData->state_global, legacySimulatorData->fr->nbv->useGpu(),
407 legacySimulatorData->fr->bMolPBC, legacySimulatorData->mdrunOptions.writeConfout,
408 opt2fn("-c", legacySimulatorData->nfile, legacySimulatorData->fnm), legacySimulatorData->inputrec,
409 legacySimulatorData->mdAtoms->mdatoms(), legacySimulatorData->top_global);
410
411 // Multi sim is turned off
412 const bool simulationsShareState = false;
413
414 energyData_ = std::make_unique<EnergyData>(
415 statePropagatorData_.get(), freeEnergyPerturbationData_.get(), legacySimulatorData->top_global,
416 legacySimulatorData->inputrec, legacySimulatorData->mdAtoms, legacySimulatorData->enerd,
417 legacySimulatorData->ekind, legacySimulatorData->constr, legacySimulatorData->fplog,
418 legacySimulatorData->fr->fcdata.get(), legacySimulatorData->mdModulesNotifier,
419 MASTER(legacySimulatorData->cr), legacySimulatorData->observablesHistory,
420 legacySimulatorData->startingBehavior, simulationsShareState);
421 }
422
build()423 ModularSimulatorAlgorithm ModularSimulatorAlgorithmBuilder::build()
424 {
425 if (algorithmHasBeenBuilt_)
426 {
427 throw SimulationAlgorithmSetupError(
428 "Tried to build ModularSimulationAlgorithm more than once.");
429 }
430 algorithmHasBeenBuilt_ = true;
431
432 // Connect propagators with thermostat / barostat
433 for (const auto& thermostatRegistration : thermostatRegistrationFunctions_)
434 {
435 for (const auto& connection : propagatorThermostatConnections_)
436 {
437 thermostatRegistration(connection);
438 }
439 }
440 for (const auto& barostatRegistration : barostatRegistrationFunctions_)
441 {
442 for (const auto& connection : propagatorBarostatConnections_)
443 {
444 barostatRegistration(connection);
445 }
446 }
447
448 ModularSimulatorAlgorithm algorithm(
449 *(legacySimulatorData_->top_global->name), legacySimulatorData_->fplog,
450 legacySimulatorData_->cr, legacySimulatorData_->mdlog, legacySimulatorData_->mdrunOptions,
451 legacySimulatorData_->inputrec, legacySimulatorData_->nrnb, legacySimulatorData_->wcycle,
452 legacySimulatorData_->fr, legacySimulatorData_->walltime_accounting);
453 registerWithInfrastructureAndSignallers(algorithm.signalHelper_.get());
454 algorithm.statePropagatorData_ = std::move(statePropagatorData_);
455 algorithm.energyData_ = std::move(energyData_);
456 algorithm.freeEnergyPerturbationData_ = std::move(freeEnergyPerturbationData_);
457 algorithm.signals_ = std::move(signals_);
458
459 // Multi sim is turned off
460 const bool simulationsShareState = false;
461
462 // Build stop handler
463 algorithm.stopHandler_ = legacySimulatorData_->stopHandlerBuilder->getStopHandlerMD(
464 compat::not_null<SimulationSignal*>(
465 &(*globalCommunicationHelper_.simulationSignals())[eglsSTOPCOND]),
466 simulationsShareState, MASTER(legacySimulatorData_->cr),
467 legacySimulatorData_->inputrec->nstlist, legacySimulatorData_->mdrunOptions.reproducible,
468 globalCommunicationHelper_.nstglobalcomm(), legacySimulatorData_->mdrunOptions.maximumHoursToRun,
469 legacySimulatorData_->inputrec->nstlist == 0, legacySimulatorData_->fplog,
470 algorithm.stophandlerCurrentStep_, algorithm.stophandlerIsNSStep_,
471 legacySimulatorData_->walltime_accounting);
472
473 // Build reset handler
474 const bool simulationsShareResetCounters = false;
475 algorithm.resetHandler_ = std::make_unique<ResetHandler>(
476 compat::make_not_null<SimulationSignal*>(
477 &(*globalCommunicationHelper_.simulationSignals())[eglsRESETCOUNTERS]),
478 simulationsShareResetCounters, legacySimulatorData_->inputrec->nsteps,
479 MASTER(legacySimulatorData_->cr), legacySimulatorData_->mdrunOptions.timingOptions.resetHalfway,
480 legacySimulatorData_->mdrunOptions.maximumHoursToRun, legacySimulatorData_->mdlog,
481 legacySimulatorData_->wcycle, legacySimulatorData_->walltime_accounting);
482
483 // Build topology holder
484 algorithm.topologyHolder_ = topologyHolderBuilder_.build(
485 *legacySimulatorData_->top_global, legacySimulatorData_->cr,
486 legacySimulatorData_->inputrec, legacySimulatorData_->fr, legacySimulatorData_->mdAtoms,
487 legacySimulatorData_->constr, legacySimulatorData_->vsite);
488
489 // Build PME load balance helper
490 if (PmeLoadBalanceHelper::doPmeLoadBalancing(legacySimulatorData_->mdrunOptions,
491 legacySimulatorData_->inputrec,
492 legacySimulatorData_->fr))
493 {
494 algorithm.pmeLoadBalanceHelper_ = std::make_unique<PmeLoadBalanceHelper>(
495 legacySimulatorData_->mdrunOptions.verbose, algorithm.statePropagatorData_.get(),
496 legacySimulatorData_->fplog, legacySimulatorData_->cr, legacySimulatorData_->mdlog,
497 legacySimulatorData_->inputrec, legacySimulatorData_->wcycle, legacySimulatorData_->fr);
498 registerWithInfrastructureAndSignallers(algorithm.pmeLoadBalanceHelper_.get());
499 }
500 // Build domdec helper
501 if (DOMAINDECOMP(legacySimulatorData_->cr))
502 {
503 algorithm.domDecHelper_ = std::make_unique<DomDecHelper>(
504 legacySimulatorData_->mdrunOptions.verbose,
505 legacySimulatorData_->mdrunOptions.verboseStepPrintInterval,
506 algorithm.statePropagatorData_.get(), algorithm.freeEnergyPerturbationData_.get(),
507 algorithm.topologyHolder_.get(),
508 globalCommunicationHelper_.moveCheckBondedInteractionsCallback(),
509 globalCommunicationHelper_.nstglobalcomm(), legacySimulatorData_->fplog,
510 legacySimulatorData_->cr, legacySimulatorData_->mdlog, legacySimulatorData_->constr,
511 legacySimulatorData_->inputrec, legacySimulatorData_->mdAtoms, legacySimulatorData_->nrnb,
512 legacySimulatorData_->wcycle, legacySimulatorData_->fr, legacySimulatorData_->vsite,
513 legacySimulatorData_->imdSession, legacySimulatorData_->pull_work);
514 registerWithInfrastructureAndSignallers(algorithm.domDecHelper_.get());
515 }
516
517 // Build trajectory element
518 auto trajectoryElement = trajectoryElementBuilder_.build(
519 legacySimulatorData_->fplog, legacySimulatorData_->nfile, legacySimulatorData_->fnm,
520 legacySimulatorData_->mdrunOptions, legacySimulatorData_->cr,
521 legacySimulatorData_->outputProvider, legacySimulatorData_->mdModulesNotifier,
522 legacySimulatorData_->inputrec, legacySimulatorData_->top_global,
523 legacySimulatorData_->oenv, legacySimulatorData_->wcycle,
524 legacySimulatorData_->startingBehavior, simulationsShareState);
525 registerWithInfrastructureAndSignallers(trajectoryElement.get());
526
527 // Build free energy element
528 std::unique_ptr<FreeEnergyPerturbationData::Element> freeEnergyPerturbationElement = nullptr;
529 if (algorithm.freeEnergyPerturbationData_)
530 {
531 freeEnergyPerturbationElement = std::make_unique<FreeEnergyPerturbationData::Element>(
532 algorithm.freeEnergyPerturbationData_.get(),
533 legacySimulatorData_->inputrec->fepvals->delta_lambda);
534 registerWithInfrastructureAndSignallers(freeEnergyPerturbationElement.get());
535 }
536
537 // Build checkpoint helper (do this last so everyone else can be a checkpoint client!)
538 {
539 checkpointHelperBuilder_.setCheckpointHandler(std::make_unique<CheckpointHandler>(
540 compat::make_not_null<SimulationSignal*>(&(*algorithm.signals_)[eglsCHKPT]),
541 simulationsShareState, legacySimulatorData_->inputrec->nstlist == 0,
542 MASTER(legacySimulatorData_->cr), legacySimulatorData_->mdrunOptions.writeConfout,
543 legacySimulatorData_->mdrunOptions.checkpointOptions.period));
544 algorithm.checkpointHelper_ = checkpointHelperBuilder_.build(
545 legacySimulatorData_->inputrec->init_step, trajectoryElement.get(),
546 legacySimulatorData_->fplog, legacySimulatorData_->cr,
547 legacySimulatorData_->observablesHistory, legacySimulatorData_->walltime_accounting,
548 legacySimulatorData_->state_global, legacySimulatorData_->mdrunOptions.writeConfout);
549 registerWithInfrastructureAndSignallers(algorithm.checkpointHelper_.get());
550 }
551
552 // Build signallers
553 {
554 /* Signallers need to be called in an exact order. Some signallers are clients
555 * of other signallers, which requires the clients signallers to be called
556 * _after_ any signaller they are registered to - otherwise, they couldn't
557 * adapt their behavior to the information they got signalled.
558 *
559 * Signallers being clients of other signallers require registration.
560 * That registration happens during construction, which in turn means that
561 * we want to construct the signallers in the reverse order of their later
562 * call order.
563 *
564 * For the above reasons, the `addSignaller` lambda defined below emplaces
565 * added signallers at the beginning of the signaller list, which will yield
566 * a signaller list which is inverse to the build order (and hence equal to
567 * the intended call order).
568 */
569 auto addSignaller = [this, &algorithm](auto signaller) {
570 registerWithInfrastructureAndSignallers(signaller.get());
571 algorithm.signallerList_.emplace(algorithm.signallerList_.begin(), std::move(signaller));
572 };
573 const auto* inputrec = legacySimulatorData_->inputrec;
574 addSignaller(energySignallerBuilder_.build(
575 inputrec->nstcalcenergy, inputrec->fepvals->nstdhdl, inputrec->nstpcouple));
576 addSignaller(trajectorySignallerBuilder_.build(
577 inputrec->nstxout, inputrec->nstvout, inputrec->nstfout,
578 inputrec->nstxout_compressed, trajectoryElement->tngBoxOut(),
579 trajectoryElement->tngLambdaOut(), trajectoryElement->tngBoxOutCompressed(),
580 trajectoryElement->tngLambdaOutCompressed(), inputrec->nstenergy));
581 addSignaller(loggingSignallerBuilder_.build(inputrec->nstlog, inputrec->init_step,
582 legacySimulatorData_->startingBehavior));
583 addSignaller(lastStepSignallerBuilder_.build(inputrec->nsteps, inputrec->init_step,
584 algorithm.stopHandler_.get()));
585 addSignaller(neighborSearchSignallerBuilder_.build(inputrec->nstlist, inputrec->init_step,
586 inputrec->init_t));
587 }
588
589 // Create element list
590 // Checkpoint helper needs to be in the call list (as first element!) to react to last step
591 algorithm.elementCallList_.emplace_back(algorithm.checkpointHelper_.get());
592 // Next, update the free energy lambda vector if needed
593 if (freeEnergyPerturbationElement)
594 {
595 algorithm.elementsOwnershipList_.emplace_back(std::move(freeEnergyPerturbationElement));
596 algorithm.elementCallList_.emplace_back(algorithm.elementsOwnershipList_.back().get());
597 }
598 // Then, move the built algorithm
599 algorithm.elementsOwnershipList_.insert(algorithm.elementsOwnershipList_.end(),
600 std::make_move_iterator(elements_.begin()),
601 std::make_move_iterator(elements_.end()));
602 algorithm.elementCallList_.insert(algorithm.elementCallList_.end(),
603 std::make_move_iterator(callList_.begin()),
604 std::make_move_iterator(callList_.end()));
605 // Finally, all trajectory writing is happening after the step
606 // (relevant data was stored by elements through energy signaller)
607 algorithm.elementsOwnershipList_.emplace_back(std::move(trajectoryElement));
608 algorithm.elementCallList_.emplace_back(algorithm.elementsOwnershipList_.back().get());
609
610 algorithm.setup();
611 return algorithm;
612 }
613
addElementToSimulatorAlgorithm(std::unique_ptr<ISimulatorElement> element)614 ISimulatorElement* ModularSimulatorAlgorithmBuilder::addElementToSimulatorAlgorithm(
615 std::unique_ptr<ISimulatorElement> element)
616 {
617 elements_.emplace_back(std::move(element));
618 return elements_.back().get();
619 }
620
elementExists(const ISimulatorElement * element) const621 bool ModularSimulatorAlgorithmBuilder::elementExists(const ISimulatorElement* element) const
622 {
623 // Check whether element exists in element list
624 if (std::any_of(elements_.begin(), elements_.end(),
625 [element](auto& existingElement) { return element == existingElement.get(); }))
626 {
627 return true;
628 }
629 // Check whether element exists in other places controlled by *this
630 return (statePropagatorData_->element() == element || energyData_->element() == element
631 || (freeEnergyPerturbationData_ && freeEnergyPerturbationData_->element() == element));
632 }
633
addElementToSetupTeardownList(ISimulatorElement * element)634 void ModularSimulatorAlgorithmBuilder::addElementToSetupTeardownList(ISimulatorElement* element)
635 {
636 // Add element if it's not already in the list
637 if (std::find(setupAndTeardownList_.begin(), setupAndTeardownList_.end(), element)
638 == setupAndTeardownList_.end())
639 {
640 setupAndTeardownList_.emplace_back(element);
641 }
642 }
643
registerLastStepCallback()644 std::optional<SignallerCallback> ModularSimulatorAlgorithm::SignalHelper::registerLastStepCallback()
645 {
646 return [this](Step step, Time gmx_unused time) { this->lastStep_ = step; };
647 }
648
registerNSCallback()649 std::optional<SignallerCallback> ModularSimulatorAlgorithm::SignalHelper::registerNSCallback()
650 {
651 return [this](Step step, Time gmx_unused time) { this->nextNSStep_ = step; };
652 }
653
GlobalCommunicationHelper(int nstglobalcomm,SimulationSignals * simulationSignals)654 GlobalCommunicationHelper::GlobalCommunicationHelper(int nstglobalcomm, SimulationSignals* simulationSignals) :
655 nstglobalcomm_(nstglobalcomm),
656 simulationSignals_(simulationSignals)
657 {
658 }
659
nstglobalcomm() const660 int GlobalCommunicationHelper::nstglobalcomm() const
661 {
662 return nstglobalcomm_;
663 }
664
simulationSignals()665 SimulationSignals* GlobalCommunicationHelper::simulationSignals()
666 {
667 return simulationSignals_;
668 }
669
setCheckBondedInteractionsCallback(CheckBondedInteractionsCallback callback)670 void GlobalCommunicationHelper::setCheckBondedInteractionsCallback(CheckBondedInteractionsCallback callback)
671 {
672 checkBondedInteractionsCallback_ = std::move(callback);
673 }
674
moveCheckBondedInteractionsCallback()675 CheckBondedInteractionsCallback GlobalCommunicationHelper::moveCheckBondedInteractionsCallback()
676 {
677 if (!checkBondedInteractionsCallback_)
678 {
679 throw SimulationAlgorithmSetupError(
680 "Requested CheckBondedInteractionsCallback before it was set.");
681 }
682 return *checkBondedInteractionsCallback_;
683 }
684
ModularSimulatorAlgorithmBuilderHelper(ModularSimulatorAlgorithmBuilder * builder)685 ModularSimulatorAlgorithmBuilderHelper::ModularSimulatorAlgorithmBuilderHelper(
686 ModularSimulatorAlgorithmBuilder* builder) :
687 builder_(builder)
688 {
689 }
690
storeElement(std::unique_ptr<ISimulatorElement> element)691 ISimulatorElement* ModularSimulatorAlgorithmBuilderHelper::storeElement(std::unique_ptr<ISimulatorElement> element)
692 {
693 return builder_->addElementToSimulatorAlgorithm(std::move(element));
694 }
695
elementIsStored(const ISimulatorElement * element) const696 bool ModularSimulatorAlgorithmBuilderHelper::elementIsStored(const ISimulatorElement* element) const
697 {
698 return builder_->elementExists(element);
699 }
700
getStoredValue(const std::string & key) const701 std::optional<std::any> ModularSimulatorAlgorithmBuilderHelper::getStoredValue(const std::string& key) const
702 {
703 const auto iter = values_.find(key);
704 if (iter == values_.end())
705 {
706 return std::nullopt;
707 }
708 else
709 {
710 return iter->second;
711 }
712 }
713
registerThermostat(std::function<void (const PropagatorThermostatConnection &)> registrationFunction)714 void ModularSimulatorAlgorithmBuilderHelper::registerThermostat(
715 std::function<void(const PropagatorThermostatConnection&)> registrationFunction)
716 {
717 builder_->thermostatRegistrationFunctions_.emplace_back(std::move(registrationFunction));
718 }
719
registerBarostat(std::function<void (const PropagatorBarostatConnection &)> registrationFunction)720 void ModularSimulatorAlgorithmBuilderHelper::registerBarostat(
721 std::function<void(const PropagatorBarostatConnection&)> registrationFunction)
722 {
723 builder_->barostatRegistrationFunctions_.emplace_back(std::move(registrationFunction));
724 }
725
registerWithThermostat(PropagatorThermostatConnection connectionData)726 void ModularSimulatorAlgorithmBuilderHelper::registerWithThermostat(PropagatorThermostatConnection connectionData)
727 {
728 builder_->propagatorThermostatConnections_.emplace_back(std::move(connectionData));
729 }
730
registerWithBarostat(PropagatorBarostatConnection connectionData)731 void ModularSimulatorAlgorithmBuilderHelper::registerWithBarostat(PropagatorBarostatConnection connectionData)
732 {
733 builder_->propagatorBarostatConnections_.emplace_back(std::move(connectionData));
734 }
735
736
737 } // namespace gmx
738