2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2018,2019,2020, 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 * \brief Declares the multi-simulation support routines.
39 * \author Mark Abraham <mark.j.abraham@gmail.com>
41 * \ingroup module_mdrunutility
43 #ifndef GMX_MDRUNUTILITY_MULTISIM_H
44 #define GMX_MDRUNUTILITY_MULTISIM_H
50 #include "gromacs/utility/arrayref.h"
51 #include "gromacs/utility/gmxmpi.h"
53 struct gmx_multisim_t;
56 * \brief Builder function for gmx_multisim_t
58 * \param[in] worldComm MPI communicator to split when
59 * multi-simulation is requested.
60 * \param[in] multidirs Strings naming the subdirectories when
61 * multi-simulation is requested, otherwise empty
63 * Splits \c worldComm into \c multidirs.size() separate
64 * simulations, if >1, and creates a communication structure
65 * between the master ranks of these simulations.
67 * Valid to call regardless of build configuration, but \c
68 * multidirs must be empty unless a real MPI build is used.
70 * \throws NotImplementedError when \c multidirs is non-empty unless using real MPI is true
71 * \throws NotImplementedError when \c multidirs has exactly one element
72 * \throws InconsistentInputError when the number of MPI ranks is not a multiple of the number of \c multidirs
73 * \throws FileIOError when the simulation cannot change to the working directory in \c multidirs
75 std::unique_ptr<gmx_multisim_t> buildMultiSimulation(MPI_Comm worldComm,
76 gmx::ArrayRef<const std::string> multidirs);
79 * \brief Coordinate multi-simulation resources for mdrun
81 * \todo Change this to class
85 /*! \brief Constructor
87 * \param[in] numSimulations The number of simulations in the MPI world.
88 * \param[in] simulationIndex The index of this simulation in the set of simulations.
89 * \param[in] mastersComm On master ranks, the communicator among master ranks;
90 * otherwise MPI_COMM_NULL.
91 * \param[in] simulationComm The communicator among ranks of this simulation.
93 * Assumes ownership of the communicators if they are neither
94 * MPI_COMM_WORLD nor MPI_COMM_NULL. If so, upon destruction will
95 * call MPI_Comm_free on them.
97 gmx_multisim_t(int numSimulations, int simulationIndex, MPI_Comm mastersComm, MPI_Comm simulationComm);
101 //! The number of simulations in the set of multi-simulations
102 int numSimulations_ = 1;
103 //! The index of the simulation that owns this object within the set
104 int simulationIndex_ = 0;
105 //! The MPI communicator between master ranks of simulations, valid only on master ranks.
106 MPI_Comm mastersComm_ = MPI_COMM_NULL;
107 //! The MPI communicator between ranks of this simulation.
108 MPI_Comm simulationComm_ = MPI_COMM_NULL;
109 /*! \brief Communication buffers needed if MPI_IN_PLACE isn't supported
111 * Other types could be added as needed.
113 * These vectors are unused when MPI_IN_PLACE is available
114 * and could be removed with preprocessing (or perhaps
115 * templating) or simply requiring MPI 2.0 (the standard
116 * introduced in 1997). However, the additional cache pressure
117 * introduced by the extra size of this type is not of great
118 * concern, since we have at most one per MPI rank.
119 * See issue #3591. */
121 std::vector<int> intBuffer_;
122 std::vector<int64_t> int64Buffer_;
126 //! Calculate the sum over the simulations of an array of ints
127 void gmx_sumi_sim(int nr, int r[], const gmx_multisim_t* ms);
129 //! Calculate the sum over the simulations of an array of large ints
130 void gmx_sumli_sim(int nr, int64_t r[], const gmx_multisim_t* ms);
132 //! Calculate the sum over the simulations of an array of floats
133 void gmx_sumf_sim(int nr, float r[], const gmx_multisim_t* ms);
135 //! Calculate the sum over the simulations of an array of doubles
136 void gmx_sumd_sim(int nr, double r[], const gmx_multisim_t* ms);
138 /*! \brief Return a vector containing the gathered values of \c
139 * localValue found on the master rank of each simulation. */
140 std::vector<int> gatherIntFromMultiSimulation(const gmx_multisim_t* ms, int localValue);
142 /*! \brief Check if val is the same on all simulations for a mdrun
145 * The string name is used to print to the log file and in a fatal error
146 * if the val's don't match. If bQuiet is true and the check passes,
147 * no output is written. */
148 void check_multi_int(FILE* log, const gmx_multisim_t* ms, int val, const char* name, gmx_bool bQuiet);
149 /*! \copydoc check_multi_int() */
150 void check_multi_int64(FILE* log, const gmx_multisim_t* ms, int64_t val, const char* name, gmx_bool bQuiet);
153 //! Convenience define for sum of reals
154 # define gmx_sum_sim gmx_sumd_sim
156 //! Convenience define for sum of reals
157 # define gmx_sum_sim gmx_sumf_sim
160 //! Are we doing multiple independent simulations?
161 static bool inline isMultiSim(const gmx_multisim_t* ms)
163 return ms != nullptr;
166 /*! \brief Return whether this rank is the master rank of a
167 * simulation, using \c ms (if it is valid) and otherwise \c
169 bool findIsSimulationMasterRank(const gmx_multisim_t* ms, MPI_Comm communicator);
171 //! Are we the master simulation of a possible multi-simulation?
172 bool isMasterSim(const gmx_multisim_t* ms);
174 /*! \brief Are we the master rank (of the master simulation, for a multi-sim).
176 * This rank prints the remaining run time etc. */
177 bool isMasterSimMasterRank(const gmx_multisim_t* ms, bool isMaster);