md.pair¶
Overview
Buckingham pair force. |
|
DLVO colloidal interaction. |
|
Dissipative Particle Dynamics. |
|
Dissipative Particle Dynamics with the LJ conservative force. |
|
DPD Conservative pair force. |
|
Ewald pair force. |
|
Expanded Gaussian pair force. |
|
Expanded Lennard-Jones pair force. |
|
Expanded Mie pair force. |
|
Force-shifted Lennard-Jones pair force. |
|
Fourier pair force. |
|
Gaussian pair force. |
|
Lennard-Jones pair force. |
|
Lennard-Jones 12-8 pair force. |
|
Lennard-Jones 8-4 pair force. |
|
Lennard-Jones-Gauss pair potential. |
|
Mie pair force. |
|
Morse pair force. |
|
Moliere pair force. |
|
Oscillating pair force. |
|
Base class pair force. |
|
Onsager reaction field pair force. |
|
Tabulated pair force. |
|
Pair potential model for globular proteins. |
|
Yukawa pair force. |
|
ZBL pair force. |
Details
Pair Potentials for molecular dynamics.
Pair force classes apply a force and virial on every particle in the simulation state commensurate with the potential energy:
where \(\vec{r}_{ij} = \mathrm{minimum\_image}(\vec{r}_j - \vec{r}_i)\).
Pair
applies a short range cutoff using a hoomd.md.nlist.NeighborList
for
performance and assumes that both \(U(r)\) and its derivatives are 0 when
\(r_{ij} \ge r_\mathrm{cut}\). Pair
also ignores particle pairs that are
excluded in the neighbor list.
Specifically, the force \(\vec{F}\) on each pair of particles \(i,j\) is:
where the cutoff radius \(r_{\mathrm{cut}}\) is given by Pair.r_cut
.
Tip
Set Pair.r_cut
to 0 to skip computations for non-interacting pairs.
Pair
splits half the energy from each pair interaction onto particles
\(i\) and \(j\):
and similarly for virials.
Shifting/smoothing mode
The function \(U_\mathrm{pair}(r)\) depends on the chosen form of the pair
potential \(U(r)\) (by the Pair
subclass) and the mode (Pair.mode
):
where \(S(r)\) is the XPLOR smoothing function:
where \(r_{\mathrm{on}}\) is given by Pair.r_on
.
The XPLOR smoothing function \(S(r)\) ensures that both the potential energy and the force going smoothly to 0 at \(r = r_{\mathrm{cut}}\), reducing the rate of energy drift in long simulations. \(r_{\mathrm{on}}\) controls the point at which the smoothing starts. Set it to modify only the tail of the potential. The WCA potential and it’s first derivative already go smoothly to 0 at the cutoff, so there is no need to apply the smoothing function. In such mixed systems, set \(r_{\mathrm{on}}\) to a value greater than \(r_{\mathrm{cut}}\) for those pairs that interact via WCA in order to enable shifting of the WCA potential to 0 at the cutoff.
Tail correction
Some pair potentials can optionally apply isotropic integrated long range tail
corrections when the tail_correction
parameter is True
. These
corrections are only valid when the shifting/smoothing mode is set to
"none"
. Following Sun 1998, the
pressure and energy corrections \(\Delta P\) and \(\Delta E\) are given
by:
and
where \(n\) is the number of unique particle types in the system, \(\rho_i\) is the number density of particles of type \(i\) in the system, \(U_{ij}(r)\) is the pair potential between particles of type \(i\) and \(j\), and \(N_i\) is the number of particles of type \(i\) in the system. These expressions assume that the radial pair distribution functions \(g_{ij}(r)\) are unity at the cutoff and beyond.
The pressure shift \(\Delta P\) appears in the additional virial term
\(W_\mathrm{additional}\) (Force.additional_virial
) and the energy shift
appears in the additional energy \(U_\mathrm{additional}\)
(Force.additional_energy
).
Warning
The value of the tail corrections depends on the number of each type of particle in the system, and these are precomputed when the pair potential object is initialized. If the number of any of the types of particles changes, the tail corrections will yield invalid results.
Anisotropic potentials
For anisotropic potentials see hoomd.md.pair.aniso
- class hoomd.md.pair.Buckingham(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Buckingham pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
Buckingham
computes the Buckingham pair force on every particle in the simulation state with:\[U(r) = A \exp\left(-\frac{r}{\rho}\right) - \frac{C}{r^6}\]Example:
nl = nlist.Cell() buck = pair.Buckingham(nl, default_r_cut=3.0) buck.params[('A', 'A')] = {'A': 2.0, 'rho'=0.5, 'C': 1.0} buck.params[('A', 'B')] = dict(A=1.0, rho=1.0, C=1.0) buck.params[('B', 'B')] = dict(A=2.0, rho=2.0, C=2.0)
- class hoomd.md.pair.DLVO(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
DLVO colloidal interaction.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
name (str) – Name of the force instance.
mode (str) – Energy shifting mode.
DLVO
computes the DLVO dispersion and electrostatic interaction pair force on every particle in the simulation state with:\[\begin{split}V_{\mathrm{DLVO}}(r) = &- \frac{A}{6} \left[ \frac{2a_1a_2}{r^2 - (a_1+a_2)^2} + \frac{2a_1a_2}{r^2 - (a_1-a_2)^2} \\ + \log \left( \frac{r^2 - (a_1+a_2)^2}{r^2 - (a_1-a_2)^2} \right) \right] \\ & + \frac{a_1 a_2}{a_1+a_2} Z e^{-\kappa(r - (a_1+a_2))}\end{split}\]where \(a_1\) is the radius of first particle in the pair, \(a_2\) is the radius of second particle in the pair, \(A\) is the Hamaker constant, \(Z\) is proportional to the surface electric potential, and \(\kappa\) is the screening parameter.
The first term corresponds to the attractive van der Waals interaction with and the second term to the repulsive double-layer interaction between two spherical surfaces. See “Intermolecular and Surface Forces” Israelachvili 2011, pp. 317.
Example:
nl = hoomd.md.nlist.Cell() dlvo = hoomd.md.pair.DLVO(nlist=nl) dlvo.params[('A', 'A')] = dict(A=1.0, kappa=1.0, Z=2, a1=1, a2=1) dlvo.params[('A', 'B')] = dict(A=2.0, kappa=0.5, Z=3, a1=1, a2=3) dlvo.params[('B', 'B')] = dict(A=2.0, kappa=0.5, Z=3, a1=3, a2=3)
- params¶
The potential parameters. The dictionary has the following keys:
A
(float
, required) - Hamaker constant \(A\) \([\mathrm{energy}]\)a1
(float
, required) - Radius of first particle \(a_1\) \([\mathrm{length}]\)a2
(float
, required) - Radius of second particle \(a_2\) \([\mathrm{length}]\)kappa
(float
, required) - screening parameter \(\kappa\) \([\mathrm{length}^{-1}]\)Z
surface electric potential (float
, required) - \(Z\) \([\mathrm{energy} \cdot \mathrm{length}^{-1}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.DPD(nlist, kT, default_r_cut=None)¶
Bases:
Pair
Dissipative Particle Dynamics.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list
kT (
hoomd.variant
orfloat
) – Temperature of thermostat \([\mathrm{energy}]\).default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
DPD
computes the DPD pair force on every particle in the simulation state. DPD includes a an interaction potential, pairwise drag force, and pairwise random force. See Groot and Warren 1997:\[F = F_{\mathrm{C}}(r) + F_{\mathrm{R,ij}}(r_{ij}) + F_{\mathrm{D,ij}}(v_{ij})\]where
\[\begin{split}F_{\mathrm{C}}(r) &= A \cdot w(r_{ij}), \\ F_{\mathrm{R, ij}}(r_{ij}) &= - \theta_{ij}\sqrt{3} \sqrt{\frac{2k_b\gamma T}{\Delta t}}\cdot w(r_{ij}), \\ F_{\mathrm{D, ij}}(r_{ij}) &= - \gamma w^2(r_{ij})\left( \hat r_{ij} \circ v_{ij} \right), \\ w(r_{ij}) &= \begin{cases} \left( 1 - r/r_{\mathrm{cut}} \right) & r < r_{\mathrm{cut}} \\ 0 & r \ge r_{\mathrm{cut}} \\ \end{cases},\end{split}\]\(\hat r_{ij}\) is a normalized vector from particle i to particle j, \(v_{ij} = v_i - v_j\), and \(\theta_{ij}\) is a uniformly distributed random number in the range \([-1, 1]\).
C. L. Phillips et. al. 2011 describes the DPD implementation details. Cite it if you utilize the DPD functionality in your work.
DPD
does not implement any energy shift / smoothing modes due to the function of the force.To use the DPD thermostat, apply the
hoomd.md.methods.ConstantVolume
orhoomd.md.methods.ConstantPressure
integration method without thermostats along withDPD
forces. Use of the DPD thermostat pair force with other integrators will result in nonphysical behavior. To useDPD
with a different conservative potential than \(F_C\), set A to zero and define the conservative pair force separately.Example:
nl = nlist.Cell() dpd = pair.DPD(nlist=nl, kT=1.0, default_r_cut=1.0) dpd.params[('A', 'A')] = dict(A=25.0, gamma=4.5) dpd.params[('A', 'B')] = dict(A=40.0, gamma=4.5) dpd.params[('B', 'B')] = dict(A=25.0, gamma=4.5) dpd.params[(['A', 'B'], ['C', 'D'])] = dict(A=40.0, gamma=4.5)
- class hoomd.md.pair.DPDConservative(nlist, default_r_cut=None)¶
Bases:
Pair
DPD Conservative pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
DPDConservative
computes the conservative part of theDPD
pair force on every particle in the simulation state with:\[U(r) = A \cdot \left( r_{\mathrm{cut}} - r \right) - \frac{1}{2} \cdot \frac{A}{r_{\mathrm{cut}}} \cdot \left(r_{\mathrm{cut}}^2 - r^2 \right).\]DPDConservative
does not implement any energy shift / smoothing modes due to the function of the force.Example:
nl = nlist.Cell() dpdc = pair.DPDConservative(nlist=nl, default_r_cut=3.0) dpdc.params[('A', 'A')] = dict(A=1.0) dpdc.params[('A', 'B')] = dict(A=2.0, r_cut = 1.0) dpdc.params[(['A', 'B'], ['C', 'D'])] = dict(A=3.0)
- params¶
The potential parameters. The dictionary has the following keys:
A
(float
, required) - \(A\) \([\mathrm{force}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.DPDLJ(nlist, kT, default_r_cut=None, mode='none')¶
Bases:
Pair
Dissipative Particle Dynamics with the LJ conservative force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
kT (
hoomd.variant
orfloat
) – Temperature of thermostat \([\mathrm{energy}]\).default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
mode (str) – Energy shifting mode.
DPDLJ
computes theDPD
thermostat combined with theLJ
pair force on every particle in the simulation state with:\[\begin{split}F &= F_{\mathrm{C}}(r) + F_{\mathrm{R,ij}}(r_{ij}) + F_{\mathrm{D,ij}}(v_{ij}), \\ F_{\mathrm{C}}(r) &= \partial U / \partial r, \\ F_{\mathrm{R, ij}}(r_{ij}) &= - \theta_{ij}\sqrt{3} \sqrt{\frac{2k_b\gamma T}{\Delta t}}\cdot w(r_{ij}), \\ F_{\mathrm{D, ij}}(r_{ij}) &= - \gamma w^2(r_{ij}) \left( \hat r_{ij} \circ v_{ij} \right), \\ U(r) &= 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right], \\ w(r_{ij}) &= \begin{cases} \left( 1 - r/r_{\mathrm{cut}} \right) & r < r_{\mathrm{cut}} \\ 0 & r \ge r_{\mathrm{cut}} \\ \end{cases},\end{split}\]\(\hat r_{ij}\) is a normalized vector from particle i to particle j, \(v_{ij} = v_i - v_j\), and \(\theta_{ij}\) is a uniformly distributed random number in the range [-1, 1].
C. L. Phillips et. al. 2011 describes the DPD implementation details. Cite it if you utilize the DPD functionality in your work.
To use the DPD thermostat, apply the
hoomd.md.methods.ConstantVolume
orhoomd.md.methods.ConstantPressure
integration method without thermostat along withDPD
forces. Use of the DPD thermostat pair force with other integrators will result in nonphysical behavior.Example:
nl = nlist.Cell() dpdlj = pair.DPDLJ(nlist=nl, kT=1.0, default_r_cut=2.5) dpdlj.params[('A', 'A')] = dict(epsilon=1.0, sigma=1.0, gamma=4.5) dpdlj.params[(['A', 'B'], ['C', 'D'])] = dict( epsilon=3.0, sigma=1.0, gamma=1.2) dpdlj.r_cut[('B', 'B')] = 2.0**(1.0/6.0)
- params¶
The DPDLJ potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - \(\varepsilon\) \([\mathrm{energy}]\)sigma
(float
, required) - \(\sigma\) \([\mathrm{length}]\)gamma
(float
, required) - \(\gamma\) \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.Ewald(nlist, default_r_cut=None)¶
Bases:
Pair
Ewald pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
Ewald
computes the Ewald pair force on every particle in the simulation state:\[U(r) = q_i q_j \left[\mathrm{erfc}\left(\kappa r + \frac{\alpha}{2\kappa}\right) \exp(\alpha r) + \mathrm{erfc}\left(\kappa r - \frac{\alpha}{2 \kappa}\right) \exp(-\alpha r)\right]\]Call
md.long_range.pppm.make_pppm_coulomb_forces
to create an instance ofEwald
andmd.long_range.pppm.Coulomb
that together implement the PPPM method for electrostatics.Example:
nl = nlist.Cell() ewald = pair.Ewald(default_r_cut=3.0, nlist=nl) ewald.params[('A', 'A')] = dict(kappa=1.0, alpha=1.5) ewald.r_cut[('A', 'B')] = 3.0
- params¶
The Ewald potential parameters. The dictionary has the following keys:
kappa
(float
, required) - Splitting parameter \(\kappa\) \([\mathrm{length}^{-1}]\)alpha
(float
, required) - Debye screening length \(\alpha\) \([\mathrm{length}^{-1}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.ExpandedGaussian(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Expanded Gaussian pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
ExpandedGaussian
computes the radially-shifted Gaussian pair force should on every particle in the simulation state:\[U(r) = \varepsilon \exp \left( -\frac{1}{2} \left( \frac{r-\Delta}{\sigma} \right)^2 \right)\]Example:
nl = nlist.Cell() expanded_gauss = pair.ExpandedGaussian(default_r_cut=3.0, nlist=nl) expanded_gauss.params[('A', 'A')] = dict(epsilon=1.0, sigma=1.0, delta=0.5) expanded_gauss.r_cut[('A', 'B')] = 3.0
- params¶
The expanded Gaussian potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - energy parameter \(\varepsilon\) \([\mathrm{energy}]\)sigma
(float
, required) - particle size \(\sigma > 0\) \([\mathrm{length}]\)delta
(float
, required) - shift distance \(\delta\) \([\mathrm{length}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.ExpandedLJ(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Expanded Lennard-Jones pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting mode.
ExpandedLJ
computes the radially-shifted Lennard-Jones pair force on every particle in the simulation state:\[U(r) = 4 \varepsilon \left[ \left( \frac{\sigma}{r - \Delta} \right)^{12} - \left( \frac{\sigma}{r - \Delta} \right)^{6} \right]\]Example:
nl = nlist.Cell() expanded_lj = pair.ExpandedLJ(default_r_cut=3.0, nlist=nl) expanded_lj.params[('A', 'A')] = dict(epsilon=1.0, sigma=1.0, delta=1.0) expanded_lj.params[('A', 'B')] = dict( epsilon=2.0, sigma=1.0, delta=0.75) expanded_lj.params[('B', 'B')] = dict(epsilon=1.0, sigma=1.0, delta=0.5)
- params¶
The potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - energy parameter \(\varepsilon\) \([\mathrm{energy}]\)sigma
(float
, required) - particle size \(\sigma\) \([\mathrm{length}]\)delta
(float
, required) - radial shift \(\Delta\) \([\mathrm{length}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.ExpandedMie(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Expanded Mie pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
ExpandedMie
computes the radially shifted Mie pair force on every particle in the simulation state:\[U(r) = \left( \frac{n}{n-m} \right) {\left( \frac{n}{m} \right)}^{\frac{m}{n-m}} \varepsilon \left[ \left( \frac{\sigma}{r-\Delta} \right)^{n} - \left( \frac {\sigma}{r-\Delta} \right)^{m} \right]\]Example:
nl = nlist.Cell() expanded_mie = pair.ExpandedMie(nlist=nl, default_r_cut=3.0) mie.params[('A', 'B')] = { "epsilon": 1.0, "sigma": 1.0, "n": 12, "m": 6, "delta": 0.5} expanded_mie.r_cut[('A', 'B')] = 2**(1.0 / 6.0) expanded_mie.params[(['A', 'B'], ['C', 'D'])] = { "epsilon": 1.5, "sigma": 2.0, "n": 12, "m": 6, "delta": 0.5}
- params¶
The Expanded Mie potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - \(\epsilon\) \([\mathrm{energy}]\).sigma
(float
, required) - \(\sigma\) \([\mathrm{length}]\).n
(float
, required) - \(n\) \([\mathrm{dimensionless}]\).m
(float
, required) - \(m\) \([\mathrm{dimensionless}]\).delta
(float
, required) - \(\Delta\) \([\mathrm{length}]\).
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.ForceShiftedLJ(nlist, default_r_cut=None)¶
Bases:
Pair
Force-shifted Lennard-Jones pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
ForceShiftedLJ
computes the modified Lennard-Jones pair force on every particle in the simulation state.\[U(r) = 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right] + \Delta V(r)\]\[\Delta V(r) = -(r - r_{\mathrm{cut}}) \frac{\partial V_{\mathrm{LJ}}}{\partial r}(r_{\mathrm{cut}})\]The force differs from the one calculated by
LJ
by the subtraction of the value of the force at \(r_{\mathrm{cut}}\), such that the force smoothly goes to zero at the cut-off. The potential is modified by a linear function. See Toxvaerd et. al. 2011 for a discussion of this potential.Example:
nl = nlist.Cell() fslj = pair.ForceShiftedLJ(nlist=nl, default_r_cut=1.5) fslj.params[('A', 'A')] = dict(epsilon=1.0, sigma=1.0)
- class hoomd.md.pair.Fourier(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Fourier pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
Fourier
computes the Fourier pair force on every particle in the simulation state:\[U(r) = \frac{1}{r^{12}} + \frac{1}{r^2}\sum_{n=1}^4 \left[ a_n cos \left( \frac{n \pi r}{r_{cut}} \right) + b_n sin \left( \frac{n \pi r}{r_{cut}} \right) \right]\]where
\[\begin{split}a_1 &= \sum_{n=2}^4 (-1)^n a_n \\ b_1 &= \sum_{n=2}^4 n (-1)^n b_n \\\end{split}\]enforce \(U(r_\mathrm{cut}) = 0\).
Example:
nl = nlist.Cell() fourier = pair.Fourier(default_r_cut=3.0, nlist=nl) fourier.params[('A', 'A')] = dict(a=[a2,a3,a4], b=[b2,b3,b4])
- params¶
The Fourier potential parameters. The dictionary has the following keys:
a
(float
, required) - array of 3 values corresponding to a2, a3 and a4 in the Fourier series \([\mathrm{dimensionless}]\)b
(float
, required) - array of 3 values corresponding to b2, b3 and b4 in the Fourier series \([\mathrm{dimensionless}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.Gaussian(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Gaussian pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
Gaussian
computes the Gaussian pair force should on every particle in the simulation state:\[U(r) = \varepsilon \exp \left( -\frac{1}{2} \left( \frac{r}{\sigma} \right)^2 \right)\]Example:
nl = nlist.Cell() gauss = pair.Gaussian(default_r_cut=3.0, nlist=nl) gauss.params[('A', 'A')] = dict(epsilon=1.0, sigma=1.0) gauss.r_cut[('A', 'B')] = 3.0
- params¶
The Gauss potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - energy parameter \(\varepsilon\) \([\mathrm{energy}]\)sigma
(float
, required) - particle size \(\sigma > 0\) \([\mathrm{length}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.LJ(nlist, default_r_cut=None, default_r_on=0.0, mode='none', tail_correction=False)¶
Bases:
Pair
Lennard-Jones pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
tail_correction (bool) – Whether to apply the isotropic integrated long range tail correction.
LJ
computes the Lennard-Jones pair force on every particle in the simulation state.\[U(r) = 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right]\]Example:
nl = nlist.Cell() lj = pair.LJ(nl, default_r_cut=3.0) lj.params[('A', 'A')] = {'sigma': 1.0, 'epsilon': 1.0} lj.r_cut[('A', 'B')] = 3.0
- class hoomd.md.pair.LJ0804(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Lennard-Jones 8-4 pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
LJ0804
computes the Lennard-Jones 8-4 pair force on every particle in the simulation state:\[U(r) = 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{8} - \left( \frac{\sigma}{r} \right)^{4} \right]\]Example:
nl = nlist.Cell() lj0804 = pair.LJ0804(nl, default_r_cut=3.0) lj0804.params[('A', 'A')] = {'sigma': 1.0, 'epsilon': 1.0} lj0804.params[('A', 'B')] = dict(epsilon=2.0, sigma=1.0) lj0804.r_cut[('A', 'B')] = 3.0
- class hoomd.md.pair.LJ1208(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Lennard-Jones 12-8 pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
LJ1208
computes the Lennard-Jones 12-8 pair force on every particle in the simulation state.Example:
nl = nlist.Cell() lj1208 = pair.LJ1208(nl, default_r_cut=3.0) lj1208.params[('A', 'A')] = {'sigma': 1.0, 'epsilon': 1.0} lj1208.params[('A', 'B')] = dict(epsilon=2.0, sigma=1.0)
\[U(r) = 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{8} \right]\]
- class hoomd.md.pair.LJGauss(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Lennard-Jones-Gauss pair potential.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
LJGauss
computes the Lennard-Jones Gauss force on all particles in the simulation state:\[U(r) = 1\ [\mathrm{energy}] \cdot \left[ \left ( \frac{1\ [\mathrm{length}]}{r} \right)^{12} - 2\ \left(\frac{1 [\mathrm{length}]}{r} \right)^{6} \right] - \epsilon \exp \left[- \frac{\left(r - r_{0}\right)^{2}}{2 \sigma^{2}} \right]\]- params¶
The potential parameters. The dictionary has the following keys:
Example:
nl = hoomd.md.nlist.Cell() ljg = pair.LJGauss(nl) ljg.params[('A', 'A')] = dict(epsilon=1.0, sigma=0.02, r0=1.6) ljg.params[('A', 'B')] = {'epsilon' : 2.0, 'sigma' : 0.02, 'r0' : 1.6} ljg.params[('A', 'B')] = {'epsilon' : 2.0, 'sigma' : 0.02, 'r0' : 1.6}
Added in version 3.1.0.
- class hoomd.md.pair.Mie(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Mie pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
Mie
computes the Mie pair force on every particle in the simulation state.\[U(r) = \left( \frac{n}{n-m} \right) {\left( \frac{n}{m} \right)}^{\frac{m}{n-m}} \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{n} - \left( \frac{\sigma}{r} \right)^{m} \right]\]Example:
nl = nlist.Cell() mie = pair.Mie(nlist=nl, default_r_cut=3.0) mie.params[('A', 'A')] = dict(epsilon=1.0, sigma=1.0, n=12, m=6) mie.r_cut[('A', 'A')] = 2**(1.0/6.0) mie.r_on[('A', 'A')] = 2.0 mie.params[(['A', 'B'], ['C', 'D'])] = dict(epsilon=1.5, sigma=2.0)
- params¶
The potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - \(\varepsilon\) \([\mathrm{energy}]\)sigma
(float
, required) - \(\sigma\) \([\mathrm{length}]\)n
(float
, required) - \(n\) \([\mathrm{dimensionless}]\)m
(float
, required) - \(m\) \([\mathrm{dimensionless}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.Moliere(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Moliere pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
Moliere
computes the Moliere pair force on every particle in the simulation state:\[\begin{split}U(r) = \frac{Z_i Z_j e^2}{4 \pi \epsilon_0 r_{ij}} \left[ 0.35 \exp \left( -0.3 \frac{r_{ij}}{a_F} \right) + \\ 0.55 \exp \left( -1.2 \frac{r_{ij}}{a_F} \right) + 0.10 \exp \left( -6.0 \frac{r_{ij}}{a_F} \right) \right]\end{split}\]Where each parameter is defined as:
\(Z_i\) - Z_i - Atomic number of species i \([\mathrm{dimensionless}]\)
\(Z_j\) - Z_j - Atomic number of species j \([\mathrm{dimensionless}]\)
\(e\) - elementary_charge - The elementary charge \([\mathrm{charge}]\)
\(a_F = \frac{0.8853 a_0}{\left( \sqrt{Z_i} + \sqrt{Z_j} \right)^{2/3}}\), where \(a_0\) is the Bohr radius \([\mathrm{length}]\)
Example:
nl = nlist.Cell() moliere = pair.Moliere(default_r_cut = 3.0, nlist=nl) Zi = 54 Zj = 7 e = 1 a0 = 1 aF = 0.8853 * a0 / (np.sqrt(Zi) + np.sqrt(Zj))**(2/3) moliere.params[('A', 'B')] = dict(qi=Zi*e, qj=Zj*e, aF=aF)
- params¶
The potential parameters. The dictionary has the following keys:
qi
(float
, required) - \(q_i = Z_i \frac{e}{\sqrt{4 \pi \epsilon_0}}\) \([\mathrm{charge}]\)qj
(float
, required) - \(q_j = Z_j \frac{e}{\sqrt{4 \pi \epsilon_0}}\) \([\mathrm{charge}]\)aF
(float
, required) - \(a_F = \frac{0.8853 a_0}{\left( \sqrt{Z_i} + \sqrt{Z_j} \right)^{2/3}}\) \([\mathrm{length}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.Morse(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Morse pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
Morse
computes the Morse pair force on every particle in the simulation state:\[U(r) = D_0 \left[ \exp \left(-2\alpha\left( r-r_0\right)\right) -2\exp \left(-\alpha\left(r-r_0\right) \right) \right]\]Example:
nl = nlist.Cell() morse = pair.Morse(default_r_cut=3.0, nlist=nl) morse.params[('A', 'A')] = dict(D0=1.0, alpha=3.0, r0=1.0) morse.r_cut[('A', 'B')] = 3.0
- params¶
The potential parameters. The dictionary has the following keys:
D0
(float
, required) - depth of the potential at its minimum \(D_0\) \([\mathrm{energy}]\)alpha
(float
, required) - the width of the potential well \(\alpha\) \([\mathrm{length}^{-1}]\)r0
(float
, required) - position of the minimum \(r_0\) \([\mathrm{length}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.OPP(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Oscillating pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
OPP
computes the oscillating pair force on all particles in the simulation state:\[U(r) = C_1 r^{-\eta_1} + C_2 r^{-\eta_2} \cos{\left(k r - \phi\right)}\]The potential was introduced in Marek Mihalkovič and C. L. Henley 2012.
Example:
nl = nlist.Cell() opp = pair.OPP(nl, default_r_cut=3.0) opp.params[('A', 'A')] = { 'C1': 1., 'C2': 1., 'eta1': 15, 'eta2': 3, 'k': 1.0, 'phi': 3.14} opp.r_cut[('A', 'B')] = 3.0
- params¶
The OPP potential parameters. The dictionary has the following keys:
C1
(float
, required) - Energy scale of the first term \(C_1\) \([\mathrm{energy}]\)C2
(float
, required) - Energy scale of the second term \(C_2\) \([\mathrm{energy}]\)eta1
(float
, required) - The inverse power to take \(r\) to in the first term, \(\eta_1\) \([\mathrm{dimensionless}]\).eta2
(float
, required) - The inverse power to take \(r\) to in the second term \(\eta_2\) \([\mathrm{dimensionless}]\).k
(float
, required) - oscillation frequency \(k\) \([\mathrm{length}^{-1}]\)phi
(float
, required) - potential phase shift \(\phi\) \([\mathrm{dimensionless}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.Pair(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Force
Base class pair force.
Pair
is the base class for all pair forces.Warning
This class should not be instantiated by users. The class can be used for
isinstance
orissubclass
checks.- r_cut¶
Cuttoff radius beyond which the energy and force are 0 \([\mathrm{length}]\). Optional: defaults to the value
default_r_cut
specified on construction.Type:
TypeParameter
[tuple
[particle_type
,particle_type
],float
])
- r_on¶
Radius at which the smoothing modification to the potential starts \([\mathrm{length}]\). Optional: defaults to the value
default_r_on
specified on construction.Type:
TypeParameter
[tuple
[particle_type
,particle_type
],float
])
- nlist¶
Neighbor list used to compute the pair force.
- compute_energy(tags1, tags2)¶
Compute the energy between two sets of particles.
- Parameters:
tags1 (
ndarray<int32>
) – a numpy array of particle tags in the first group.tags2 (
ndarray<int32>
) – a numpy array of particle tags in the second group.
\[U = \sum_{i \in \mathrm{tags1}, j \in \mathrm{tags2}} V_{ij}(r)\]where \(V_{ij}(r)\) is the pairwise energy between two particles \(i\) and \(j\).
Assumed properties of the sets tags1 and tags2 are:
tags1 and tags2 are disjoint
all elements in tags1 and tags2 are unique
tags1 and tags2 are contiguous numpy arrays of dtype int32
None of these properties are validated.
Examples:
tags=numpy.linspace(0,N-1,1, dtype=numpy.int32) # computes the energy between even and odd particles U = mypair.compute_energy(tags1=numpy.array(tags[0:N:2]), tags2=numpy.array(tags[1:N:2]))
- class hoomd.md.pair.ReactionField(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Onsager reaction field pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
ReactionField
computes the Onsager reaction field pair force on all particles in the simulation state.Reaction field electrostatics is an approximation to the screened electrostatic interaction, which assumes that the medium can be treated as an electrostatic continuum of dielectric constant \(\epsilon_{RF}\) outside the cutoff sphere of radius \(r_{\mathrm{cut}}\). See: Barker et. al. 1973.
By default (
use_charge=False
), the reaction field potential ignores the particle charges. Two parameters, \(\varepsilon\) and \(\epsilon_{RF}\) are needed. If \(\epsilon_{RF}\) is specified as zero, it will represent infinity:\[U(r) = \varepsilon \left[ \frac{1}{r} + \frac{(\epsilon_{RF}-1) r^2}{(2 \epsilon_{RF} + 1) r_c^3} \right]\]When use_charge is set to
True
, the following formula is evaluated instead:\[U(r) = q_i q_j \varepsilon \left[ \frac{1}{r} + \frac{(\epsilon_{RF}-1) r^2}{(2 \epsilon_{RF} + 1) r_c^3} \right]\]where \(q_i\) and \(q_j\) are the charges of the particle pair.
Example:
nl = nlist.Cell() reaction_field = pair.reaction_field(nl, default_r_cut=3.0) reaction_field.params[('A', 'B')] = dict(epsilon=1.0, eps_rf=1.0) reaction_field.params[('B', 'B')] = dict( epsilon=1.0, eps_rf=0.0, use_charge=True)
- params¶
The potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - \(\varepsilon\) \([\mathrm{energy} \cdot \mathrm{length}]\)eps_rf
(float
, required) - \(\epsilon_{RF}\) \([\mathrm{dimensionless}]\)use_charge
(bool
, optional) - evaluate pair force using particle charges (default: False)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.TWF(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Pair potential model for globular proteins.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
TWF
computes the Ten-wolde Frenkel potential on all particles in the simulation state:\[U(r) = \frac{4 \epsilon}{\alpha^2} {\left[ {\left(\frac{\sigma^2}{r^2} - 1 \right)}^{-6} - \alpha {\left(\frac{\sigma^2}{r^2} - 1 \right)}^{-3}\right]}\]The potential was introdcued in Pieter Rein ten Wolde and Daan Frenkel 1997.
Example:
nl = nlist.Cell() twf = hoomd.md.pair.TWF(nl, default_r_cut=3.0) twf.params[('A', 'A')] = {'sigma': 1.0, 'epsilon': 1.0, 'alpha': 50.0} twf.r_cut[('A', 'B')] = 3.0
- params¶
The LJ potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - energy parameter \(\varepsilon\) \([\mathrm{energy}]\)sigma
(float
, required) - particle size \(\sigma\) \([\mathrm{length}]\)alpha
(float
, required) - controls well-width \(\alpha\) \([\mathrm{dimensionless}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.Table(nlist, default_r_cut=None)¶
Bases:
Pair
Tabulated pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
Table
computes the tabulated pair force on every particle in the simulation state.The force \(\vec{F}\) is:
\[\begin{split}\vec{F}(\vec{r}) = \begin{cases} 0 & r < r_{\mathrm{min}} \\ F_\mathrm{table}(r)\hat{r} & r_{\mathrm{min}} \le r < r_{\mathrm{cut}} \\ 0 & r \ge r_{\mathrm{cut}} \\ \end{cases}\end{split}\]and the potential \(U(r)\) is:
\[\begin{split}U(r) = \begin{cases} 0 & r < r_{\mathrm{min}} \\ U_\mathrm{table}(r) & r_{\mathrm{min}} \le r < r_{\mathrm{cut}} \\ 0 & r \ge r_{\mathrm{cut}} \\ \end{cases}\end{split}\]where \(\vec{r}\) is the vector pointing from one particle to the other in the pair,
r_min
is defined inparams
, andr_cut
is defined inPair.r_cut
.Provide \(F_\mathrm{table}(r)\) and \(U_\mathrm{table}(r)\) on evenly spaced grid points points between \(r_{\mathrm{min}}\) and \(r_{\mathrm{cut}}\).
Table
linearly interpolates values when \(r\) lies between grid points and between the last grid point and \(r=r_{\mathrm{cut}}\). The force must be specificed commensurate with the potential: \(F = -\frac{\partial U}{\partial r}\).Table
does not support energy shifting or smoothing modes.Note
For potentials that diverge near r=0, to set r_min to a non-zero value.
Note
The implicitly defined \(r\) values are those that would be returned by
numpy.linspace(r_min, r_cut, len(U), endpoint=False)
.Tip
Define non-interacting potentials with:
table.params[(type1, type2)] = dict(r_min=0, U=[0], F=[0]) table.r_cut[(type1, type2)] = 0
There must be at least one element in U and F, and the
r_cut
value of 0 disables the interaction entirely.- params¶
The potential parameters. The dictionary has the following keys:
r_min
(float
, required) - the minimum distance to apply the tabulated potential, corresponding to the first element of the energy and force arrays \([\mathrm{length}]\).U
((N,)numpy.ndarray
offloat
, required) - the tabulated energy values \([\mathrm{energy}]\).F
((N,)numpy.ndarray
offloat
, required) - the tabulated force values \([\mathrm{force}]\). Must have the same length asU
.
- Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.Yukawa(nlist, default_r_cut=None, default_r_on=0.0, mode='none')¶
Bases:
Pair
Yukawa pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
mode (str) – Energy shifting/smoothing mode.
Yukawa
computes the Yukawa pair force son every particle in the simulation state:\[U(r) = \varepsilon \frac{ \exp \left( -\kappa r \right) }{r}\]Example:
nl = nlist.Cell() yukawa = pair.Yukawa(default_r_cut=3.0, nlist=nl) yukawa.params[('A', 'A')] = dict(epsilon=1.0, kappa=1.0) yukawa.r_cut[('A', 'B')] = 3.0
- params¶
The Yukawa potential parameters. The dictionary has the following keys:
epsilon
(float
, required) - energy parameter \(\varepsilon\) \([\mathrm{energy}]\)kappa
(float
, required) - scaling parameter \(\kappa\) \([\mathrm{length}^{-1}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
- class hoomd.md.pair.ZBL(nlist, default_r_cut=None, default_r_on=0.0)¶
Bases:
Pair
ZBL pair force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list.
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).
ZBL
computes the Ziegler-Biersack-Littmark pair force on every particle in the simulation state:\[\begin{split}U(r) = \frac{Z_i Z_j e^2}{4 \pi \epsilon_0 r_{ij}} \left[ 0.1818 \exp \left( -3.2 \frac{r_{ij}}{a_F} \right) \right. \\ + 0.5099 \exp \left( -0.9423 \frac{r_{ij}}{a_F} \right) \\ + 0.2802 \exp \left( -0.4029 \frac{r_{ij}}{a_F} \right) \\ + \left. 0.02817 \exp \left( -0.2016 \frac{r_{ij}}{a_F} \right) \right]\end{split}\]Where each parameter is defined as:
\(Z_i\) - Z_i - Atomic number of species i \([\mathrm{dimensionless}]\)
\(Z_j\) - Z_j - Atomic number of species j \([\mathrm{dimensionless}]\)
\(e\) - elementary_charge - The elementary charge \([\mathrm{charge}]\)
\(a_F = \frac{0.8853 a_0}{ Z_i^{0.23} + Z_j^{0.23} }\), where \(a_0\) is the Bohr radius \([\mathrm{length}]\)
Example:
nl = nlist.Cell() zbl = pair.ZBL(default_r_cut=3.0, nlist=nl) Zi = 54 Zj = 7 e = 1 a0 = 1 aF = 0.8853 * a0 / (Zi**(0.23) + Zj**(0.23)) zbl.params[('A', 'B')] = dict(qi=Zi*e, qj=Zj*e, aF=aF)
- params¶
The ZBL potential parameters. The dictionary has the following keys:
q_i
(float
, required) - \(q_i=Z_i \frac{e}{\sqrt{4 \pi \epsilon_0}}\) \([\mathrm{charge}]\)q_j
(float
, required) - \(q_j=Z_j \frac{e}{\sqrt{4 \pi \epsilon_0}}\) \([\mathrm{charge}]\)a_F
(float
, required) - \(a_F = \frac{0.8853 a_0}{ Z_i^{0.23} + Z_j^{0.23} }\) \([\mathrm{length}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
Modules