2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2019,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.
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.
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.
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.
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.
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.
35 /*! \libinternal \file
37 * Declares gmx::MDModulesNotifiers.
39 * \author Christian Blau <blau@kth.se>
41 * \ingroup module_utility
44 #ifndef GMX_UTILITY_MDMODULESNOTIFIERS_H
45 #define GMX_UTILITY_MDMODULESNOTIFIERS_H
50 #include "gromacs/utility/mdmodulesnotifier.h"
54 enum class PbcType : int;
59 class KeyValueTreeObject;
60 class KeyValueTreeObjectBuilder;
61 class LocalAtomSetManager;
62 class IndexGroupsAndNames;
63 class SeparatePmeRanksPermitted;
64 struct MDModulesCheckpointReadingDataOnMaster;
65 struct MDModulesCheckpointReadingBroadcast;
66 struct MDModulesWriteCheckpointData;
68 /*! \libinternal \brief Check if module outputs energy to a specific field.
70 * Ensures that energy is output for this module.
72 struct MDModulesEnergyOutputToDensityFittingRequestChecker
74 //! Trigger output to density fitting energy field
75 bool energyOutputToDensityFitting_ = false;
79 * \brief Collect errors for the energy calculation frequency.
81 * Collect errors regarding energy calculation frequencies as strings that then
82 * may be used to issue errors.
84 * \note The mdp option "nstcalcenergy" is altered after reading the .mdp input
85 * and only used in certain integrators, thus this class is to be used
86 * only after all these operations are done.
88 class EnergyCalculationFrequencyErrors
91 //! Construct by setting the energy calculation frequency
92 EnergyCalculationFrequencyErrors(int64_t energyCalculationIntervalInSteps) :
93 energyCalculationIntervalInSteps_(energyCalculationIntervalInSteps)
96 //! Return the number of steps of an energy calculation interval
97 std::int64_t energyCalculationIntervalInSteps() const
99 return energyCalculationIntervalInSteps_;
101 //! Collect error messages
102 void addError(const std::string& errorMessage) { errorMessages_.push_back(errorMessage); }
103 //! Return error messages
104 const std::vector<std::string>& errorMessages() const { return errorMessages_; }
107 //! The frequency of energy calculations
108 const std::int64_t energyCalculationIntervalInSteps_;
109 //! The error messages
110 std::vector<std::string> errorMessages_;
113 /*! \libinternal \brief Provides the simulation time step in ps.
115 struct SimulationTimeStep
118 const double delta_t;
122 * \brief Group of notifers to organize that MDModules
123 * can receive callbacks they subscribe to.
125 * MDModules use members of this struct to subscribe to notifications
126 * of particular events. When the event occurs, the callback provided
127 * by a particular MDModule will be passed a parameter of the
128 * particular type they are interested in.
130 * Typically, during the setup phase, modules subscribe to notifiers
131 * that interest them by passing callbacks that expect a single parameter
132 * that describes the event. These are stored for later use. See the
133 * sequence diagram that follows:
138 modules [label = "mdModules:\nMDModules"],
139 notifiers [label="notifiers\nMDModulesNotifiers"],
140 notifier [label="exampleNotifier:\nBuildMDModulesNotifier\n<EventX, EventY>::type"],
141 moduleA [label="moduleA"],
142 moduleB [label="moduleB"],
143 moduleC [label="moduleC"];
145 modules box moduleC [label = "mdModules creates and owns moduleA, moduleB, and moduleC"];
146 modules =>> notifiers [label="creates"];
147 notifiers =>> notifier [label="creates"];
148 notifier =>> notifiers [label="returns"];
149 notifiers =>> modules [label="returns"];
151 modules =>> moduleA [label="provides notifiers"];
152 moduleA =>> moduleA [label="unpacks\nnotifiers.exampleNotifier"];
153 moduleA =>> notifier [label="subscribes with\ncallback(EventX&)"];
154 notifier =>> notifier [label="records subscription\nto EventX"];
155 moduleA =>> notifier [label="subscribes with\ncallback(EventY&)"];
156 notifier =>> notifier [label="records subscription\nto EventY"];
157 moduleA =>> modules [label="returns"];
159 modules =>> moduleB [label="provides notifiers"];
160 moduleB =>> moduleB [label="unpacks\nnotifiers.exampleNotifier"];
161 moduleA =>> notifier [label="subscribes with\ncallback(EventY&)"];
162 notifier =>> notifier [label="records subscription\nto EventY"];
163 moduleB =>> modules [label="returns"];
165 modules =>> moduleC [label="provides notifiers"];
166 moduleC =>> moduleC [label="unpacks and keeps\nnotifiers.exampleNotifier"];
167 moduleC =>> modules [label="returns"];
171 * When the event occurs later on, the stored callbacks are used to
172 * allow the modules to react. See the following sequence diagram,
173 * which assumes that exampleNotifier was configured as in the
174 * previous sequence diagram.
180 moduleC [label="moduleC"],
181 notifier [label="exampleNotifier:\nBuildMDModulesNotifier\n<EventX, EventY>::type"],
182 moduleA [label="moduleA"],
183 moduleB [label="moduleB"];
185 moduleC box moduleB [label = "Later, when ModuleC is doing work"];
186 moduleC =>> moduleC [label="generates EventX"];
187 moduleC =>> moduleC [label="generates EventY"];
188 moduleC =>> notifier [label="calls notify(eventX)"];
189 notifier =>> moduleA [label="calls callback(eventX)"];
190 moduleA =>> moduleA [label="reacts to eventX"];
191 moduleA =>> notifier [label="returns"];
193 notifier =>> moduleC [label="returns"];
194 moduleC =>> notifier [label="calls notify(eventY)"];
195 notifier =>> moduleA [label="calls callback(eventY)"];
196 moduleA =>> moduleA [label="reacts to eventY"];
197 moduleA =>> notifier [label="returns"];
198 notifier =>> moduleB [label="calls callback(eventY)"];
199 moduleB =>> moduleB [label="reacts to eventY"];
200 moduleB =>> notifier [label="returns"];
201 notifier =>> moduleC [label="returns"];
204 * The template arguments to the members of this struct are the
205 * parameters passed to the callback functions, one type per
206 * callback. Arguments passed as pointers are always meant to be
207 * modified, but never meant to be stored (in line with the policy
211 struct MDModulesNotifiers
213 /*! \brief Handles subscribing and calling pre-processing callback functions.
215 * EnergyCalculationFrequencyErrors* allows modules to check if they match
216 * their required calculation frequency
217 * and add their error message if needed
218 * to the collected error messages
219 * IndexGroupsAndNames provides modules with atom indices and their names
220 * KeyValueTreeObjectBuilder enables writing of module internal data to
223 BuildMDModulesNotifier<EnergyCalculationFrequencyErrors*, IndexGroupsAndNames, KeyValueTreeObjectBuilder>::type preProcessingNotifier_;
225 /*! \brief Handles subscribing and calling checkpointing callback functions.
227 * MDModulesCheckpointReadingDataOnMaster provides modules with their
228 * checkpointed data on the master
229 * node and checkpoint file version
230 * MDModulesCheckpointReadingBroadcast provides modules with a communicator
231 * and the checkpoint file version to
232 * distribute their data
233 * MDModulesWriteCheckpointData provides the modules with a key-value-tree
234 * builder to store their checkpoint data and
235 * the checkpoint file version
237 BuildMDModulesNotifier<MDModulesCheckpointReadingDataOnMaster, MDModulesCheckpointReadingBroadcast, MDModulesWriteCheckpointData>::type
238 checkpointingNotifier_;
240 /*! \brief Handles subscribing and calling callbacks during simulation setup.
242 * const KeyValueTreeObject& provides modules with the internal data they
243 * wrote to .tpr files
244 * LocalAtomSetManager* enables modules to add atom indices to local atom sets
246 * const gmx_mtop_t& provides the topology of the system to the modules
247 * MDModulesEnergyOutputToDensityFittingRequestChecker* enables modules to
248 * report if they want to write their energy output
249 * to the density fitting field in the energy files
250 * SeparatePmeRanksPermitted* enables modules to report if they want
251 * to disable dedicated PME ranks
252 * const PbcType& provides modules with the periodic boundary condition type
253 * that is used during the simulation
254 * const SimulationTimeStep& provides modules with the simulation time-step
255 * that allows them to interconvert between step
257 * const t_commrec& provides a communicator to the modules during simulation
260 BuildMDModulesNotifier<const KeyValueTreeObject&,
261 LocalAtomSetManager*,
263 MDModulesEnergyOutputToDensityFittingRequestChecker*,
264 SeparatePmeRanksPermitted*,
266 const SimulationTimeStep&,
267 const t_commrec&>::type simulationSetupNotifier_;