2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2017,2018,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.
38 * Implements bias sharing checking functionality.
40 * \author Berk Hess <hess@kth.se>
46 #include "biassharing.h"
54 #include "gromacs/gmxlib/network.h"
55 #include "gromacs/mdrunutility/multisim.h"
56 #include "gromacs/mdtypes/awh_params.h"
57 #include "gromacs/mdtypes/commrec.h"
58 #include "gromacs/utility/arrayref.h"
59 #include "gromacs/utility/exceptions.h"
60 #include "gromacs/utility/gmxassert.h"
61 #include "gromacs/utility/stringutil.h"
69 //! Determines and returns which of the local biases are shared with who how many other simulations
70 std::multiset<int> getGlobalShareIndices(ArrayRef<const int> localShareIndices, MPI_Comm simulationMastersComm)
74 MPI_Comm_size(simulationMastersComm, &numSimulations);
76 MPI_Comm_rank(simulationMastersComm, &ourRank);
77 std::vector<int> biasCountsIn(numSimulations, 0);
78 std::vector<int> biasCounts(numSimulations, 0);
79 biasCountsIn[ourRank] = localShareIndices.size();
80 MPI_Allreduce(biasCountsIn.data(), biasCounts.data(), numSimulations, MPI_INT, MPI_SUM, simulationMastersComm);
81 // Now we need to gather the share indices to all (master) ranks.
82 // We could use MPI_Allgatherv, but thread-MPI does not support that and using
83 // MPI_Allreduce produces simpler code, so we use that.
86 for (int rank = 0; rank < numSimulations; rank++)
90 ourOffset = totNumBiases;
92 totNumBiases += biasCounts[rank];
94 // Fill a buffer with zeros and our part of sharing indices
95 std::vector<int> shareIndicesAllIn(totNumBiases, 0);
96 std::copy(localShareIndices.begin(), localShareIndices.end(), shareIndicesAllIn.begin() + ourOffset);
97 // Gather all sharing indices to all (master) ranks
98 std::vector<int> shareIndicesAll(totNumBiases);
99 MPI_Allreduce(shareIndicesAllIn.data(), shareIndicesAll.data(), totNumBiases, MPI_INT, MPI_SUM, simulationMastersComm);
101 GMX_UNUSED_VALUE(simulationMastersComm);
103 ArrayRef<const int> shareIndicesAll = localShareIndices;
106 std::multiset<int> shareIndicesSet;
107 for (int shareIndex : shareIndicesAll)
111 shareIndicesSet.insert(shareIndex);
115 return shareIndicesSet;
120 BiasSharing::BiasSharing(const AwhParams& awhParams, const t_commrec& commRecord, MPI_Comm simulationMastersComm) :
121 commRecord_(commRecord)
123 if (MASTER(&commRecord))
125 std::vector<int> localShareIndices;
126 int shareGroupPrev = 0;
127 for (int k = 0; k < awhParams.numBias(); k++)
129 const int shareGroup = awhParams.awhBiasParams()[k].shareGroup();
130 GMX_RELEASE_ASSERT(shareGroup >= 0, "Bias share group values should be >= 0");
131 localShareIndices.push_back(shareGroup);
134 if (shareGroup <= shareGroupPrev)
137 InvalidInputError("AWH biases that are shared should use increasing "
138 "share-group values"));
140 shareGroupPrev = shareGroup;
143 std::multiset<int> globalShareIndices =
144 getGlobalShareIndices(localShareIndices, simulationMastersComm);
146 int numSimulations = 1;
148 MPI_Comm_size(simulationMastersComm, &numSimulations);
150 MPI_Comm_rank(simulationMastersComm, &myRank);
153 numSharingSimulations_.resize(awhParams.numBias(), 1);
154 sharingSimulationIndices_.resize(awhParams.numBias(), 0);
155 multiSimCommPerBias_.resize(awhParams.numBias(), MPI_COMM_NULL);
157 for (int shareIndex : globalShareIndices)
159 if (globalShareIndices.count(shareIndex) > 1)
161 const auto& findBiasIndex =
162 std::find(localShareIndices.begin(), localShareIndices.end(), shareIndex);
163 const index localBiasIndex = (findBiasIndex == localShareIndices.end()
165 : findBiasIndex - localShareIndices.begin());
167 if (static_cast<int>(globalShareIndices.count(shareIndex)) == numSimulations)
169 splitComm = simulationMastersComm;
174 const int haveLocally = (localBiasIndex >= 0 ? 1 : 0);
175 MPI_Comm_split(simulationMastersComm, haveLocally, myRank, &splitComm);
176 createdCommList_.push_back(splitComm);
178 GMX_RELEASE_ASSERT(false, "Can not have sharing without MPI");
181 if (localBiasIndex >= 0)
183 numSharingSimulations_[localBiasIndex] = globalShareIndices.count(shareIndex);
185 MPI_Comm_rank(splitComm, &sharingSimulationIndices_[localBiasIndex]);
187 multiSimCommPerBias_[localBiasIndex] = splitComm;
194 if (commRecord.nnodes > 1)
196 numSharingSimulations_.resize(awhParams.numBias());
198 numSharingSimulations_.data(), numSharingSimulations_.size(), MPI_INT, 0, commRecord.mpi_comm_mygroup);
203 BiasSharing::~BiasSharing()
206 for (MPI_Comm comm : createdCommList_)
208 MPI_Comm_free(&comm);
219 std::enable_if_t<std::is_same_v<T, int>, MPI_Datatype> mpiType()
225 std::enable_if_t<std::is_same_v<T, long>, MPI_Datatype> mpiType()
231 std::enable_if_t<std::is_same_v<T, double>, MPI_Datatype> mpiType()
241 * Sum an array over all simulations on master ranks or all ranks of each simulation.
243 * This assumes the data is identical on all ranks within each simulation.
245 * \param[in,out] data The data to sum.
246 * \param[in] multiSimComm Communicator for the master ranks of sharing simulations.
247 * \param[in] broadcastWithinSimulation Broadcast the result to all ranks within the simulation
248 * \param[in] commRecord Struct for intra-simulation communication.
251 void sumOverSimulations(ArrayRef<T> data,
252 MPI_Comm multiSimComm,
253 const bool broadcastWithinSimulation,
254 const t_commrec& commRecord)
257 if (MASTER(&commRecord))
259 MPI_Allreduce(MPI_IN_PLACE, data.data(), data.size(), mpiType<T>(), MPI_SUM, multiSimComm);
261 if (broadcastWithinSimulation && commRecord.nnodes > 1)
263 gmx_bcast(data.size() * sizeof(T), data.data(), commRecord.mpi_comm_mygroup);
266 GMX_UNUSED_VALUE(data);
267 GMX_UNUSED_VALUE(commRecord);
268 GMX_UNUSED_VALUE(broadcastWithinSimulation);
269 GMX_UNUSED_VALUE(multiSimComm);
273 void BiasSharing::sumOverMasterRanks(ArrayRef<int> data, const int biasIndex) const
275 sumOverSimulations(data, multiSimCommPerBias_[biasIndex], false, commRecord_);
278 void BiasSharing::sumOverMasterRanks(ArrayRef<long> data, const int biasIndex) const
280 sumOverSimulations(data, multiSimCommPerBias_[biasIndex], false, commRecord_);
283 void BiasSharing::sum(ArrayRef<int> data, const int biasIndex) const
285 sumOverSimulations(data, multiSimCommPerBias_[biasIndex], true, commRecord_);
288 void BiasSharing::sum(ArrayRef<double> data, const int biasIndex) const
290 sumOverSimulations(data, multiSimCommPerBias_[biasIndex], true, commRecord_);
293 bool haveBiasSharingWithinSimulation(const AwhParams& awhParams)
295 bool haveSharing = false;
297 for (int k = 0; k < awhParams.numBias(); k++)
299 int shareGroup = awhParams.awhBiasParams()[k].shareGroup();
302 for (int i = k + 1; i < awhParams.numBias(); i++)
304 if (awhParams.awhBiasParams()[i].shareGroup() == shareGroup)
315 void biasesAreCompatibleForSharingBetweenSimulations(const AwhParams& awhParams,
316 ArrayRef<const size_t> pointSize,
317 const BiasSharing& biasSharing)
319 /* Check the point sizes. This is a sufficient condition for running
320 * as shared multi-sim run. No physics checks are performed here.
322 const auto& awhBiasParams = awhParams.awhBiasParams();
323 for (int b = 0; b < gmx::ssize(awhBiasParams); b++)
325 if (awhBiasParams[b].shareGroup() > 0)
327 const int numSim = biasSharing.numSharingSimulations(b);
330 // This bias is not actually shared
333 const int simIndex = biasSharing.sharingSimulationIndex(b);
334 std::vector<int> intervals(numSim * 2);
335 intervals[numSim * 0 + simIndex] = awhParams.nstSampleCoord();
336 intervals[numSim * 1 + simIndex] = awhParams.numSamplesUpdateFreeEnergy();
337 biasSharing.sumOverMasterRanks(intervals, b);
338 for (int sim = 1; sim < numSim; sim++)
340 if (intervals[sim] != intervals[0])
342 GMX_THROW(InvalidInputError(
343 "All simulations should have the same AWH sample interval"));
345 if (intervals[numSim + sim] != intervals[numSim])
348 InvalidInputError("All simulations should have the same AWH "
349 "free-energy update interval"));
353 std::vector<long> pointSizes(numSim);
354 pointSizes[simIndex] = pointSize[b];
355 biasSharing.sumOverMasterRanks(pointSizes, b);
356 for (int sim = 1; sim < numSim; sim++)
358 if (pointSizes[sim] != pointSizes[0])
360 GMX_THROW(InvalidInputError(
361 gmx::formatString("Shared AWH bias %d has different grid sizes in "
362 "different simulations\n",