From 9cf1bd23abe8ac9a4e2638845d34ae1b1ccf19b7 Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Thu, 21 Oct 2021 08:51:00 +0000 Subject: [PATCH] Improve AWH bias sharing documentation --- docs/reference-manual/special/awh.rst | 11 +++++++++ docs/user-guide/mdp-options.rst | 14 +++++------ .../applied_forces/awh/biassharing.cpp | 12 +++++----- src/gromacs/applied_forces/awh/biassharing.h | 24 ++++++++++++------- src/gromacs/applied_forces/awh/biasstate.cpp | 13 +++++----- .../applied_forces/awh/tests/biassharing.cpp | 2 +- 6 files changed, 48 insertions(+), 28 deletions(-) diff --git a/docs/reference-manual/special/awh.rst b/docs/reference-manual/special/awh.rst index f9e30ea273..1cc7e935db 100644 --- a/docs/reference-manual/special/awh.rst +++ b/docs/reference-manual/special/awh.rst @@ -498,6 +498,17 @@ equal to the length of the sampling interval, the sampling interval is considered covered when at least one walker has independently traversed the sampling interval. +In practice biases are shared by setting :mdp:`awh-share-multisim` to true +and :mdp:`awh1-share-group` (for bias 1) to a non-zero value. Here, bias 1 +will be shared between simulations that have the same share group value. +Sharing can be different for bias 1, 2, etc. (although there are +few use cases where this is useful). Technically there are no restrictions +on sharing, apart from that biases that are shared need to have the same +number of grid points and the update intervals should match. +Note that biases can not be shared within a simulation. +The latter could be useful, especially for multimeric proteins, but this +is more difficult to implement. + .. _awhreweight: Reweighting and combining biased data diff --git a/docs/user-guide/mdp-options.rst b/docs/user-guide/mdp-options.rst index c9ce3a20a0..d0e990632a 100644 --- a/docs/user-guide/mdp-options.rst +++ b/docs/user-guide/mdp-options.rst @@ -2095,15 +2095,15 @@ AWH adaptive biasing .. mdp-value:: positive - Share the bias and PMF estimates within and/or between simulations. - Within a simulation, the bias will be shared between biases that have the - same :mdp:`awh1-share-group` index (note that the current code does not support this). - With :mdp-value:`awh-share-multisim=yes` and - :ref:`gmx mdrun` option ``-multidir`` the bias will also be shared across simulations. + Share the bias and PMF estimates between simulations. This currently + only works between biases with the same index. Note that currently + sharing within a single simulation is not supported. + The bias will be shared across simulations that specify the same + value for :mdp:`awh1-share-group`. To enable this, use + :mdp-value:`awh-share-multisim=yes` and the :ref:`gmx mdrun` option + ``-multidir``. Sharing may increase convergence initially, although the starting configurations can be critical, especially when sharing between many biases. - Currently, the share value should increase with increasing bias index - (or be 0). .. mdp:: awh1-ndim diff --git a/src/gromacs/applied_forces/awh/biassharing.cpp b/src/gromacs/applied_forces/awh/biassharing.cpp index 176d0c350d..8c519a5640 100644 --- a/src/gromacs/applied_forces/awh/biassharing.cpp +++ b/src/gromacs/applied_forces/awh/biassharing.cpp @@ -270,22 +270,22 @@ void sumOverSimulations(ArrayRef data, #endif // GMX_MPI } -void BiasSharing::sumOverMasterRanks(ArrayRef data, const int biasIndex) const +void BiasSharing::sumOverSharingMasterRanks(ArrayRef data, const int biasIndex) const { sumOverSimulations(data, multiSimCommPerBias_[biasIndex], false, commRecord_); } -void BiasSharing::sumOverMasterRanks(ArrayRef data, const int biasIndex) const +void BiasSharing::sumOverSharingMasterRanks(ArrayRef data, const int biasIndex) const { sumOverSimulations(data, multiSimCommPerBias_[biasIndex], false, commRecord_); } -void BiasSharing::sum(ArrayRef data, const int biasIndex) const +void BiasSharing::sumOverSharingSimulations(ArrayRef data, const int biasIndex) const { sumOverSimulations(data, multiSimCommPerBias_[biasIndex], true, commRecord_); } -void BiasSharing::sum(ArrayRef data, const int biasIndex) const +void BiasSharing::sumOverSharingSimulations(ArrayRef data, const int biasIndex) const { sumOverSimulations(data, multiSimCommPerBias_[biasIndex], true, commRecord_); } @@ -334,7 +334,7 @@ void biasesAreCompatibleForSharingBetweenSimulations(const AwhParams& awhP std::vector intervals(numSim * 2); intervals[numSim * 0 + simIndex] = awhParams.nstSampleCoord(); intervals[numSim * 1 + simIndex] = awhParams.numSamplesUpdateFreeEnergy(); - biasSharing.sumOverMasterRanks(intervals, b); + biasSharing.sumOverSharingMasterRanks(intervals, b); for (int sim = 1; sim < numSim; sim++) { if (intervals[sim] != intervals[0]) @@ -352,7 +352,7 @@ void biasesAreCompatibleForSharingBetweenSimulations(const AwhParams& awhP std::vector pointSizes(numSim); pointSizes[simIndex] = pointSize[b]; - biasSharing.sumOverMasterRanks(pointSizes, b); + biasSharing.sumOverSharingMasterRanks(pointSizes, b); for (int sim = 1; sim < numSim; sim++) { if (pointSizes[sim] != pointSizes[0]) diff --git a/src/gromacs/applied_forces/awh/biassharing.h b/src/gromacs/applied_forces/awh/biassharing.h index e73ec5118a..5f571d73c2 100644 --- a/src/gromacs/applied_forces/awh/biassharing.h +++ b/src/gromacs/applied_forces/awh/biassharing.h @@ -85,17 +85,25 @@ public: //! Returns the index of our simulation in the simulations sharing bias \p biasIndex int sharingSimulationIndex(int biasIndex) const { return sharingSimulationIndices_[biasIndex]; } - //! Sums data over the master ranks of all simulations sharing bias \p biasIndex - void sumOverMasterRanks(ArrayRef data, int biasIndex) const; + //! Sums data of type int over the master ranks of all simulations sharing bias \p biasIndex + void sumOverSharingMasterRanks(ArrayRef data, int biasIndex) const; - //! Sums data over the master ranks of all simulations sharing bias \p biasIndex - void sumOverMasterRanks(ArrayRef data, int biasIndex) const; + //! Sums data of type long over the master ranks of all simulations sharing bias \p biasIndex + void sumOverSharingMasterRanks(ArrayRef data, int biasIndex) const; - //! Sums data over all simulations sharing bias \p biasIndex and broadcasts to all ranks within the simulation - void sum(ArrayRef data, int biasIndex) const; + /*! \brief Sums data of type int over all simulations sharing bias \p biasIndex + * + * The summing is performed over the master ranks of the simulations sharing bias \p biasIndex + * and the result is broadcasted to the other ranks within each simulation. + */ + void sumOverSharingSimulations(ArrayRef data, int biasIndex) const; - //! Sums data over all simulations sharing bias \p biasIndex and broadcasts to all ranks within the simulation - void sum(ArrayRef data, int biasIndex) const; + /*! \brief Sums data of type long over all simulations sharing bias \p biasIndex + * + * The summing is performed over the master ranks of the simulations sharing bias \p biasIndex + * and the result is broadcasted to the other ranks within each simulation. + */ + void sumOverSharingSimulations(ArrayRef data, int biasIndex) const; private: //! The number of simulations sharing for each bias diff --git a/src/gromacs/applied_forces/awh/biasstate.cpp b/src/gromacs/applied_forces/awh/biasstate.cpp index 3a43f29db1..7c8b1d0e2b 100644 --- a/src/gromacs/applied_forces/awh/biasstate.cpp +++ b/src/gromacs/applied_forces/awh/biasstate.cpp @@ -124,7 +124,7 @@ void sumPmf(gmx::ArrayRef pointState, int numSharedUpdate, const Bia buffer[i] = pointState[i].inTargetRegion() ? std::exp(pointState[i].logPmfSum()) : 0; } - biasSharing->sum(gmx::ArrayRef(buffer), biasIndex); + biasSharing->sumOverSharingSimulations(gmx::ArrayRef(buffer), biasIndex); /* Take log again to get (non-normalized) PMF */ double normFac = 1.0 / numSharedUpdate; @@ -671,7 +671,7 @@ void mergeSharedUpdateLists(std::vector* updateList, } /* Sum over the sims to get all the flagged points */ - biasSharing.sum(arrayRefFromArray(numUpdatesOfPoint.data(), numPoints), biasIndex); + biasSharing.sumOverSharingSimulations(arrayRefFromArray(numUpdatesOfPoint.data(), numPoints), biasIndex); /* Collect the indices of the flagged points in place. The resulting array will be the merged update list.*/ updateList->clear(); @@ -790,8 +790,8 @@ void sumHistograms(gmx::ArrayRef pointState, coordVisits[localIndex] = ps.numVisitsIteration(); } - biasSharing->sum(gmx::ArrayRef(weightSum), biasIndex); - biasSharing->sum(gmx::ArrayRef(coordVisits), biasIndex); + biasSharing->sumOverSharingSimulations(gmx::ArrayRef(weightSum), biasIndex); + biasSharing->sumOverSharingSimulations(gmx::ArrayRef(coordVisits), biasIndex); /* Transfer back the result */ for (size_t localIndex = 0; localIndex < localUpdateList.size(); localIndex++) @@ -1025,8 +1025,9 @@ bool BiasState::isSamplingRegionCovered(const BiasParams& params, /* For multiple dimensions this may not be the best way to do it. */ for (int d = 0; d < grid.numDimensions(); d++) { - biasSharing_->sum(gmx::arrayRefFromArray(checkDim[d].covered.data(), grid.axis(d).numPoints()), - params.biasIndex); + biasSharing_->sumOverSharingSimulations( + gmx::arrayRefFromArray(checkDim[d].covered.data(), grid.axis(d).numPoints()), + params.biasIndex); } } diff --git a/src/gromacs/applied_forces/awh/tests/biassharing.cpp b/src/gromacs/applied_forces/awh/tests/biassharing.cpp index 70f75e41cb..c77a7a503f 100644 --- a/src/gromacs/applied_forces/awh/tests/biassharing.cpp +++ b/src/gromacs/applied_forces/awh/tests/biassharing.cpp @@ -101,7 +101,7 @@ void parallelTestFunction(const void gmx_unused* dummy) const std::array input{ 1, 2, 4, 8 }; std::array buffer{ input[myRank] }; - biasSharing.sumOverMasterRanks(buffer, 0); + biasSharing.sumOverSharingMasterRanks(buffer, 0); int expectedSum = 0; for (int i = 0; i < c_numSharingBiases; i++) { -- 2.22.0