md.methods#

Overview

Brownian

Brownian dynamics.

ConstantPressure

Constant pressure dynamics.

ConstantVolume

Constant volume dynamics.

DisplacementCapped

Newtonian dynamics with a cap on the maximum displacement per time step.

Langevin

Langevin dynamics.

Method

Base class integration method.

OverdampedViscous

Overdamped viscous dynamics.

Thermostatted

Base class for thermostatted integrators.

Details

Integration methods for molecular dynamics.

Integration methods work with hoomd.md.Integrator to define the equations of motion for the system. Each individual method applies the given equations of motion to a subset of particles.

Thermostatted methods

Thermostatted methods require usage of a thermostat, see hoomd.md.methods.thermostats.

Integration methods with constraints

For methods that constrain motion to a manifold see hoomd.md.methods.rattle.

class hoomd.md.methods.Brownian(filter, kT, default_gamma=1.0, default_gamma_r=(1.0, 1.0, 1.0))#

Bases: Method

Brownian dynamics.

Parameters:
  • filter (hoomd.filter.filter_like) – Subset of particles to apply this method to.

  • kT (hoomd.variant.variant_like) – Temperature of the simulation \([\mathrm{energy}]\).

  • default_gamma (float) – Default drag coefficient for all particle types \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\).

  • default_gamma_r ([float, float, float]) – Default rotational drag coefficient tensor for all particles \([\mathrm{time}^{-1}]\).

Brownian integrates particles forward in time according to the overdamped Langevin equations of motion, sometimes called Brownian dynamics or the diffusive limit.

The translational degrees of freedom follow:

\[ \begin{align}\begin{aligned}\frac{d\vec{r}}{dt} &= \frac{\vec{F}_\mathrm{C} + \vec{F}_\mathrm{R}}{\gamma},\\\langle \vec{F}_\mathrm{R} \rangle &= 0,\\\langle |\vec{F}_\mathrm{R}|^2 \rangle &= 2 d k T \gamma / \delta t,\\\langle \vec{v}(t) \rangle &= 0,\\\langle |\vec{v}(t)|^2 \rangle &= d k T / m,\end{aligned}\end{align} \]

where \(\vec{F}_\mathrm{C} = \vec{F}_\mathrm{net}\) is the net force on the particle from all forces (hoomd.md.Integrator.forces) and constraints (hoomd.md.Integrator.constraints), \(\gamma\) is the translational drag coefficient (gamma), \(\vec{F}_\mathrm{R}\) is a uniform random force, \(\vec{v}\) is the particle’s velocity, and \(d\) is the dimensionality of the system. The magnitude of the random force is chosen via the fluctuation-dissipation theorem to be consistent with the specified drag and temperature, \(kT\).

About axes where \(I^i > 0\), the rotational degrees of freedom follow:

\[ \begin{align}\begin{aligned}\frac{d\mathbf{q}}{dt} &= \frac{\vec{\tau}_\mathrm{C} + \vec{\tau}_\mathrm{R}}{\gamma_r},\\\langle \vec{\tau}_\mathrm{R} \rangle &= 0,\\\langle \tau_\mathrm{R}^i \cdot \tau_\mathrm{R}^i \rangle &= 2 k T \gamma_r^i / \delta t,\\\langle \vec{L}(t) \rangle &= 0,\\\langle L^i(t) \cdot L^i(t) \rangle &= k T \cdot I^i,\end{aligned}\end{align} \]

where \(\vec{\tau}_\mathrm{C} = \vec{\tau}_\mathrm{net}\), \(\gamma_r^i\) is the i-th component of the rotational drag coefficient (gamma_r), \(\tau_\mathrm{R}^i\) is a component of the uniform random the torque, \(L^i\) is the i-th component of the particle’s angular momentum and \(I^i\) is the i-th component of the particle’s moment of inertia. The magnitude of the random torque is chosen via the fluctuation-dissipation theorem to be consistent with the specified drag and temperature, \(kT\).

Brownian uses the numerical integration method from I. Snook 2007, The Langevin and Generalised Langevin Approach to the Dynamics of Atomic, Polymeric and Colloidal Systems, section 6.2.5, with the exception that \(\vec{F}_\mathrm{R}\) is drawn from a uniform random number distribution.

Warning

This numerical method has errors in \(O(\delta t)\), which is much larger than the errors of most other integration methods which are in \(O(\delta t^2)\). As a consequence, expect to use much smaller values of \(\delta t\) with Brownian compared to e.g. Langevin or ConstantVolume.

In Brownian dynamics, particle velocities and angular momenta are completely decoupled from positions. At each time step, Brownian draws a new velocity distribution consistent with the current set temperature so that hoomd.md.compute.ThermodynamicQuantities will report appropriate temperatures and pressures when logged or used by other methods.

The attributes gamma and gamma_r set the translational and rotational damping coefficients, respectivley, by particle type.

Example:

brownian = hoomd.md.methods.Brownian(filter=hoomd.filter.All(), kT=1.5)
simulation.operations.integrator.methods = [brownian]
filter#

Subset of particles to apply this method to.

Type:

hoomd.filter.filter_like

kT#

Temperature of the simulation \([\mathrm{energy}]\).

Examples:

brownian.kT = 1.0
brownian.kT = hoomd.variant.Ramp(A=2.0,
                                 B=1.0,
                                 t_start=0,
                                 t_ramp=1_000_000)
Type:

hoomd.variant.Variant

gamma#

The drag coefficient for each particle type \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\).

Example:

brownian.gamma['A'] = 0.5
Type:

TypeParameter[ particle type, float ]

gamma_r#

The rotational drag coefficient tensor for each particle type \([\mathrm{time}^{-1}]\).

Example:

brownian.gamma_r['A'] = [1.0, 2.0, 3.0]
Type:

TypeParameter[particle type,[float, float , float]]

class hoomd.md.methods.ConstantPressure(filter, S, tauS, couple, thermostat=None, box_dof=[True, True, True, False, False, False], rescale_all=False, gamma=0.0)#

Bases: Thermostatted

Constant pressure dynamics.

Parameters:
  • filter (hoomd.filter.filter_like) – Subset of particles on which to apply this method.

  • thermostat (hoomd.md.methods.thermostats.Thermostat) – Thermostat to control temperature. Setting this to None samples a constant enthalpy (NPH) integration.

  • S (tuple[variant.variant_like, ...] or variant.variant_like) –

    Stress component set points for the barostat.

    In Voigt notation: \([S_{xx}, S_{yy}, S_{zz}, S_{yz}, S_{xz}, S_{xy}]\) \([\mathrm{pressure}]\). In case of isotropic pressure P use S = p to imply (\([p, p, p, 0, 0, 0]\)).

  • tauS (float) – Coupling constant for the barostat \([\mathrm{time}]\).

  • couple (str) – Couplings of diagonal elements of the stress tensor. One of “none”, “xy”, “xz”,”yz”, or “xyz”.

  • box_dof (list [ bool ]) – Box degrees of freedom with six boolean elements in the order x, y, z, xy, xz, yz. Defaults to [True,True,True,False,False,False]). When True, rescale corresponding lengths or tilt factors and components of particle coordinates and velocities.

  • rescale_all (bool) – When True, rescale all particles, not only those selected by the filter. Defaults to False.

  • gamma (float) – Friction constant for the box degrees of freedom. Defaults to 0 \([\mathrm{time}^{-1}]\).

ConstantPressure integrates translational and rotational degrees of freedom of the system held at constant pressure with a barostat. The barostat introduces additional degrees of freedom in the Hamiltonian that couple with box parameters. Use a thermostat to model an isothermal-isobaric (NPT) ensemble. Use no thermostat (thermostat = None) to model a isoenthalpic-isobaric (NPH) ensemble.

The barostat tensor is \(\nu_{\mathrm{ij}}\). Access these quantities using barostat_dof.

By default, ConstantPressure performs integration in a cubic box under hydrostatic pressure by simultaneously rescaling the lengths Lx, Ly and Lz of the simulation box by the same factors. Set the couplings and/or box degrees of freedom to change this default.

Couplings define which diagonal elements of the pressure tensor \(P_{\alpha,\beta}\) should be averaged over, so that the corresponding box lengths are rescaled by the same amount.

Valid couplings are:

  • 'none' (all box lengths are updated independently)

  • 'xy' (Lx and Ly are coupled)

  • 'xz' (Lx and Lz are coupled)

  • 'yz' (Ly and Lz are coupled)

  • 'xyz' (Lx, Ly, and Lz are coupled)

The degrees of freedom of the box set which lengths and tilt factors of the box should be updated, and how particle coordinates and velocities should be rescaled. The box_dof tuple controls the way the box is rescaled and updated. The first three elements box_dof[:3] controls whether the x, y, and z box lengths are rescaled and updated, respectively. The last three entries box_dof[3:] control the rescaling or the tilt factors xy, xz, and yz. All options also appropriately rescale particle coordinates and velocities.

By default, the x, y, and z degrees of freedom are updated. [True,True,True,False,False,False]

Note

When any of the diagonal x, y, z degrees of freedom is not being integrated, pressure tensor components along that direction are not considered for the remaining degrees of freedom.

ConstantPressure numerically integrates the equations of motion using the symplectic Martyna-Tobias-Klein integrator with a Langevin piston. The equation of motion of box dimensions is given by:

\[ \begin{align}\begin{aligned}\frac{d^2 L}{dt^2} &= V W^{-1} (S - S_{ext}) - \gamma \frac{dL}{dt} + R(t)\\\langle R \rangle &= 0\\\langle |R|^2 \rangle &= 2 \gamma kT \delta t W^{-1}\end{aligned}\end{align} \]

Where \(\gamma\) is the friction on the barostat piston, which damps unphysical volume oscillations at the cost of non-deterministic integration, and \(R\) is a random force, chosen appropriately for the coupled degrees of freedom.

Note

The barostat coupling constant tauS should be set within a reasonable range to avoid abrupt fluctuations in the box volume and to avoid long time to equilibration. The recommended value for most systems is \(\tau_S = 1000 \delta t\).

Note

If \(\gamma\) is used, its value should be chosen so that the system is near critical damping. A good initial guess is \(\gamma \approx 2 \tau_S^{-1}\). A value too high will result in long relaxation times.

Note

Set gamma = 0 to obtain the same MTK equations of motion used in HOOMD-blue releases prior to 4.0.0.

Examples:

NPH integrator with cubic symmetry:

nph = hoomd.md.methods.ConstantPressure(filter=hoomd.filter.All(),
                                        tauS=1.0,
                                        S=2.0,
                                        couple="xyz")
simulation.operations.integrator.methods = [nph]

NPT integrator with cubic symmetry:

npt = hoomd.md.methods.ConstantPressure(
    filter=hoomd.filter.All(),
    tauS=1.0,
    S=2.0,
    couple="xyz",
    thermostat=hoomd.md.methods.thermostats.Bussi(kT=1.5))
simulation.operations.integrator.methods = [npt]

NPT integrator with tetragonal symmetry:

npt = hoomd.md.methods.ConstantPressure(
    filter=hoomd.filter.All(),
    tauS = 1.0,
    S=2.0,
    couple="xy",
    thermostat=hoomd.md.methods.thermostats.Bussi(kT=1.5))
simulation.operations.integrator.methods = [npt]

NPT integrator with orthorhombic symmetry:

npt = hoomd.md.methods.ConstantPressure(
    filter=hoomd.filter.All(),
    tauS = 1.0,
    S=2.0,
    couple="none",
    thermostat=hoomd.md.methods.thermostats.Bussi(kT=1.5))
simulation.operations.integrator.methods = [npt]

NPT integrator with triclinic symmetry:

npt = hoomd.md.methods.ConstantPressure(
    filter=hoomd.filter.All(),
    tauS = 1.0,
    S=2.0,
    couple="none",
    box_dof=[True, True, True, True, True, True],
    thermostat=hoomd.md.methods.thermostats.Bussi(kT=1.5))
simulation.operations.integrator.methods = [npt]
filter#

Subset of particles on which to apply this method.

Type:

hoomd.filter.filter_like

thermostat#

Temperature control for the integrator.

Type:

hoomd.md.methods.thermostats.Thermostat

S#

Stress components set point for the barostat. In Voigt notation, \([S_{xx}, S_{yy}, S_{zz}, S_{yz}, S_{xz}, S_{xy}]\) \([\mathrm{pressure}]\).

Examples:

npt.S = 4.0
npt.S = hoomd.variant.Ramp(A=1.0,
                           B=2.0,
                           t_start=0,
                           t_ramp=1_000_000)
Type:

tuple[hoomd.variant.Variant,…]

tauS#

Coupling constant for the barostat \([\mathrm{time}]\).

Example:

npt.tauS = 2.0
Type:

float

couple#

Couplings of diagonal elements of the stress tensor, can be ‘none’, ‘xy’, ‘xz’, ‘yz’, or ‘xyz’.

Example:

npt.couple = 'none'
Type:

str

box_dof#

Box degrees of freedom with six boolean elements in the order [x, y, z, xy, xz, yz].

Example:

npt.box_dof = [False, False, True, False, False, False]
Type:

list[bool]

rescale_all#

When True, rescale all particles, not only those selected by the filter.

Example:

npt.rescale_all = True
Type:

bool

gamma#

Friction constant for the box degrees of freedom \([\mathrm{time^{-1}}]\).

Type:

float

barostat_dof#

Additional degrees of freedom for the barostat (\(\nu_{xx}\), \(\nu_{xy}\), \(\nu_{xz}\), \(\nu_{yy}\), \(\nu_{yz}\), \(\nu_{zz}\))

Save and restore the barostat degrees of freedom when continuing simulations:

Examples:

Save before exiting:

numpy.save(file=path / 'barostat_dof.npy',
           arr=npt.barostat_dof)

Load when continuing:

npt = hoomd.md.methods.ConstantPressure(
    filter=hoomd.filter.All(),
    tauS=1.0,
    S=2.0,
    couple="xyz",
    thermostat=hoomd.md.methods.thermostats.Bussi(kT=1.5))
simulation.operations.integrator.methods = [npt]

npt.barostat_dof = numpy.load(file=path / 'barostat_dof.npy')
Type:

tuple[float, float, float, float, float, float]

property barostat_energy#

Energy the barostat contributes to the Hamiltonian \([\mathrm{energy}]\).

Example:

logger.add(obj=npt, quantities=['barostat_energy'])

(Loggable: category=”scalar”)

thermalize_barostat_dof()#

Set the thermostat and barostat momenta to random values.

thermalize_barostat_dof sets random values for the the barostat momentum \(\nu_{\mathrm{ij}}\).

Important

You must call Simulation.run before thermalize_barostat_dof.

simulation.run(0)
npt.thermalize_barostat_dof()
class hoomd.md.methods.ConstantVolume(filter, thermostat=None)#

Bases: Thermostatted

Constant volume dynamics.

Parameters:

ConstantVolume numerically integrates the translational degrees of freedom using Velocity-Verlet and the rotational degrees of freedom with a scheme based on Kamberaj 2005.

When set, the thermostat rescales the particle velocities to model a canonical (NVT) ensemble. Use no thermostat (thermostat = None) to perform constant energy integration.

Examples:

NVE integration:

nve = hoomd.md.methods.ConstantVolume(filter=hoomd.filter.All())
simulation.operations.integrator.methods = [nve]

NVT integration:

nvt = hoomd.md.methods.ConstantVolume(
    filter=hoomd.filter.All(),
    thermostat=hoomd.md.methods.thermostats.Bussi(kT=1.5))
simulation.operations.integrator.methods = [nvt]
filter#

Subset of particles on which to apply this method.

Type:

hoomd.filter.filter_like

thermostat#

Temperature control for the integrator.

Examples:

nvt.thermostat.kT = 1.0
nvt.thermostat = hoomd.md.methods.thermostats.Bussi(kT=0.5)
Type:

hoomd.md.methods.thermostats.Thermostat

class hoomd.md.methods.DisplacementCapped(filter, maximum_displacement: Variant | float)#

Bases: ConstantVolume

Newtonian dynamics with a cap on the maximum displacement per time step.

The method limits particle motion to a maximum displacement allowed each time step which may be helpful to relax a high energy initial condition.

Warning

This method does not conserve energy or momentum.

Parameters:

DisplacementCapped integrates integrates translational and rotational degrees of freedom using modified microcanoncial dynamics. See NVE for the basis of the algorithm.

Example:

displacement_capped = hoomd.md.methods.DisplacementCapped(
    filter=hoomd.filter.All(),
    maximum_displacement=1e-3)
simulation.operations.integrator.methods = [displacement_capped]
filter#

Subset of particles on which to apply this method.

Type:

hoomd.filter.filter_like

maximum_displacement#

The maximum displacement allowed for a particular timestep \([\mathrm{length}]\).

displacement_capped.maximum_displacement = 1e-5
Type:

hoomd.variant.variant_like

class hoomd.md.methods.Langevin(filter, kT, tally_reservoir_energy=False, default_gamma=1.0, default_gamma_r=(1.0, 1.0, 1.0))#

Bases: Method

Langevin dynamics.

Parameters:
  • filter (hoomd.filter.filter_like) – Subset of particles to apply this method to.

  • kT (hoomd.variant.variant_like) – Temperature of the simulation \([\mathrm{energy}]\).

  • tally_reservoir_energy (bool) – When True, track the energy exchange between the thermal reservoir and the particles. Defaults to False \([\mathrm{energy}]\).

  • default_gamma (float) – Default drag coefficient for all particle types \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\).

  • default_gamma_r ([float, float, float]) – Default rotational drag coefficient tensor for all particles \([\mathrm{time}^{-1}]\).

Langevin integrates particles forward in time according to the Langevin equations of motion, modelling a canonical ensemble (NVT).

The translational degrees of freedom follow:

\[ \begin{align}\begin{aligned}m \frac{d\vec{v}}{dt} &= \vec{F}_\mathrm{C} - \gamma \cdot \vec{v} + \vec{F}_\mathrm{R}\\\langle \vec{F}_\mathrm{R} \rangle &= 0\\\langle |\vec{F}_\mathrm{R}|^2 \rangle &= 2 d kT \gamma / \delta t\end{aligned}\end{align} \]

where \(\vec{F}_\mathrm{C}\) is the force on the particle from all potentials and constraint forces, \(\gamma\) is the drag coefficient, \(\vec{v}\) is the particle’s velocity, \(\vec{F}_\mathrm{R}\) is a uniform random force, and \(d\) is the dimensionality of the system (2 or 3). The magnitude of the random force is chosen via the fluctuation-dissipation theorem to be consistent with the specified drag and temperature, \(kT\).

About axes where \(I^i > 0\), the rotational degrees of freedom follow:

\[ \begin{align}\begin{aligned}I \frac{d\vec{L}}{dt} &= \vec{\tau}_\mathrm{C} - \gamma_r \cdot \vec{L} + \vec{\tau}_\mathrm{R}\\\langle \vec{\tau}_\mathrm{R} \rangle &= 0,\\\langle \tau_\mathrm{R}^i \cdot \tau_\mathrm{R}^i \rangle &= 2 k T \gamma_r^i / \delta t,\end{aligned}\end{align} \]

where \(\vec{\tau}_\mathrm{C} = \vec{\tau}_\mathrm{net}\), \(\gamma_r^i\) is the i-th component of the rotational drag coefficient (gamma_r), \(\tau_\mathrm{R}^i\) is a component of the uniform random the torque, \(\vec{L}\) is the particle’s angular momentum and \(I\) is the the particle’s moment of inertia. The magnitude of the random torque is chosen via the fluctuation-dissipation theorem to be consistent with the specified drag and temperature, \(kT\).

Langevin numerically integrates the translational degrees of freedom using Velocity-Verlet and the rotational degrees of freedom with a scheme based on Kamberaj 2005.

The attributes gamma and gamma_r set the translational and rotational damping coefficients, respectivley, by particle type.

Example:

langevin = hoomd.md.methods.Langevin(filter=hoomd.filter.All(), kT=1.5)
simulation.operations.integrator.methods = [langevin]
filter#

Subset of particles to apply this method to.

Type:

hoomd.filter.filter_like

kT#

Temperature of the simulation \([\mathrm{energy}]\).

Examples:

langevin.kT = 1.0
langevin.kT = hoomd.variant.Ramp(A=2.0,
                                 B=1.0,
                                 t_start=0,
                                 t_ramp=1_000_000)
Type:

hoomd.variant.Variant

tally_reservoir_energy#

When True, track the energy exchange between the thermal reservoir and the particles.

Example:

langevin.tally_reservoir_energy = True
Type:

bool

gamma#

The drag coefficient for each particle type \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\).

Example:

langevin.gamma['A'] = 0.5
Type:

TypeParameter[ particle type, float ]

gamma_r#

The rotational drag coefficient tensor for each particle type \([\mathrm{time}^{-1}]\).

Example:

langevin.gamma_r['A'] = [1.0, 2.0, 3.0]
Type:

TypeParameter[particle type,[float, float , float]]

property reservoir_energy#

Energy absorbed by the reservoir \([\mathrm{energy}]\).

Set tally_reservoir_energy to True to track the reservoir energy.

Example:

langevin.tally_reservoir_energy = True
logger.add(obj=langevin, quantities=['reservoir_energy'])

Warning

When continuing a simulation, the energy of the reservoir will be reset to zero.

(Loggable: category=”scalar”)

class hoomd.md.methods.Method#

Bases: AutotunedObject

Base class integration method.

Provides common methods for all subclasses.

Note

Users should use the subclasses and not instantiate Method directly.

class hoomd.md.methods.OverdampedViscous(filter, default_gamma=1.0, default_gamma_r=(1.0, 1.0, 1.0))#

Bases: Method

Overdamped viscous dynamics.

Parameters:
  • filter (hoomd.filter.filter_like) – Subset of particles to apply this method to.

  • default_gamma (float) – Default drag coefficient for all particle types \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\).

  • default_gamma_r ([float, float, float]) – Default rotational drag coefficient tensor for all particles \([\mathrm{time}^{-1}]\).

OverdampedViscous integrates particles forward in time following Newtonian dynamics in the overdamped limit where there is no inertial term. (in the limit that the mass \(m\) and moment of inertia \(I\) go to 0):

\[ \begin{align}\begin{aligned}\frac{d\vec{r}}{dt} &= \vec{v}\\\vec{v(t)} &= \frac{\vec{F}_\mathrm{C}}{\gamma}\\\frac{d\mathbf{q}}{dt} &= \vec{\tau}\\\tau^i &= \frac{\tau_\mathrm{C}^i}{\gamma_r^i}\end{aligned}\end{align} \]

where \(\vec{F}_\mathrm{C} = \vec{F}_\mathrm{net}\) is the net force on the particle from all forces (hoomd.md.Integrator.forces) and constraints (hoomd.md.Integrator.constraints), \(\gamma\) is the translational drag coefficient (gamma), \(\vec{v}\) is the particle’s velocity, \(\tau_\mathrm{C}^i\) is the i-th component of the net torque from all forces and constraints, and \(\gamma_r^i\) is the i-th component of the rotational drag coefficient (gamma_r).

The attributes gamma and gamma_r set the translational and rotational damping coefficients, respectivley, by particle type.

Warning

This numerical method has errors in \(O(\delta t)\), which is much larger than the errors of most other integration methods which are in \(O(\delta t^2)\). As a consequence, expect to use much smaller values of \(\delta t\) with Brownian compared to e.g. Langevin or ConstantVolume.

Tip

OverdampedViscous can be used to simulate systems of athermal active matter.

Note

OverdampedViscous models systems in the limit that \(m\) and moment of inertia \(I\) go to 0. However, you must still set non-zero moments of inertia to enable the integration of rotational degrees of freedom.

Example:

overdamped_viscous = hoomd.md.methods.OverdampedViscous(
    filter=hoomd.filter.All())
simulation.operations.integrator.methods = [overdamped_viscous]
filter#

Subset of particles to apply this method to.

Type:

hoomd.filter.filter_like

gamma#

The drag coefficient for each particle type \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\).

Example:

overdamped_viscous.gamma['A'] = 0.5
Type:

TypeParameter[ particle type, float ]

gamma_r#

The rotational drag coefficient tensor for each particle type \([\mathrm{time}^{-1}]\).

Example:

overdamped_viscous.gamma_r['A'] = [1.0, 2.0, 3.0]
Type:

TypeParameter[particle type,[float, float , float]]

class hoomd.md.methods.Thermostatted#

Bases: Method

Base class for thermostatted integrators.

Provides a common interface for all methods using thermostats

Note

Users should use the subclasses and not instantiate Thermostatted directly.

Modules