md.methods¶
Overview
Base class integration method. |
|
Applies the Berendsen thermostat. |
|
Brownian dynamics. |
|
Newtonian dynamics with a cap on the maximum displacement per time step. |
|
Langevin dynamics. |
|
Constant pressure, constant enthalpy dynamics. |
|
Constant pressure, constant temperature dynamics. |
|
Constant volume, constant energy dynamics. |
|
Constant volume, constant temperature dynamics. |
|
Overdamped viscous dynamics. |
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.
Integration methods with constraints
For methods that constrain motion to a manifold see hoomd.md.methods.rattle
- class hoomd.md.methods.Berendsen(filter, kT, tau)¶
Bases:
Method
Applies the Berendsen thermostat.
- Parameters
filter (hoomd.filter.filter_like) – Subset of particles to apply this method to.
kT (hoomd.variant.variant_like) – Temperature of the simulation. \([energy]\)
tau (float) – Time constant of thermostat. \([time]\)
Berendsen
rescales the velocities of all particles on each time step. The rescaling is performed so that the difference in the current temperature from the set point decays exponentially: Berendsen et. al. 1984.\[\frac{dT_\mathrm{cur}}{dt} = \frac{T - T_\mathrm{cur}}{\tau}\]Attention
Berendsen
does not function with MPI parallel simulations.Attention
Berendsen
does not integrate rotational degrees of freedom.Examples:
berendsen = hoomd.md.methods.Berendsen( filter=hoomd.filter.All(), kT=0.2, tau=10.0) integrator = hoomd.md.Integrator( dt=0.001, methods=[berendsen], forces=[lj])
- filter¶
Subset of particles to apply this method to.
- kT¶
Temperature of the simulation. \([energy]\)
- class hoomd.md.methods.Brownian(filter, kT, alpha=None, 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}]\).
alpha (float) – When set, use \(\alpha d_i\) for the drag coefficient where \(d_i\) is particle diameter \([\mathrm{mass} \cdot \mathrm{length}^{-1} \cdot \mathrm{time}^{-1}]\). Defaults to
None
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. It integrates both the translational and rotational degrees of freedom.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, \(T\).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, \(T\).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.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 thathoomd.md.compute.ThermodynamicQuantities
will report appropriate temperatures and pressures when logged or used by other methods.Brownian dynamics neglects the acceleration term in the Langevin equation. This assumption is valid when overdamped: \(\frac{m}{\gamma} \ll \delta t\). Use
Langevin
if your system is not overdamped.You can set \(\gamma\) in two ways:
Specify \(\alpha\) which scales the particle diameter to \(\gamma = \alpha d_i\).
After the method object is created, specify the attribute
gamma
andgamma_r
for rotational damping or random torque to assign them directly, with independent values for each particle type in the system.
Examples:
brownian = hoomd.md.methods.Brownian(filter=hoomd.filter.All(), kT=0.2, alpha=1.0) integrator = hoomd.md.Integrator(dt=0.001, methods=[brownian], forces=[lj])
Examples of using
gamma
andgamma_r
:brownian = hoomd.md.methods.Brownian(filter=hoomd.filter.All(), kT=0.2) brownian.gamma.default = 2.0 brownian.gamma_r.default = [1.0, 2.0, 3.0]
- filter¶
Subset of particles to apply this method to.
- kT¶
Temperature of the simulation \([\mathrm{energy}]\).
- alpha¶
When set, use \(\alpha d_i\) for the drag coefficient where \(d_i\) is particle diameter \([\mathrm{mass} \cdot \mathrm{length}^{-1} \cdot \mathrm{time}^{-1}]\).
- Type
- class hoomd.md.methods.DisplacementCapped(filter, maximum_displacement: Union[Variant, float])¶
Bases:
NVE
Newtonian dynamics with a cap on the maximum displacement per time step.
The method employs a maximum displacement allowed each time step. This method can be helpful to relax a system with too much overlaps without “blowing up” the system.
Warning
This method does not conserve energy or momentum.
- Parameters
filter (hoomd.filter.filter_like) – Subset of particles on which to apply this method.
maximum_displacement (hoomd.variant.variant_like) – The maximum displacement allowed for a particular timestep \([\mathrm{length}]\).
DisplacementCapped
integrates integrates translational and rotational degrees of freedom using modified microcanoncial dynamics. SeeNVE
for the basis of the algorithm.Examples:
relaxer = hoomd.md.methods.DisplacementCapped( filter=hoomd.filter.All(), maximum_displacement=1e-3) integrator = hoomd.md.Integrator( dt=0.005, methods=[relaxer], forces=[lj])
- filter¶
Subset of particles on which to apply this method.
- maximum_displacement¶
The maximum displacement allowed for a particular timestep \([\mathrm{length}]\).
- class hoomd.md.methods.Langevin(filter, kT, alpha=None, 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}]\).
alpha (float) – When set, use \(\alpha d_i\) for the drag coefficient where \(d_i\) is particle diameter \([\mathrm{mass} \cdot \mathrm{length}^{-1} \cdot \mathrm{time}^{-1}]\). Defaults to None.
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.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, \(T\).
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, \(T\).Langevin
numerically integrates the translational degrees of freedom using Velocity-Verlet and the rotational degrees of freedom with a scheme based on Kamberaj 2005.Langevin dynamics includes the acceleration term in the Langevin equation. This assumption is valid when underdamped: \(\frac{m}{\gamma} \gg \delta t\). Use
Brownian
if your system is not underdamped.You can set \(\gamma\) in two ways:
Specify \(\alpha\) which scales the particle diameter to \(\gamma = \alpha d_i\).
After the method object is created, specify the attribute
gamma
andgamma_r
for rotational damping or random torque to assign them directly, with independent values for each particle type in the system.
Examples:
langevin = hoomd.md.methods.Langevin(filter=hoomd.filter.All(), kT=0.2, alpha=1.0) integrator = hoomd.md.Integrator(dt=0.001, methods=[langevin], forces=[lj])
Examples of using
gamma
andgamma_r
:langevin = hoomd.md.methods.Langevin(filter=hoomd.filter.All(), kT=0.2) langevin.gamma.default = 2.0 langevin.gamma_r.default = [1.0,2.0,3.0]
Warning
When restarting a simulation, the energy of the reservoir will be reset to zero.
- filter¶
Subset of particles to apply this method to.
- kT¶
Temperature of the simulation \([\mathrm{energy}]\).
- alpha¶
When set, use \(\alpha d_i\) for the drag coefficient where \(d_i\) is particle diameter \([\mathrm{mass} \cdot \mathrm{length}^{-1} \cdot \mathrm{time}^{-1}]\). Defaults to None.
- Type
- tally_reservoir_energy¶
When True, track the energy exchange between the thermal reservoir and the particles. \([\mathrm{energy}]\).
- Type
- gamma¶
The drag coefficient for each particle type. Used when
alpha
isNone
. \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\).- Type
TypeParameter[
particle type
,float
]
- gamma_r¶
The rotational drag coefficient tensor for each particle type \([\mathrm{time}^{-1}]\).
- property reservoir_energy¶
Energy absorbed by the reservoir \([\mathrm{energy}]\).
Set
tally_reservoir_energy
toTrue
to track the reservoir energy.(
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.NPH(filter, S, tauS, couple, box_dof=(True, True, True, False, False, False), rescale_all=False, gamma=0.0)¶
Bases:
Method
Constant pressure, constant enthalpy dynamics.
- Parameters
filter (hoomd.filter.filter_like) – Subset of particles on which to apply this method.
S (tuple[hoomd.variant.variant_like, ...] or hoomd.variant.variant_like) –
Stress components set point 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 (\([p, p, p, 0, 0, 0]\)), use
S = p
.tauS (float) – Coupling constant for the barostat \([\mathrm{time}]\).
couple (str) – Couplings of diagonal elements of the stress tensor, can be “none”, “xy”, “xz”,”yz”, or “all”, default to “all”.
box_dof (
tuple
[bool
]) – Box degrees of freedom with six boolean elements corresponding to x, y, z, xy, xz, yz, each. Default to [True,True,True,False,False,False]). If turned on to True, rescale corresponding lengths or tilt factors and components of particle coordinates and velocities.rescale_all (bool) – if True, rescale all particles, not only those in the group, Default to False.
gamma (float) – Dimensionless damping factor for the box degrees of freedom, Default to 0.
NPH
integrates translational and rotational degrees of freedom forward in time in the Isoenthalpic-isobaric ensemble. The barostat is introduced as additional degrees of freedom in the Hamiltonian that couple with the box parameters.The barostat tensor is \(\nu_{\mathrm{ij}}\). Access these quantities
barostat_dof
.See also
Except for the thermostat,
NPH
shares parameters withNPT
. SeeNPT
for descriptions of the coupling and other barostat parameters.Examples:
dt = 0.005 tauS = 1000 * dt nph = hoomd.md.methods.NPH(filter=hoomd.filter.All(), tauS=tauS, S=2.0) # orthorhombic symmetry nph = hoomd.md.methods.NPH(filter=hoomd.filter.All(), tauS=tauS, S=2.0, couple="none") # tetragonal symmetry nph = hoomd.md.methods.NPH(filter=hoomd.filter.All(), tauS=tauS, S=2.0, couple="xy") # triclinic symmetry nph = hoomd.md.methods.NPH(filter=hoomd.filter.All(), tauS=tauS, S=2.0, couple="none", rescale_all=True) integrator = hoomd.md.Integrator(dt=dt, methods=[nph], forces=[lj])
- filter¶
Subset of particles on which to apply this method.
- S¶
Stress components set point for the barostat totalling 6 components. In Voigt notation, \([S_{xx}, S_{yy}, S_{zz}, S_{yz}, S_{xz}, S_{xy}]\) \([\mathrm{pressure}]\). Stress can be reset after method object is created. For example, an isotopic pressure can be set by
nph.S = 4.
- Type
- couple¶
Couplings of diagonal elements of the stress tensor, can be “none”, “xy”, “xz”,”yz”, or “all”.
- Type
- box_dof¶
Box degrees of freedom with six boolean elements corresponding to x, y, z, xy, xz, yz, each.
- barostat_dof¶
Additional degrees of freedom for the barostat (\(\nu_{xx}\), \(\nu_{xy}\), \(\nu_{xz}\), \(\nu_{yy}\), \(\nu_{yz}\), \(\nu_{zz}\))
- property barostat_energy¶
Energy the barostat contributes to the Hamiltonian \([\mathrm{energy}]\).
(
Loggable
: category=”scalar”)
- thermalize_barostat_dof()¶
Set the barostat momentum to random values.
thermalize_barostat_dof
sets a random value for the barostat \(\nu_{\mathrm{ij}}\). Callthermalize_barostat_dof
to set a new random state for the barostat.Important
You must call
Simulation.run
beforethermalize_barostat_dof
. Callrun(steps=0)
to prepare a newly createdhoomd.Simulation
.See also
- class hoomd.md.methods.NPT(filter, kT, tau, S, tauS, couple, box_dof=[True, True, True, False, False, False], rescale_all=False, gamma=0.0)¶
Bases:
Method
Constant pressure, constant temperature dynamics.
- Parameters
filter (hoomd.filter.filter_like) – Subset of particles on which to apply this method.
kT (hoomd.variant.variant_like) – Temperature set point for the thermostat \([\mathrm{energy}]\).
tau (float) – Coupling constant for the thermostat \([\mathrm{time}]\).
S (tuple[hoomd.variant.variant_like, ...] or hoomd.variant.variant_like) –
Stress components set point 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 (\([p, p, p, 0, 0, 0]\)), use
S = p
.tauS (float) – Coupling constant for the barostat \([\mathrm{time}]\).
couple (str) – Couplings of diagonal elements of the stress tensor, can be “none”, “xy”, “xz”,”yz”, or “xyz”.
box_dof (
list
[bool
]) – Box degrees of freedom with six boolean elements corresponding to x, y, z, xy, xz, yz, each. Default to [True,True,True,False,False,False]). If turned on to True, rescale corresponding lengths or tilt factors and components of particle coordinates and velocities.rescale_all (bool) – if True, rescale all particles, not only those in the group, Default to False.
gamma (float) – Dimensionless damping factor for the box degrees of freedom, Default to 0.
NPT
integrates integrates translational and rotational degrees of freedom in the Isothermal-isobaric ensemble. The thermostat and barostat are introduced as additional degrees of freedom in the Hamiltonian that couple with the particle velocities and angular momenta and the box parameters.The translational thermostat has a momentum \(\xi\) and position \(\eta\). The rotational thermostat has momentum \(\xi_{\mathrm{rot}}\) and position \(\eta_\mathrm{rot}\). The barostat tensor is \(\nu_{\mathrm{ij}}\). Access these quantities using
translational_thermostat_dof
,rotational_thermostat_dof
, andbarostat_dof
.By default,
NPT
performs integration in a cubic box under hydrostatic pressure by simultaneously rescaling the lengths Lx, Ly and Lz of the simulation box. Set the integration mode to change this default.The integration mode is defined by a set of couplings and by specifying the box degrees of freedom that are put under barostat control. 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)
Degrees of freedom of the box specify 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 elementsbox_dof[:3]
controls whether the x, y, and z box lengths are rescaled and updated, respectively. The last three entriesbox_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
If 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.
For example:
Specifying all couplings and x, y, and z degrees of freedom amounts to cubic symmetry (default)
Specifying xy couplings and x, y, and z degrees of freedom amounts to tetragonal symmetry.
Specifying no couplings and all degrees of freedom amounts to a fully deformable triclinic unit cell
NPT
numerically integrates the equations of motion using the symplectic Martyna-Tobias-Klein equations of motion for NPT. For optimal stability, the update equations leave the phase-space measure invariant and are manifestly time-reversible.See also
Note
The coupling constant
tau
should be set within a reasonable range to avoid abrupt fluctuations in the kinetic temperature and to avoid long time to equilibration. The recommended value for most systems is \(\tau = 100 \delta t\).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 recommend value for most systems is \(\tau_S = 1000 \delta t\).Important
Ensure that your initial condition includes non-zero particle velocities and angular momenta (when appropriate). The coupling between the thermostat and the velocities and angular momenta occurs via multiplication, so
NPT
cannot convert a zero velocity into a non-zero one except through particle collisions.Examples:
npt = hoomd.md.methods.NPT(filter=hoomd.filter.All(), tau=1.0, kT=0.65, tauS = 1.2, S=2.0, couple="xyz") # orthorhombic symmetry npt = hoomd.md.methods.NPT(filter=hoomd.filter.All(), tau=1.0, kT=0.65, tauS = 1.2, S=2.0, couple="none") # tetragonal symmetry npt = hoomd.md.methods.NPT(filter=hoomd.filter.All(), tau=1.0, kT=0.65, tauS = 1.2, S=2.0, couple="xy") # triclinic symmetry npt = hoomd.md.methods.NPT(filter=hoomd.filter.All(), tau=1.0, kT=0.65, tauS = 1.2, S=2.0, couple="none", rescale_all=True) integrator = hoomd.md.Integrator(dt=0.005, methods=[npt], forces=[lj])
- filter¶
Subset of particles on which to apply this method.
- kT¶
Temperature set point for the thermostat \([\mathrm{energy}]\).
- 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}]\). Stress can be reset after the method object is created. For example, an isotropic pressure can be set by
npt.S = 4.
- Type
- couple¶
Couplings of diagonal elements of the stress tensor, can be “none”, “xy”, “xz”,”yz”, or “xyz”.
- Type
- box_dof¶
Box degrees of freedom with six boolean elements corresponding to x, y, z, xy, xz, yz, each.
- translational_thermostat_dof¶
Additional degrees of freedom for the translational thermostat (\(\xi\), \(\eta\))
- rotational_thermostat_dof¶
Additional degrees of freedom for the rotational thermostat (\(\xi_\mathrm{rot}\), \(\eta_\mathrm{rot}\))
- barostat_dof¶
Additional degrees of freedom for the barostat (\(\nu_{xx}\), \(\nu_{xy}\), \(\nu_{xz}\), \(\nu_{yy}\), \(\nu_{yz}\), \(\nu_{zz}\))
- property barostat_energy¶
Energy the barostat contributes to the Hamiltonian \([\mathrm{energy}]\).
(
Loggable
: category=”scalar”)
- thermalize_thermostat_and_barostat_dof()¶
Set the thermostat and barostat momenta to random values.
thermalize_thermostat_and_barostat_dof
sets a random value for the momentum \(\xi\) and the barostat \(\nu_{\mathrm{ij}}\). WhenIntegrator.integrate_rotational_dof
isTrue
, it also sets a random value for the rotational thermostat momentum \(\xi_{\mathrm{rot}}\). Callthermalize_thermostat_and_barostat_dof
to set a new random state for the thermostat and barostat.Important
You must call
Simulation.run
beforethermalize_thermostat_and_barostat_dof
. Callrun(steps=0)
to prepare a newly createdhoomd.Simulation
.See also
- class hoomd.md.methods.NVE(filter)¶
Bases:
Method
Constant volume, constant energy dynamics.
- Parameters
filter (hoomd.filter.filter_like) – Subset of particles on which to apply this method.
NVE
integrates integrates translational and rotational degrees of freedom in the microcanonical ensemble. The equations of motion are derived from the hamiltonian:\[H = U + K_\mathrm{translational} + K_\mathrm{rotational}\]NVE
numerically integrates the translational degrees of freedom using Velocity-Verlet and the rotational degrees of freedom with a scheme based on Kamberaj 2005.Examples:
nve = hoomd.md.methods.NVE(filter=hoomd.filter.All()) integrator = hoomd.md.Integrator(dt=0.005, methods=[nve], forces=[lj])
- filter¶
Subset of particles on which to apply this method.
- class hoomd.md.methods.NVT(filter, kT, tau)¶
Bases:
Method
Constant volume, constant temperature dynamics.
- Parameters
filter (hoomd.filter.filter_like) – Subset of particles on which to apply this method.
kT (hoomd.variant.variant_like) – Temperature set point for the Nosé-Hoover thermostat \([\mathrm{energy}]\).
tau (float) – Coupling constant for the Nosé-Hoover thermostat \([\mathrm{time}]\).
NVT
integrates integrates translational and rotational degrees of freedom in the canonical ensemble using the Nosé-Hoover thermostat. The thermostat is introduced as additional degrees of freedom in the Hamiltonian that couple with the velocities and angular momenta of the particles.The translational thermostat has a momentum \(\xi\) and position \(\eta\). The rotational thermostat has momentum \(\xi_{\mathrm{rot}}\) and position \(\eta_\mathrm{rot}\). Access these quantities using
translational_thermostat_dof
androtational_thermostat_dof
.NVT
numerically integrates the equations of motion using the symplectic Martyna-Tobias-Klein formalism described refs. G. J. Martyna, D. J. Tobias, M. L. Klein 1994 and J. Cao, G. J. Martyna 1996.Note
The coupling constant
tau
should be set within a reasonable range to avoid abrupt fluctuations in the kinetic temperature and to avoid long time to equilibration. The recommended value for most systems is \(\tau = 100 \delta t\).Important
Ensure that your initial condition includes non-zero particle velocities and angular momenta (when appropriate). The coupling between the thermostat and the velocities and angular momenta occurs via multiplication, so
NVT
cannot convert a zero velocity into a non-zero one except through particle collisions.Examples:
nvt=hoomd.md.methods.NVT(filter=hoomd.filter.All(), kT=1.0, tau=0.5) integrator = hoomd.md.Integrator(dt=0.005, methods=[nvt], forces=[lj])
- filter¶
Subset of particles on which to apply this method.
- kT¶
Temperature set point for the Nosé-Hoover thermostat \([\mathrm{energy}]\).
- translational_thermostat_dof¶
Additional degrees of freedom for the translational thermostat (\(\xi\), \(\eta\))
- rotational_thermostat_dof¶
Additional degrees of freedom for the rotational thermostat (\(\xi_\mathrm{rot}\), \(\eta_\mathrm{rot}\))
- thermalize_thermostat_dof()¶
Set the thermostat momenta to random values.
thermalize_thermostat_dof
sets a random value for the momentum \(\xi\). WhenIntegrator.integrate_rotational_dof
isTrue
, it also sets a random value for the rotational thermostat momentum \(\xi_{\mathrm{rot}}\). Callthermalize_thermostat_dof
to set a new random state for the thermostat.Important
You must call
Simulation.run
beforethermalize_thermostat_dof
. Callrun(steps=0)
to prepare a newly createdhoomd.Simulation
.See also
- class hoomd.md.methods.OverdampedViscous(filter, alpha=None, 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.
alpha (float) – When set, use \(\alpha d_i\) for the drag coefficient where \(d_i\) is particle diameter \([\mathrm{mass} \cdot \mathrm{length}^{-1} \cdot \mathrm{time}^{-1}]\). Defaults to
None
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, \(d\) is the dimensionality of the system, \(\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
).You can set \(\gamma\) in two ways:
Specify \(\alpha\) which scales the particle diameter to \(\gamma = \alpha d_i\).
After the method object is created, specify the attribute
gamma
andgamma_r
for rotational damping or random torque to assign them directly, with independent values for each particle type in the system.
Tip
OverdampedViscous
can be used to simulate systems of athermal active matter, such as athermal Active Brownian Particles.Note
Even though
OverdampedViscous
models systems in the limit that \(m\) and moment of inertia \(I\) go to 0, you must still set non-zero moments of inertia to enable the integration of rotational degrees of freedom.Examples:
odv = hoomd.md.methods.OverdampedViscous(filter=hoomd.filter.All()) odv.gamma.default = 2.0 odv.gamma_r.default = [1.0, 2.0, 3.0]
- filter¶
Subset of particles to apply this method to.
- alpha¶
When set, use \(\alpha d_i\) for the drag coefficient where \(d_i\) is particle diameter \([\mathrm{mass} \cdot \mathrm{length}^{-1} \cdot \mathrm{time}^{-1}]\).
- Type
Modules