md.compute¶
Overview
Compute harmonic averaged thermodynamic properties of particles. |
|
Compute thermodynamic properties of a subset of the system. |
Details
Compute properties of molecular dynamics simulations.
The MD compute classes compute instantaneous properties of the simulation state
and provide results as loggable quantities for use with hoomd.logging.Logger
or by direct access via the Python API.
- class hoomd.md.compute.HarmonicAveragedThermodynamicQuantities(filter, kT, harmonic_pressure=0)¶
Bases:
Compute
Compute harmonic averaged thermodynamic properties of particles.
- Parameters:
filter (hoomd.filter.filter_like) – Particle filter to compute thermodynamic properties for.
kT (float) – Temperature of the system \([\mathrm{energy}]\).
harmonic_pressure (float) – Harmonic contribution to the pressure \([\mathrm{pressure}]\). If omitted, the HMA pressure can still be computed, but will be similar in precision to the conventional pressure.
HarmonicAveragedThermodynamicQuantities
acts on a given subset of particles and calculates harmonically mapped average (HMA) properties of those particles when requested. HMA computes properties more precisely (with less variance) for atomic crystals in NVT simulations. The presence of diffusion (vacancy hopping, etc.) will prevent HMA from providing improvement. HMA tracks displacements from the lattice positions, which are saved either during first call toSimulation.run
or when the compute is first added to the simulation, whichever occurs last.Note
HarmonicAveragedThermodynamicQuantities
is an implementation of the methods section of Sabry G. Moustafa, Andrew J. Schultz, and David A. Kofke. (2015). “Very fast averaging of thermal properties of crystals by molecular simulation”. Phys. Rev. E 92, 043303 doi:10.1103/PhysRevE.92.043303Examples:
hma = hoomd.compute.HarmonicAveragedThermodynamicQuantities( filter=hoomd.filter.Type('A'), kT=1.0)
- filter¶
Subset of particles compute thermodynamic properties for.
- Type:
- class hoomd.md.compute.ThermodynamicQuantities(filter)¶
Bases:
Compute
Compute thermodynamic properties of a subset of the system.
- Parameters:
filter (
hoomd.filter
) – Particle filter to compute thermodynamic properties for.
ThermodynamicQuantities
acts on a subset of particles in the system and calculates thermodynamic properties of those particles. Add aThermodynamicQuantities
instance to a logger to save these quantities to a file, seehoomd.logging.Logger
for more details.Note
For compatibility with
hoomd.md.constrain.Rigid
,ThermodynamicQuantities
performs all sums \(\sum_{i \in \mathrm{filter}}\) over free particles and rigid body centers - ignoring constituent particles to avoid double counting.Examples:
f = filter.Type('A') compute.ThermodynamicQuantities(filter=f)
- property degrees_of_freedom¶
Number of degrees of freedom in the subset \(N_{\mathrm{dof}}\).
\[N_\mathrm{dof} = N_\mathrm{dof, translational} + N_\mathrm{dof, rotational}\]See also
hoomd.State.update_group_dof
describes when \(N_{\mathrm{dof}}\) is updated.(
Loggable
: category=”scalar”)
- property kinetic_energy¶
Total kinetic energy \(K\) of the subset \([\mathrm{energy}]\).
\[K = K_\mathrm{rotational} + K_\mathrm{translational}\](
Loggable
: category=”scalar”)
- property kinetic_temperature¶
Instantaneous thermal energy \(kT_k\) of the subset \([\mathrm{energy}]\).
\[kT_k = 2 \cdot \frac{K}{N_{\mathrm{dof}}}\](
Loggable
: category=”scalar”)
- property potential_energy¶
Potential energy \(U\) that the subset contributes to the system \([\mathrm{energy}]\).
\[U = U_{\mathrm{net},\mathrm{additional}} + \sum_{i \in \mathrm{filter}} U_{\mathrm{net},i},\]where the net energy terms are computed by
hoomd.md.Integrator
over all of the forces inhoomd.md.Integrator.forces
.(
Loggable
: category=”scalar”)
- property pressure¶
Instantaneous pressure \(P\) of the subset \([\mathrm{pressure}]\).
\[P = \frac{ 2 \cdot K_{\mathrm{translational}} + W_\mathrm{isotropic} }{D \cdot V},\]where \(D\) is the dimensionality of the system, \(V\) is the volume of the simulation box (or area in 2D), and \(W_\mathrm{isotropic}\) is the isotropic virial:
\[\begin{split}W_\mathrm{isotropic} = & \left( W_{\mathrm{net},\mathrm{additional}}^{xx} + W_{\mathrm{net},\mathrm{additional}}^{yy} + W_{\mathrm{net},\mathrm{additional}}^{zz} \right) \\ + & \sum_{i \in \mathrm{filter}} \left( W_\mathrm{{net},i}^{xx} + W_\mathrm{{net},i}^{yy} + W_\mathrm{{net},i}^{zz} \right)\end{split}\]where the net virial terms are computed by
hoomd.md.Integrator
over all of the forces inhoomd.md.Integrator.forces
and \(W^{zz}=0\) in 2D simulations.(
Loggable
: category=”scalar”)
- property pressure_tensor¶
Instantaneous pressure tensor of the subset \([\mathrm{pressure}]\).
The six components of the pressure tensor are given in the order: (\(P^{xx}\), \(P^{xy}\), \(P^{xz}\), \(P^{yy}\), \(P^{yz}\), \(P^{zz}\)):
\[P^{kl} = \frac{1}{V} \left( W_{\mathrm{net},\mathrm{additional}}^{kl} + \sum_{i \in \mathrm{filter}} m_i \cdot v_i^k \cdot v_i^l + W_{\mathrm{net},i}^{kl} \right),\]where the net virial terms are computed by
hoomd.md.Integrator
over all of the forces inhoomd.md.Integrator.forces
, \(v_i^k\) is the k-th component of the velocity of particle \(i\) and \(V\) is the total simulation box volume (or area in 2D).(
Loggable
: category=”sequence”)
- property rotational_degrees_of_freedom¶
Number of rotational degrees of freedom in the subset \(N_\mathrm{dof, rotational}\).
Integration methods (
hoomd.md.methods
) determine the number of degrees of freedom they give to each particle. Each integration method operates on a subset of the system \(\mathrm{filter}_m\) that may be distinct from the subset from the subset given toThermodynamicQuantities
.When
hoomd.md.Integrator.integrate_rotational_dof
isFalse
, \(N_\mathrm{dof, rotational} = 0\). When it isTrue
, the given degrees of freedom depend on the dimensionality of the system.In 2D:
\[N_\mathrm{dof, rotational} = \sum_{m \in \mathrm{methods}} \; \sum_{i \in \mathrm{filter} \cup \mathrm{filter}_m} \left[ I_{z,i} > 0 \right]\]where \(I\) is the particles’s moment of inertia.
In 3D:
\[N_\mathrm{dof, rotational} = \sum_{m \in \mathrm{methods}} \; \sum_{i \in \mathrm{filter} \cup \mathrm{filter}_m} \left[ I_{x,i} > 0 \right] + \left[ I_{y,i} > 0 \right] + \left[ I_{z,i} > 0 \right]\]See also
hoomd.State.update_group_dof
describes when \(N_{\mathrm{dof, rotational}}\) is updated.(
Loggable
: category=”scalar”)
- property rotational_kinetic_energy¶
Rotational kinetic energy \(K_\mathrm{rotational}\) of the subset \([\mathrm{energy}]\).
\[ \begin{align}\begin{aligned}\begin{split}K_\mathrm{rotational,d} = \frac{1}{2} \sum_{i \in \mathrm{filter}} \begin{cases} \frac{L_{d,i}^2}{I_{d,i}} & I^d_i > 0 \\ 0 & I^d_i = 0 \end{cases}\end{split}\\K_\mathrm{rotational} = K_\mathrm{rotational,x} + K_\mathrm{rotational,y} + K_\mathrm{rotational,z}\end{aligned}\end{align} \]\(I\) is the moment of inertia and \(L\) is the angular momentum in the (diagonal) reference frame of the particle.
(
Loggable
: category=”scalar”)
- property translational_degrees_of_freedom¶
Number of translational degrees of freedom in the subset \(N_{\mathrm{dof, translational}}\).
When using a single integration method on all particles that is momentum conserving, the center of mass motion is conserved and the number of translational degrees of freedom is:
\[N_\mathrm{dof, translational} = DN - D\frac{N}{N_\mathrm{particles}} - N_\mathrm{constraints}(\mathrm{filter})\]where \(D\) is the dimensionality of the system and \(N_\mathrm{constraints}(\mathrm{filter})\) is the number of degrees of freedom removed by constraints (
hoomd.md.constrain
) in the subset. The fraction \(\frac{N}{N_\mathrm{particles}}\) distributes the momentum conservation constraint evenly whenThermodynamicQuantities
is applied to multiple subsets.Note
When using rigid bodies (
hoomd.md.constrain.Rigid
), \(N\) is the number of rigid body centers plus free particles selected by the filter and \(N_\mathrm{particles}\) is total number of rigid body centers plus free particles in the whole system.When
hoomd.md.Integrator.rigid
is not set, \(N\) is the total number of particles selected by the filter and \(N_\mathrm{particles}\) is the total number of particles in the system, regardless of theirbody
value.When using multiple integration methods, a single integration method on fewer than all particles, or a single integration method that is not momentum conserving,
hoomd.md.Integrator
assumes that linear momentum is not conserved and counts the center of mass motion in the degrees of freedom:\[N_{\mathrm{dof, translational}} = DN - N_\mathrm{constraints}(\mathrm{filter})\](
Loggable
: category=”scalar”)