#include <memory>
#include "gromacs/domdec/localatomsetmanager.h"
+#include "gromacs/fileio/checkpoint.h"
#include "gromacs/fileio/mrcdensitymap.h"
#include "gromacs/math/multidimarray.h"
+#include "gromacs/mdlib/broadcaststructs.h"
+#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/imdmodule.h"
#include "gromacs/utility/classhelpers.h"
#include "gromacs/utility/exceptions.h"
* - constructing local atom sets in the simulation parameter setup
* by taking a LocalAtomSetManager * as parameter
* - the type of periodic boundary conditions that are used
- * by taking a PeriodicBoundaryConditionTypeEnum as parameter
+ * by taking a PeriodicBoundaryConditionType as parameter
+ * - the writing of checkpoint data
+ * by taking a MdModulesWriteCheckpointData as parameter
+ * - the reading of checkpoint data
+ * by taking a MdModulesCheckpointReadingDataOnMaster as parameter
+ * - the broadcasting of checkpoint data
+ * by taking MdModulesCheckpointReadingBroadcast as parameter
*/
explicit DensityFitting(MdModulesNotifier *notifier)
{
this->setEnergyOutputRequest(energyOutputRequest);
};
notifier->notifier_.subscribe(requestEnergyOutput);
+
+ // writing checkpoint data
+ const auto checkpointDataWriting = [this](MdModulesWriteCheckpointData checkpointData) {
+ this->writeCheckpointData(checkpointData);
+ };
+ notifier->notifier_.subscribe(checkpointDataWriting);
+
+ // reading checkpoint data
+ const auto checkpointDataReading = [this](MdModulesCheckpointReadingDataOnMaster checkpointData) {
+ this->readCheckpointDataOnMaster(checkpointData);
+ };
+ notifier->notifier_.subscribe(checkpointDataReading);
+
+ // broadcasting checkpoint data
+ const auto checkpointDataBroadcast = [this](MdModulesCheckpointReadingBroadcast checkpointData) {
+ this->broadcastCheckpointData(checkpointData);
+ };
+ notifier->notifier_.subscribe(checkpointDataBroadcast);
}
//! From IMDModule
densityFittingSimulationParameters_.referenceDensity(),
densityFittingSimulationParameters_.transformationToDensityLattice(),
densityFittingSimulationParameters_.localAtomSet(),
- densityFittingSimulationParameters_.periodicBoundaryConditionType());
+ densityFittingSimulationParameters_.periodicBoundaryConditionType(),
+ densityFittingState_);
forceProviders->addForceProvider(forceProvider_.get());
}
}
energyOutputRequest->energyOutputToDensityFitting_ = densityFittingOptions_.active();
}
+ /*! \brief Write internal density fitting data to checkpoint file.
+ * \param[in] checkpointWriting enables writing to the Key-Value-Tree
+ * that is used for storing the checkpoint
+ * information
+ */
+ void writeCheckpointData(MdModulesWriteCheckpointData checkpointWriting)
+ {
+ if (densityFittingOptions_.active())
+ {
+ const DensityFittingForceProviderState &state =
+ forceProvider_->state();
+ checkpointWriting.builder_.addValue<std::int64_t>(
+ DensityFittingModuleInfo::name_ + "-stepsSinceLastCalculation",
+ state.stepsSinceLastCalculation_);
+ }
+ }
+
+ /*! \brief Read the internal parameters from the checkpoint file on master
+ * \param[in] checkpointReading holding the checkpoint information
+ */
+ void readCheckpointDataOnMaster(MdModulesCheckpointReadingDataOnMaster checkpointReading)
+ {
+ if (densityFittingOptions_.active())
+ {
+ if (checkpointReading.checkpointedData_.keyExists(
+ DensityFittingModuleInfo::name_ + "-stepsSinceLastCalculation"))
+ {
+ densityFittingState_.stepsSinceLastCalculation_ =
+ checkpointReading.checkpointedData_[
+ DensityFittingModuleInfo::name_ + "-stepsSinceLastCalculation"].cast<std::int64_t>();
+ }
+ }
+ }
+
+ /*! \brief Broadcast the internal parameters.
+ * \param[in] checkpointBroadcast containing the communication record to
+ * broadcast the checkpoint information
+ */
+ void broadcastCheckpointData(MdModulesCheckpointReadingBroadcast checkpointBroadcast)
+ {
+ if (densityFittingOptions_.active())
+ {
+ if (PAR(&(checkpointBroadcast.cr_)))
+ {
+ block_bc(&(checkpointBroadcast.cr_), densityFittingState_.stepsSinceLastCalculation_);
+ }
+ }
+ }
+
private:
//! The output provider
DensityFittingOutputProvider densityFittingOutputProvider_;
/*! \brief Parameters for density fitting that become available at
* simulation setup time.
*/
- DensityFittingSimulationParameterSetup densityFittingSimulationParameters_;
+ DensityFittingSimulationParameterSetup densityFittingSimulationParameters_;
+ //! The internal parameters of density fitting force provider
+ DensityFittingForceProviderState densityFittingState_;
GMX_DISALLOW_COPY_AND_ASSIGN(DensityFitting);
};
basic_mdspan<const float, dynamicExtents3D> referenceDensity,
const TranslateAndScale &transformationToDensityLattice,
const LocalAtomSet &localAtomSet,
- int pbcType);
+ int pbcType,
+ const DensityFittingForceProviderState &state);
~Impl();
void calculateForces(const ForceProviderInput &forceProviderInput, ForceProviderOutput *forceProviderOutput);
+
+ DensityFittingForceProviderState state();
+
private:
const DensityFittingParameters ¶meters_;
+ DensityFittingForceProviderState state_;
LocalAtomSet localAtomSet_;
GaussianSpreadKernelParameters::Shape spreadKernel_;
TranslateAndScale transformationToDensityLattice_;
RVec referenceDensityCenter_;
int pbcType_;
+ real forceConstantScale_;
+
};
DensityFittingForceProvider::Impl::~Impl() = default;
basic_mdspan<const float, dynamicExtents3D> referenceDensity,
const TranslateAndScale &transformationToDensityLattice,
const LocalAtomSet &localAtomSet,
- int pbcType) :
+ int pbcType,
+ const DensityFittingForceProviderState &state) :
parameters_(parameters),
+ state_(state),
localAtomSet_(localAtomSet),
spreadKernel_(makeSpreadKernel(parameters_.gaussianTransformSpreadingWidth_,
parameters_.gaussianTransformSpreadingRangeInMultiplesOfWidth_,
transformedCoordinates_(localAtomSet_.numAtomsLocal()),
amplitudeLookup_(parameters_.amplitudeLookupMethod_),
transformationToDensityLattice_(transformationToDensityLattice),
- pbcType_(pbcType)
+ pbcType_(pbcType),
+ forceConstantScale_(parameters_.calculationIntervalInSteps_)
{
referenceDensityCenter_ = {
real(referenceDensity.extent(XX))/2,
void DensityFittingForceProvider::Impl::calculateForces(const ForceProviderInput &forceProviderInput,
ForceProviderOutput *forceProviderOutput)
{
+ // do nothing but count number of steps when not in density fitting step
+ if (state_.stepsSinceLastCalculation_ % parameters_.calculationIntervalInSteps_ != 0)
+ {
+ ++(state_.stepsSinceLastCalculation_);
+ return;
+ }
+
+ state_.stepsSinceLastCalculation_ = 0;
+
// do nothing if there are no density fitting atoms on this node
if (localAtomSet_.numAtomsLocal() == 0)
{
for (const auto localAtomIndex : localAtomSet_.localIndex())
{
forceProviderOutput->forceWithVirial_.force_[localAtomIndex] +=
- parameters_.forceConstant_ * *densityForceIterator;
+ forceConstantScale_ * parameters_.forceConstant_ * *densityForceIterator;
++densityForceIterator;
}
forceProviderOutput->enerd_.term[F_DENSITYFITTING] += energy;
}
+DensityFittingForceProviderState
+DensityFittingForceProvider::Impl::state()
+{
+ return state_;
+}
+
/********************************************************************
* DensityFittingForceProvider
*/
basic_mdspan<const float, dynamicExtents3D> referenceDensity,
const TranslateAndScale &transformationToDensityLattice,
const LocalAtomSet &localAtomSet,
- int pbcType)
- : impl_(new Impl(parameters, referenceDensity, transformationToDensityLattice, localAtomSet, pbcType))
+ int pbcType,
+ const DensityFittingForceProviderState &state)
+ : impl_(new Impl(parameters, referenceDensity, transformationToDensityLattice, localAtomSet, pbcType, state))
{}
void DensityFittingForceProvider::calculateForces(const ForceProviderInput &forceProviderInput,
impl_->calculateForces(forceProviderInput, forceProviderOutput);
}
+DensityFittingForceProviderState DensityFittingForceProvider::state()
+{
+ return impl_->state();
+}
+
} // namespace gmx
struct DensityFittingParameters;
+/*! \internal
+ * \brief Parameters defining the internal density fitting force provider state.
+ */
+struct DensityFittingForceProviderState
+{
+ /*! \brief The steps since the last force calculation.
+ * Used if density fitting is to be calculated every N steps.
+ */
+ std::int64_t stepsSinceLastCalculation_ = 0;
+};
+
/*! \internal \brief
* Implements IForceProvider for density-fitting forces.
*/
basic_mdspan<const float, dynamicExtents3D> referenceDensity,
const TranslateAndScale &transformationToDensityLattice,
const LocalAtomSet &localAtomSet,
- int pbcType);
+ int pbcType,
+ const DensityFittingForceProviderState &state);
~DensityFittingForceProvider();
- /* \brief Calculate forces that maximise goodness-of-fit with a reference density map
+ /*!\brief Calculate forces that maximise goodness-of-fit with a reference density map.
* \param[in] forceProviderInput input for force provider
* \param[out] forceProviderOutput output for force provider
*/
void calculateForces(const ForceProviderInput &forceProviderInput,
ForceProviderOutput *forceProviderOutput) override;
+
+ //! Return the state of the forceprovider.
+ DensityFittingForceProviderState state();
+
private:
class Impl;
PrivateImplPointer<Impl> impl_;
densityfittingMdpTransformFromString<real>(rules, &fromStdString<real>, c_gaussianTransformSpreadingWidthTag_);
densityfittingMdpTransformFromString<real>(rules, &fromStdString<real>, c_gaussianTransformSpreadingRangeInMultiplesOfWidthTag_);
densityfittingMdpTransformFromString<std::string>(rules, stringIdentityTransform, c_referenceDensityFileNameTag_);
+ densityfittingMdpTransformFromString<std::int64_t>(rules, &fromStdString<std::int64_t>, c_everyNStepsTag_);
}
void DensityFittingOptions::buildMdpOutput(KeyValueTreeObjectBuilder *builder) const
addDensityFittingMdpOutputValue(builder, parameters_.gaussianTransformSpreadingWidth_, c_gaussianTransformSpreadingWidthTag_);
addDensityFittingMdpOutputValue(builder, parameters_.gaussianTransformSpreadingRangeInMultiplesOfWidth_, c_gaussianTransformSpreadingRangeInMultiplesOfWidthTag_);
addDensityFittingMdpOutputValue(builder, referenceDensityFileName_, c_referenceDensityFileNameTag_);
+ addDensityFittingMdpOutputValue(builder, parameters_.calculationIntervalInSteps_, c_everyNStepsTag_);
}
}
section.addOption(RealOption(c_gaussianTransformSpreadingWidthTag_.c_str()).store(¶meters_.gaussianTransformSpreadingWidth_));
section.addOption(RealOption(c_gaussianTransformSpreadingRangeInMultiplesOfWidthTag_.c_str()).store(¶meters_.gaussianTransformSpreadingRangeInMultiplesOfWidth_));
section.addOption(StringOption(c_referenceDensityFileNameTag_.c_str()).store(&referenceDensityFileName_));
+ section.addOption(Int64Option(c_everyNStepsTag_.c_str()).store(¶meters_.calculationIntervalInSteps_));
}
bool DensityFittingOptions::active() const
const DensityFittingParameters &DensityFittingOptions::buildParameters()
{
+ // the options modules does not know unsigned integers so any input of this
+ // kind is rectified here
+ if (parameters_.calculationIntervalInSteps_ < 1)
+ {
+ parameters_.calculationIntervalInSteps_ = 1;
+ }
return parameters_;
}
const std::string c_referenceDensityFileNameTag_ = "reference-density-filename";
std::string referenceDensityFileName_ = "reference.mrc";
+ const std::string c_everyNStepsTag_ = "nst";
+
DensityFittingParameters parameters_;
};
real gaussianTransformSpreadingWidth_ = 0.2;
//! The spreading range for spreading atoms onto the grid in multiples of the spreading width
real gaussianTransformSpreadingRangeInMultiplesOfWidth_ = 4.0;
+ //! Apply density fitting forces only every n-steps
+ std::int64_t calculationIntervalInSteps_ = 1;
};
/*!\brief Check if two structs holding density fitting parameters are equal.
densityFittingModule_->initForceProviders(&densityFittingForces_);
}
- private:
+ protected:
KeyValueTreeBuilder mdpValueBuilder_;
MdModulesNotifier notifier_;
EXPECT_ANY_THROW(intializeForceProviders());
}
-
} // namespace
} // namespace gmx
"density-guided-simulation-gaussian-transform-spreading-width = 0.2\n"
"density-guided-simulation-gaussian-transform-spreading-range-in-multiples-of-width = 4\n"
"density-guided-simulation-reference-density-filename = reference.mrc\n"
+ "density-guided-simulation-nst = 1\n"
};
EXPECT_EQ(expected, stream.toString());
pullHistory //!< Pull history statistics (sums since last written output)
};
-namespace gmx
-{
-
-struct MdModulesCheckpointReadingDataOnMaster
-{
- //! The data of the MdModules that is stored in the checkpoint file
- const KeyValueTreeObject &checkpointedData;
- //! The version of the read ceckpoint file
- int checkpointFileVersion_;
-};
-
-/*! \libinternal
- * \brief Provides the MdModules with the communication record to broadcast.
- */
-struct MdModulesCheckpointReadingBroadcast
-{
- //! The communication record
- const t_commrec &cr;
- //! The version of the read file version
- int checkpointFileVersion_;
-};
-
-/*! \libinternal \brief Writing the MdModules data to a checkpoint file.
- */
-struct MdModulesWriteCheckpointData
-{
- //! Builder for the Key-Value-Tree to store the MdModule checkpoint data
- KeyValueTreeObjectBuilder builder;
- //! The version of the read file version
- int checkpointFileVersion_;
-};
-} // namespace gmx
-
//! \brief Return the name of a checkpoint entry based on part and part entry
static const char *entryName(StatePart part, int ecpt)
{
#include "gromacs/math/vectypes.h"
#include "gromacs/utility/basedefinitions.h"
+#include "gromacs/utility/keyvaluetreebuilder.h"
class energyhistory_t;
struct gmx_file_position_t;
struct t_inputrec;
class t_state;
struct t_trxframe;
+
namespace gmx
{
+
struct MdModulesNotifier;
-}
+class KeyValueTreeObject;
+
+struct MdModulesCheckpointReadingDataOnMaster
+{
+ //! The data of the MdModules that is stored in the checkpoint file
+ const KeyValueTreeObject &checkpointedData_;
+ //! The version of the read ceckpoint file
+ int checkpointFileVersion_;
+};
+
+/*! \libinternal
+ * \brief Provides the MdModules with the communication record to broadcast.
+ */
+struct MdModulesCheckpointReadingBroadcast
+{
+ //! The communication record
+ const t_commrec &cr_;
+ //! The version of the read file version
+ int checkpointFileVersion_;
+};
+
+/*! \libinternal \brief Writing the MdModules data to a checkpoint file.
+ */
+struct MdModulesWriteCheckpointData
+{
+ //! Builder for the Key-Value-Tree to store the MdModule checkpoint data
+ KeyValueTreeObjectBuilder builder_;
+ //! The version of the read file version
+ int checkpointFileVersion_;
+};
+
+} // namespace gmx
/* the name of the environment variable to disable fsync failure checks with */
#define GMX_IGNORE_FSYNC_FAILURE_ENV "GMX_IGNORE_FSYNC_FAILURE"