md.pair

Overview

Buckingham

Buckingham pair force.

DLVO

DLVO colloidal interaction.

DPD

Dissipative Particle Dynamics.

DPDLJ

Dissipative Particle Dynamics with the LJ conservative force.

DPDConservative

DPD Conservative pair force.

Ewald

Ewald pair force.

ExpandedLJ

Expanded Lennard-Jones pair force.

ExpandedMie

Expanded Mie pair force.

ForceShiftedLJ

Force-shifted Lennard-Jones pair force.

Fourier

Fourier pair force.

Gauss

Gaussian pair force.

LJ

Lennard-Jones pair force.

LJ1208

Lennard-Jones 12-8 pair force.

LJ0804

Lennard-Jones 8-4 pair force.

LJGauss

Lennard-Jones-Gauss pair potential.

Mie

Mie pair force.

Morse

Morse pair force.

Moliere

Moliere pair force.

OPP

Oscillating pair force.

Pair

Base class pair force.

ReactionField

Onsager reaction field pair force.

Table

Tabulated pair force.

TWF

Pair potential model for globular proteins.

Yukawa

Yukawa pair force.

ZBL

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:

\[U_\mathrm{pair,total} = \frac{1}{2} \sum_{i=0}^\mathrm{N_particles-1} \sum_{j \ne i, (i,j) \notin \mathrm{exclusions}} U_\mathrm{pair}(r_{ij})\]

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:

\[\begin{split}\vec{F} = \begin{cases} -\nabla U_\mathrm{pair}(r) & r < r_{\mathrm{cut}} \\ 0 & r \ge r_{\mathrm{cut}} \\ \end{cases}\end{split}\]

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\):

\[U_i = \frac{1}{2} \sum_{j \ne i, (i,j) \notin \mathrm{exclusions}} U_\mathrm{pair}(r_{ij}) [r_{ij} < r_\mathrm{cut}]\]

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):

\[\begin{split}U_\mathrm{pair}(r) = \begin{cases} U_(r) & \mathrm{mode\ is\ no\_shift} \\ U(r) - U(r_{\mathrm{cut}}) & \mathrm{mode\ is\ shift} \\ S(r) \cdot U_(r) & \mathrm{mode\ is\ xplor} \land r_{\mathrm{on}} < r_{\mathrm{cut}} \\ U(r) - U(r_{\mathrm{cut}}) & \mathrm{mode\ is\ xplor} \land r_{\mathrm{on}} \ge r_{\mathrm{cut}} \end{cases}\end{split}\]

where \(S(r)\) is the XPLOR smoothing function:

\[\begin{split}S(r) = \begin{cases} 1 & r < r_{\mathrm{on}} \\ \frac{(r_{\mathrm{cut}}^2 - r^2)^2 \cdot (r_{\mathrm{cut}}^2 + 2r^2 - 3r_{\mathrm{on}}^2)}{(r_{\mathrm{cut}}^2 - r_{\mathrm{on}}^2)^3} & r_{\mathrm{on}} \le r \le r_{\mathrm{cut}} \\ 0 & r > r_{\mathrm{cut}} \\ \end{cases}\end{split}\]

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:

\[\Delta P = \frac{-2\pi}{3} \sum_{i=1}^{n} \rho_i \sum_{j=1}^{n} \rho_j \int_{r_\mathrm{cut}}^{\infty} \left( r \frac{\mathrm{d}U_{ij}(r)}{\mathrm{d}r} \right) r^2 \mathrm{d}r\]

and

\[\Delta E = 2\pi \sum_{i=1}^{n} N_i \sum_{j=1}^{n} \rho_j \int_{r_\mathrm{cut}}^{\infty} U_{ij}(r) r^2 \mathrm{d}r,\]

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)
params

The potential parameters. The dictionary has the following keys:

  • A (float, required) - \(A\) \([\mathrm{energy}]\)

  • rho (float, required) - \(\rho\) \([\mathrm{length}]\)

  • C (float, required) - \(C\) \([\mathrm{energy}]\)

Type: TypeParameter [tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]

mode

Energy shifting/smoothing mode: "none" or "shift".

Type: str

class hoomd.md.pair.DPD(nlist, kT, default_r_cut=None)

Bases: Pair

Dissipative Particle Dynamics.

Parameters

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.NVE integration method along with DPD forces. Use of the DPD thermostat pair force with other integrators will result in nonphysical behavior. To use DPD 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)
params

The force parameters. The dictionary has the following keys:

  • A (float, required) - \(A\) \([\mathrm{force}]\)

  • gamma (float, required) - \(\gamma\) \([\mathrm{mass} \cdot \mathrm{time}^{-1}]\)

Type: TypeParameter [tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none".

Type: str

class hoomd.md.pair.DPDConservative(nlist, default_r_cut=None)

Bases: Pair

DPD Conservative pair force.

Parameters

DPDConservative computes the conservative part of the DPD 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]

mode

Energy shifting/smoothing mode: "none".

Type: str

class hoomd.md.pair.DPDLJ(nlist, kT, default_r_cut=None, mode='none')

Bases: Pair

Dissipative Particle Dynamics with the LJ conservative force.

Parameters

DPDLJ computes the DPD thermostat combined with the LJ 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.NVE integration method along with DPD 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]

mode

Energy shifting/smoothing mode: "none" or "shift".

Type: str

class hoomd.md.pair.Ewald(nlist, default_r_cut=None)

Bases: Pair

Ewald pair force.

Parameters

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 of Ewald and md.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]

mode

Energy shifting/smoothing mode: "none".

Type: str

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]\]

Note

To replicate the behavior of the SLJ potential in HOOMD-blue v2, set hoomd.md.pair.Pair.r_cut to r_cut_unshifted + delta.

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

class hoomd.md.pair.ForceShiftedLJ(nlist, default_r_cut=None)

Bases: Pair

Force-shifted Lennard-Jones pair force.

Parameters

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)
params

The potential parameters. The dictionary has the following keys:

  • epsilon (float, required) - \(\varepsilon\) \([\mathrm{energy}]\)

  • sigma (float, required) - \(\sigma\) \([\mathrm{length}]\)

Type: TypeParameter [tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none".

Type: str

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

class hoomd.md.pair.Gauss(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.

Gauss 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.Gauss(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\) \([\mathrm{length}]\)

Type: TypeParameter [tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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
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}]\)

Type: TypeParameter [tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

tail_correction

Whether to apply the isotropic integrated long range tail correction.

Type: bool

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
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}]\)

Type: TypeParameter [tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]\]
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}]\)

Type: TypeParameter [tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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:

  • epsilon (float, required) - energy parameter \(\varepsilon\) \([\mathrm{energy}]\)

  • sigma (float, required) - Gaussian width \(\sigma\) \([\mathrm{length}]\)

  • r0 (float, required) - Gaussian center \(r_0\) \([\mathrm{length}]\)

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}

New 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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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 or issubclass 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])

mode

mode, optional: defaults to "none". Possible values: "none", "shift", "xplor"

Type: str

nlist

Neighbor list used to compute the pair force.

Type: hoomd.md.nlist.NeighborList

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

class hoomd.md.pair.Table(nlist, default_r_cut=None)

Bases: Pair

Tabulated pair force.

Parameters

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 in params, and r_cut is defined in Pair.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 of float, required) - the tabulated energy values \([\mathrm{energy}]\).

  • F ((N,) numpy.ndarray of float, required) - the tabulated force values \([\mathrm{force}]\). Must have the same length as U.

Type

TypeParameter [ tuple [particle_type, particle_type], dict]

mode

Energy shifting/smoothing mode: "none".

Type

str

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]

mode

Energy shifting/smoothing mode: "none", "shift", or "xplor".

Type: str

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]

mode

Energy shifting/smoothing mode: "none".

Type: str

Modules