Gromacs  2024.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Typedefs | Enumerations | Functions | Variables
The modular simulator

Description

The modular simulator improves extensibility, adds Monte Carlo capabilities, promotes data locality and communication via interfaces, supports multi-stepping integrators, and paves the way for some task parallelism.

For more information, see page_modularsimulator

Todo:
Can we link to docs/doxygen/lib/modularsimulator.md?
Author
Pascal Merz pasca.nosp@m.l.me.nosp@m.rz@me.nosp@m..com

Classes

class  gmx::ElementNotFoundError
 Exception class signalling that a requested element was not found. More...
 
class  gmx::MissingElementConnectionError
 Exception class signalling that elements were not connected properly. More...
 
class  gmx::SimulationAlgorithmSetupError
 Exception class signalling that the ModularSimulatorAlgorithm was set up in an incompatible way. More...
 
class  gmx::CheckpointError
 Exception class signalling an error in reading or writing modular checkpoints. More...
 

Typedefs

typedef std::function< void()> gmx::CheckBondedInteractionsCallback
 The function type allowing to request a check of the number of bonded interactions.
 
using gmx::Step = int64_t
 Step number.
 
using gmx::Time = double
 Simulation time.
 
typedef std::function< void()> gmx::SimulatorRunFunction
 The function type that can be scheduled to be run during the simulator run.
 
typedef std::function< void(SimulatorRunFunction)> gmx::RegisterRunFunction
 The function type that allows to register run functions.
 
typedef std::function< void(Step,
Time, const
RegisterRunFunction &)> 
gmx::SchedulingFunction
 The function type scheduling run functions for a step / time using a RegisterRunFunction reference.
 
typedef std::function< void(Step,
Time)> 
gmx::SignallerCallback
 The function type that can be registered to signallers for callback.
 
typedef std::function< void(gmx_mdoutf
*, Step, Time, bool, bool)> 
gmx::ITrajectoryWriterCallback
 Function type for trajectory writing clients.
 
typedef std::function< void(Step)> gmx::PropagatorCallback
 Generic callback to the propagator.
 
typedef std::function< void()> gmx::DomDecCallback
 Callback used by the DomDecHelper object to inform clients about system re-partitioning.
 
using gmx::ReferenceTemperatureCallback = std::function< void(ArrayRef< const real >, ReferenceTemperatureChangeAlgorithm algorithm)>
 Callback updating the reference temperature.
 

Enumerations

enum  gmx::ComputeGlobalsAlgorithm { LeapFrog, VelocityVerlet }
 The different global reduction schemes we know about.
 
enum  gmx::EnergySignallerEvent { EnergyCalculationStep, VirialCalculationStep, FreeEnergyCalculationStep }
 The energy events signalled by the EnergySignaller.
 
enum  gmx::TrajectoryEvent { StateWritingStep, EnergyWritingStep }
 The trajectory writing events.
 
enum  gmx::ModularSimulatorBuilderState { AcceptingClientRegistrations, NotAcceptingClientRegistrations }
 Enum allowing builders to store whether they can accept client registrations.
 
enum  gmx::ReportPreviousStepConservedEnergy { Yes, No, Count }
 Enum describing whether an element is reporting conserved energy from the previous step.
 
enum  gmx::ScaleVelocities { PreStepOnly, PreStepAndPostStep }
 Which velocities the thermostat scales.
 
enum  gmx::IntegrationStage {
  gmx::IntegrationStage::PositionsOnly, gmx::IntegrationStage::VelocitiesOnly, gmx::IntegrationStage::LeapFrog, gmx::IntegrationStage::VelocityVerletPositionsAndVelocities,
  gmx::IntegrationStage::ScaleVelocities, gmx::IntegrationStage::ScalePositions, gmx::IntegrationStage::Count
}
 The different integration types we know about. More...
 
enum  gmx::NumPositionScalingValues { gmx::NumPositionScalingValues::None, gmx::NumPositionScalingValues::Single, gmx::NumPositionScalingValues::Multiple, gmx::NumPositionScalingValues::Count }
 Sets the number of different position scaling values. More...
 
enum  gmx::NumVelocityScalingValues { gmx::NumVelocityScalingValues::None, gmx::NumVelocityScalingValues::Single, gmx::NumVelocityScalingValues::Multiple, Count }
 Sets the number of different velocity scaling values. More...
 
enum  gmx::ParrinelloRahmanVelocityScaling { gmx::ParrinelloRahmanVelocityScaling::No, gmx::ParrinelloRahmanVelocityScaling::Diagonal, gmx::ParrinelloRahmanVelocityScaling::Anisotropic, Count }
 Describes the properties of the Parrinello-Rahman pressure scaling matrix. More...
 

Functions

template<typename... Ts>
auto gmx::checkUseModularSimulator (Ts &&...args) -> decltype(ModularSimulator::isInputCompatible(std::forward< Ts >(args)...))
 Whether or not to use the ModularSimulator. More...
 
virtual void gmx::ISimulatorElement::scheduleTask (Step, Time, const RegisterRunFunction &)=0
 Query whether element wants to run at step / time. More...
 
virtual void gmx::ISimulatorElement::elementSetup ()=0
 Method guaranteed to be called after construction, before simulator run.
 
virtual void gmx::ISimulatorElement::elementTeardown ()=0
 Method guaranteed to be called after simulator run, before deconstruction.
 
virtual gmx::ISimulatorElement::~ISimulatorElement ()=default
 Standard virtual destructor.
 
virtual void gmx::ISignaller::signal (Step, Time)=0
 Function run before every step of scheduling.
 
virtual void gmx::ISignaller::setup ()=0
 Method guaranteed to be called after construction, before simulator run.
 
virtual gmx::ISignaller::~ISignaller ()=default
 Standard virtual destructor.
 
virtual gmx::INeighborSearchSignallerClient::~INeighborSearchSignallerClient ()=default
 Standard virtual destructor.
 
virtual std::optional
< SignallerCallback > 
gmx::INeighborSearchSignallerClient::registerNSCallback ()=0
 Return callback to NeighborSearchSignaller.
 
virtual gmx::ILastStepSignallerClient::~ILastStepSignallerClient ()=default
 Standard virtual destructor.
 
virtual std::optional
< SignallerCallback > 
gmx::ILastStepSignallerClient::registerLastStepCallback ()=0
 Return callback to LastStepSignaller.
 
virtual gmx::ILoggingSignallerClient::~ILoggingSignallerClient ()=default
 Standard virtual destructor.
 
virtual std::optional
< SignallerCallback > 
gmx::ILoggingSignallerClient::registerLoggingCallback ()=0
 Return callback to LoggingSignaller.
 
virtual gmx::IEnergySignallerClient::~IEnergySignallerClient ()=default
 Standard virtual destructor.
 
virtual std::optional
< SignallerCallback > 
gmx::IEnergySignallerClient::registerEnergyCallback (EnergySignallerEvent)=0
 Return callback to EnergySignaller.
 
virtual gmx::ITrajectorySignallerClient::~ITrajectorySignallerClient ()=default
 Standard virtual destructor.
 
virtual std::optional
< SignallerCallback > 
gmx::ITrajectorySignallerClient::registerTrajectorySignallerCallback (TrajectoryEvent)=0
 Return callback to TrajectoryElement.
 
virtual gmx::ITrajectoryWriterClient::~ITrajectoryWriterClient ()=default
 Standard virtual destructor.
 
virtual void gmx::ITrajectoryWriterClient::trajectoryWriterSetup (gmx_mdoutf *outf)=0
 Setup method with valid output pointer.
 
virtual void gmx::ITrajectoryWriterClient::trajectoryWriterTeardown (gmx_mdoutf *outf)=0
 Teardown method with valid output pointer.
 
virtual std::optional
< ITrajectoryWriterCallback > 
gmx::ITrajectoryWriterClient::registerTrajectoryWriterCallback (TrajectoryEvent)=0
 Return callback to TrajectoryElement.
 
virtual gmx::ITopologyHolderClient::~ITopologyHolderClient ()=default
 Standard virtual destructor.
 
virtual void gmx::ITopologyHolderClient::setTopology (const gmx_localtop_t *)=0
 Pass pointer to new local topology.
 
virtual gmx::ICheckpointHelperClient::~ICheckpointHelperClient ()=default
 Standard virtual destructor.
 
virtual void gmx::ICheckpointHelperClient::saveCheckpointState (std::optional< WriteCheckpointData > checkpointData, const t_commrec *cr)=0
 Write checkpoint (CheckpointData object only passed on main rank)
 
virtual void gmx::ICheckpointHelperClient::restoreCheckpointState (std::optional< ReadCheckpointData > checkpointData, const t_commrec *cr)=0
 Read checkpoint (CheckpointData object only passed on main rank)
 
virtual const std::string & gmx::ICheckpointHelperClient::clientID ()=0
 Get unique client id.
 
 gmx::ElementNotFoundError::ElementNotFoundError (const ExceptionInitializer &details)
 Creates an exception object with the provided initializer/reason. More...
 
 gmx::MissingElementConnectionError::MissingElementConnectionError (const ExceptionInitializer &details)
 Creates an exception object with the provided initializer/reason. More...
 
 gmx::SimulationAlgorithmSetupError::SimulationAlgorithmSetupError (const ExceptionInitializer &details)
 Creates an exception object with the provided initializer/reason. More...
 
 gmx::CheckpointError::CheckpointError (const ExceptionInitializer &details)
 Creates an exception object with the provided initializer/reason. More...
 
 gmx::PropagatorTag::PropagatorTag (std::string_view name)
 Explicit constructor.
 
 gmx::PropagatorTag::operator const std::string & () const
 Can be used as string.
 
bool gmx::PropagatorTag::operator== (const PropagatorTag &other) const
 Equality operator.
 
bool gmx::PropagatorTag::operator!= (const PropagatorTag &other) const
 Inequality operator.
 
 gmx::TimeStep::TimeStep (real timeStep)
 Explicit constructor.
 
 gmx::TimeStep::operator const real & () const
 Can be used as underlying type.
 
 gmx::Offset::Offset (int offset)
 Explicit constructor.
 
 gmx::Offset::operator const int & () const
 Can be used as underlying type.
 
bool gmx::PropagatorConnection::hasStartVelocityScaling () const
 Whether the propagator offers start velocity scaling.
 
bool gmx::PropagatorConnection::hasEndVelocityScaling () const
 Whether the propagator offers end velocity scaling.
 
bool gmx::PropagatorConnection::hasPositionScaling () const
 Whether the propagator offers position scaling.
 
bool gmx::PropagatorConnection::hasParrinelloRahmanScaling () const
 Whether the propagator offers Parrinello-Rahman scaling.
 
virtual gmx::IDomDecHelperClient::~IDomDecHelperClient ()=default
 Standard virtual destructor.
 
virtual DomDecCallback gmx::IDomDecHelperClient::registerDomDecCallback ()=0
 Register function to be informed about system re-partitioning.
 

Variables

PropagatorTag gmx::PropagatorConnection::tag
 The tag of the creating propagator.
 
std::function< void(int,
ScaleVelocities)> 
gmx::PropagatorConnection::setNumVelocityScalingVariables
 Function object for setting velocity scaling variables.
 
std::function< void(int)> gmx::PropagatorConnection::setNumPositionScalingVariables
 Function object for setting velocity scaling variables.
 
std::function< ArrayRef< real >)> gmx::PropagatorConnection::getViewOnStartVelocityScaling
 Function object for receiving view on velocity scaling (before step)
 
std::function< ArrayRef< real >)> gmx::PropagatorConnection::getViewOnEndVelocityScaling
 Function object for receiving view on velocity scaling (after step)
 
std::function< ArrayRef< real >)> gmx::PropagatorConnection::getViewOnPositionScaling
 Function object for receiving view on position scaling.
 
std::function
< PropagatorCallback()> 
gmx::PropagatorConnection::getVelocityScalingCallback
 Function object to request callback allowing to signal a velocity scaling step.
 
std::function
< PropagatorCallback()> 
gmx::PropagatorConnection::getPositionScalingCallback
 Function object to request callback allowing to signal a position scaling step.
 
std::function< Matrix3x3 *()> gmx::PropagatorConnection::getViewOnPRScalingMatrix
 Function object for receiving view on pressure scaling matrix.
 
std::function
< PropagatorCallback()> 
gmx::PropagatorConnection::getPRScalingCallback
 Function object to request callback allowing to signal a Parrinello-Rahman scaling step.
 

Enumeration Type Documentation

enum gmx::IntegrationStage
strong

The different integration types we know about.

Enumerator
PositionsOnly 

Moves the position vector by the given time step.

VelocitiesOnly 

Moves the velocity vector by the given time step.

LeapFrog 

Manual fusion of the previous two propagators.

VelocityVerletPositionsAndVelocities 

Manual position (full dt) and velocity (half dt) fusion.

ScaleVelocities 

Only scale velocities, don't propagate.

ScalePositions 

Only scale positions, don't propagate.

Count 

The number of enum entries.

Sets the number of different position scaling values.

Enumerator
None 

No position scaling (either this step or ever)

Single 

Single scaling value (either one group or all values =1)

Multiple 

Multiple scaling values, need to use T-group indices.

Count 

The number of enum entries.

Sets the number of different velocity scaling values.

Enumerator
None 

No velocity scaling (either this step or ever)

Single 

Single T-scaling value (either one group or all values =1)

Multiple 

Multiple T-scaling values, need to use T-group indices.

Describes the properties of the Parrinello-Rahman pressure scaling matrix.

Enumerator
No 

Do not apply velocity scaling (not a PR-coupling run or step)

Diagonal 

Apply velocity scaling using a diagonal matrix.

Anisotropic 

Apply velocity scaling using a matrix with off-diagonal elements.

Function Documentation

gmx::CheckpointError::CheckpointError ( const ExceptionInitializer details)
inlineexplicit

Creates an exception object with the provided initializer/reason.

Parameters
[in]detailsInitializer for the exception.
Exceptions
std::bad_allocif out of memory.

It is possible to call this constructor either with an explicit ExceptionInitializer object (useful for more complex cases), or a simple string if only a reason string needs to be provided.

template<typename... Ts>
auto gmx::checkUseModularSimulator ( Ts &&...  args) -> decltype(ModularSimulator::isInputCompatible(std::forward<Ts>(args)...))

Whether or not to use the ModularSimulator.

GMX_DISABLE_MODULAR_SIMULATOR environment variable allows to disable modular simulator for all uses.

See ModularSimulator::isInputCompatible() for function signature.

gmx::ElementNotFoundError::ElementNotFoundError ( const ExceptionInitializer details)
inlineexplicit

Creates an exception object with the provided initializer/reason.

Parameters
[in]detailsInitializer for the exception.
Exceptions
std::bad_allocif out of memory.

It is possible to call this constructor either with an explicit ExceptionInitializer object (useful for more complex cases), or a simple string if only a reason string needs to be provided.

gmx::MissingElementConnectionError::MissingElementConnectionError ( const ExceptionInitializer details)
inlineexplicit

Creates an exception object with the provided initializer/reason.

Parameters
[in]detailsInitializer for the exception.
Exceptions
std::bad_allocif out of memory.

It is possible to call this constructor either with an explicit ExceptionInitializer object (useful for more complex cases), or a simple string if only a reason string needs to be provided.

virtual void gmx::ISimulatorElement::scheduleTask ( Step  ,
Time  ,
const RegisterRunFunction  
)
pure virtual

Query whether element wants to run at step / time.

Element can register one or more functions to be run at that step through the registration pointer.

gmx::SimulationAlgorithmSetupError::SimulationAlgorithmSetupError ( const ExceptionInitializer details)
inlineexplicit

Creates an exception object with the provided initializer/reason.

Parameters
[in]detailsInitializer for the exception.
Exceptions
std::bad_allocif out of memory.

It is possible to call this constructor either with an explicit ExceptionInitializer object (useful for more complex cases), or a simple string if only a reason string needs to be provided.