md.pair

Overview

Pair

Common pair potential documentation.

Buckingham

Buckingham pair potential.

DLVO

DLVO colloidal interaction.

DPD

Dissipative Particle Dynamics.

DPDLJ

Dissipative Particle Dynamics with a LJ conservative force.

DPDConservative

DPD Conservative pair force.

Ewald

Ewald pair potential.

ExpandedMie

Expanded Mie pair potential.

ForceShiftedLJ

Force-shifted Lennard-Jones pair potential.

Fourier

Fourier pair potential.

Gauss

Gaussian pair potential.

LJ

Lennard-Jones pair potential.

LJ1208

Lennard-Jones 12-8 pair potential.

LJ0804

Lennard-Jones 8-4 pair potential.

Mie

Mie pair potential.

Morse

Morse pair potential.

Moliere

Moliere pair potential.

OPP

Oscillating pair potential.

ReactionField

Onsager reaction field pair potential.

SLJ

Shifted Lennard-Jones pair potential.

Table

Tabulated pair potential.

TWF

Pair potential model for globular proteins.

Yukawa

Yukawa pair potential.

ZBL

ZBL pair potential.

Details

Pair Potentials for molecular dynamics.

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

Buckingham pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Buckingham pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{Buckingham}}(r) = & A \exp\left(-\frac{r}{\rho}\right) - \frac{C}{r^6}; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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

DLVO colloidal interaction.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a DLVO dispersion and electrostatic interaction should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} 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))}; & r < (r_{\mathrm{cut}} + \Delta) \\ = & 0; & r \ge (r_{\mathrm{cut}} + \Delta) \end{eqnarray*}

where \(a_i\) is the radius of particle \(i\), \(\Delta = (d_i + d_j)/2\) and \(d_i\) is the diameter of particle \(i\).

The first term corresponds to the attractive van der Waals interaction with \(A\) being the Hamaker constant, the second term to the repulsive double-layer interaction between two spherical surfaces with Z proportional to the surface electric potential. See Israelachvili 2011, pp. 317.

The DLVO potential does not need charge, but does need diameter. See SLJ for an explanation on how diameters are handled in the neighbor lists.

Due to the way that DLVO modifies the cutoff condition, it will not function properly with the xplor shifting mode. See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

params

The potential parameters. The dictionary has the following keys:

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

  • kappa (float, required) - scaling parameter \(\kappa\) \([\mathrm{length}^{-1}]\)

  • Z (float, required) - \(Z\) \([\mathrm{length}^{-1}]\)

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

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

Example:

nl = nlist.cell()
dlvo = hoomd.md.pair.DLVO(nlist=nl)
dlvo.params[('A', 'A')] = {"epsilon": 1.0, "kappa": 1.0}
dlvo.params[('A', 'B')] = {"epsilon": 2.0, "kappa": 0.5,}
dlvo.params[(['A', 'B'], ['C', 'D'])] = {"epsilon": 0.5, "kappa": 3.0}
class hoomd.md.pair.DPD(nlist, kT, default_r_cut=None, default_r_on=0.0)

Dissipative Particle Dynamics.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list

  • kT (hoomd.variant or float) – Temperature of thermostat \([\mathrm{energy}]\).

  • default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).

  • default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).

DPD specifies that a DPD pair force should be applied between every non-excluded particle pair in the simulation, including an interaction potential, pairwise drag force, and pairwise random force. See Groot and Warren 1997.

\begin{eqnarray*} F = F_{\mathrm{C}}(r) + F_{\mathrm{R,ij}}(r_{ij}) + F_{\mathrm{D,ij}}(v_{ij}) \\ \end{eqnarray*}
\begin{eqnarray*} 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) \\ \end{eqnarray*}
\begin{eqnarray*} w(r_{ij}) = &\left( 1 - r/r_{\mathrm{cut}} \right); & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

where \(\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 in HOOMD-blue. 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, an hoomd.md.methods.NVE integrator must be applied to the system and the user must specify a temperature. Use of the dpd thermostat pair force with other integrators will result in unphysical behavior. To use pair.dpd with a different conservative potential than \(F_C\), set A to zero and define the conservative pair potential separately. Note that DPD thermostats are often defined in terms of \(\sigma\) where \(\sigma = \sqrt{2k_b\gamma T}\).

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]

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, default_r_on=0.0)

DPD Conservative pair force.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list.

  • default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).

  • default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).

DPDConservative specifies the conservative part of the DPD pair potential should be applied between every non-excluded particle pair in the simulation. No thermostat (e.g. Drag Force and Random Force) is applied, as is in DPD.

\begin{eqnarray*} V_{\mathrm{DPD-C}}(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); & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

DPDConservative does not implement any energy shift / smoothing modes due to the function of the force.

params

The potential parameters. The dictionary has the following keys:

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

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

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)
class hoomd.md.pair.DPDLJ(nlist, kT, default_r_cut=None, default_r_on=0.0, mode='none')

Dissipative Particle Dynamics with a LJ conservative force.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list.

  • kT (hoomd.variant or float) – Temperature of thermostat \([\mathrm{energy}]\).

  • default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).

  • default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).

  • mode (str) – Energy shifting mode.

DPDLJ specifies that a DPD thermostat and a Lennard-Jones pair potential should be applied between every non-excluded particle pair in the simulation.

C. L. Phillips et. al. 2011 describes the DPD implementation details in HOOMD-blue. Cite it if you utilize the DPD functionality in your work.

\begin{eqnarray*} F = F_{\mathrm{C}}(r) + F_{\mathrm{R,ij}}(r_{ij}) + F_{\mathrm{D,ij}}(v_{ij}) \\ \end{eqnarray*}
\begin{eqnarray*} F_{\mathrm{C}}(r) = & \partial V_{\mathrm{LJ}} / \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) \\ \end{eqnarray*}
\begin{eqnarray*} V_{\mathrm{LJ}}(r) = & 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}
\begin{eqnarray*} w(r_{ij}) = &\left( 1 - r/r_{\mathrm{cut}} \right); & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

where \(\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].

To use the DPD thermostat, an hoomd.md.methods.NVE integrator must be applied to the system and the user must specify a temperature. Use of the dpd thermostat pair force with other integrators will result in unphysical behavior.

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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)
class hoomd.md.pair.Ewald(nlist, default_r_cut=None, default_r_on=0.0)

Ewald pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list.

  • default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).

  • default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).

Ewald specifies that a Ewald pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{ewald}}(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]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

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.

See Pair for details on how forces are calculated. Note Ewald does not support energy shifting or smoothing.

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]

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
class hoomd.md.pair.ExpandedMie(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Expanded Mie pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list.

  • default_r_cut (float) – Default cutoff radius (in distance units).

  • default_r_on (float) – Default turn-on radius (in distance units).

  • mode (str) – Energy shifting/smoothing mode.

ExpandedMie specifies that a radially shifted Mie pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{mie}}(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]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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}
class hoomd.md.pair.ForceShiftedLJ(nlist, default_r_cut=None, default_r_on=0.0)

Force-shifted Lennard-Jones pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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.

ForceShiftedLJ specifies that a modified Lennard-Jones pair force should be applied between non-excluded particle pair in the simulation. 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. This potential can be used as a substitute for LJ, when the exact analytical form of the latter is not required but a smaller cut-off radius is desired for computational efficiency. See Toxvaerd et. al. 2011 for a discussion of this potential.

\begin{eqnarray*} V(r) = & 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right] + \Delta V(r); & r < r_{\mathrm{cut}}\\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}
\[\Delta V(r) = -(r - r_{\mathrm{cut}}) \frac{\partial V_{\mathrm{LJ}}}{\partial r}(r_{\mathrm{cut}})\]

See Pair for details on how forces are calculated.

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]

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

Fourier pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Fourier pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{Fourier}}(r) = & \frac{1}{r^{12}} + \frac{1}{r^2}\sum_{n=1}^4 [a_n cos(\frac{n \pi r}{r_{cut}}) + b_n sin(\frac{n \pi r}{r_{cut}})]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*} where: \begin{eqnarray*} a_1 = \sum_{n=2}^4 (-1)^n a_n \end{eqnarray*} \begin{eqnarray*} b_1 = \sum_{n=2}^4 n (-1)^n b_n \end{eqnarray*} is calculated to enforce close to zero value at r_cut.

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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])
class hoomd.md.pair.Gauss(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Gaussian pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Gaussian pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{gauss}}(r) = & \varepsilon \exp \left[ -\frac{1}{2} \left( \frac{r}{\sigma} \right)^2 \right]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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
class hoomd.md.pair.LJ(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Lennard-Jones pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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.

LJ specifies that a Lennard-Jones pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{LJ}}(r) = & 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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

Lennard-Jones 8-4 pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Lennard-Jones 8-4 pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{LJ}}(r) = & 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{8} - \left( \frac{\sigma}{r} \right)^{4} \right]; & r < r_{\mathrm{cut}} \\ = & 0: & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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

Lennard-Jones 12-8 pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Lennard-Jones 12-8 pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{LJ}}(r) = & 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{8} \right]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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)
class hoomd.md.pair.Mie(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Mie pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Mie pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{mie}}(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]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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)
class hoomd.md.pair.Moliere(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Moliere pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Moliere type pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{Moliere}}(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]; & r < r_{\mathrm{cut}} \\ = & 0; & r > r_{\mathrm{cut}} \\ \end{eqnarray*}

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

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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)
class hoomd.md.pair.Morse(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Morse pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Morse pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{morse}}(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]; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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
class hoomd.md.pair.OPP(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Oscillating pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that an oscillating pair potential should be applied between every non-excluded particle pair in the simulation. The OPP potential can be used to model metallic interactions.

\begin{equation*} V_{\mathrm{OPP}}(r) = C_1 r^{-\eta_1} + C_2 r^{-\eta_2} \cos{\left(k r - \phi\right)} \end{equation*}

See Pair for details on how forces are calculate. Note OPP does not support energy shifting or smoothing.

The potential comes from Marek Mihalkovič and C. L. Henley 2012.

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]

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
class hoomd.md.pair.Pair(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Common pair potential documentation.

Users should not invoke Pair directly. It is a base command that provides common features to all standard pair forces. Common documentation for all pair potentials is documented here.

All pair force commands specify that a given potential energy and force be computed on all non-excluded particle pairs in the system within a short range cutoff distance \(r_{\mathrm{cut}}\).

The force \(\vec{F}\) applied between each pair of particles is:

\begin{eqnarray*} \vec{F} = & -\nabla V(r); & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

where \(\vec{r}\) is the vector pointing from one particle to the other in the pair, and \(V(r)\) is chosen by a mode switch:

\begin{eqnarray*} V(r) = & V_{\mathrm{pair}}(r); & \mathrm{mode\ is\ no\_shift} \\ = & V_{\mathrm{pair}}(r) - V_{\mathrm{pair}}(r_{\mathrm{cut}}); & \mathrm{mode\ is\ shift} \\ = & S(r) \cdot V_{\mathrm{pair}}(r); & \mathrm{mode\ is\ xplor\ and\ } r_{\mathrm{on}} < r_{\mathrm{cut}} \\ = & V_{\mathrm{pair}}(r) - V_{\mathrm{pair}}(r_{\mathrm{cut}}); & \mathrm{mode\ is\ xplor\ and\ } r_{\mathrm{on}} \ge r_{\mathrm{cut}} \end{eqnarray*}

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

\begin{eqnarray*} S(r) = & 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{eqnarray*}

and \(V_{\mathrm{pair}}(r)\) is the specific pair potential chosen by the respective command.

Enabling the XPLOR smoothing function \(S(r)\) results in 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, so it can be set to only slightly modify the tail of the potential. It is suggested that you plot your potentials with various values of \(r_{\mathrm{on}}\) in order to find a good balance between a smooth potential function and minimal modification of the original \(V_{\mathrm{pair}}(r)\). A good value for the LJ potential is \(r_{\mathrm{on}} = 2 \cdot \sigma\).

The split smoothing / shifting of the potential when the mode is xplor is designed for use in mixed WCA / LJ systems. 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.

The following coefficients must be set per unique pair of particle types. See hoomd.md.pair for information on how to set coefficients.

r_cut
r_cut \([\mathrm{length}]\), optional: defaults to the

value default_r_cut specified on construction.

Type: TypeParameter [tuple [particle_type, particle_type], float])

r_on
r_on \([\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

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]))
property nlist

Neighbor list used to compute the pair potential.

class hoomd.md.pair.ReactionField(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Onsager reaction field pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that an Onsager reaction field pair potential should be applied between every non-excluded particle pair in the simulation.

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.

\[V_{\mathrm{RF}}(r) = \varepsilon \left[ \frac{1}{r} + \frac{(\epsilon_{RF}-1) r^2}{(2 \epsilon_{RF} + 1) r_c^3} \right]\]

By default, the reaction field potential does not require charge or diameter to be set. Two parameters, \(\varepsilon\) and \(\epsilon_{RF}\) are needed. If \(\epsilon_{RF}\) is specified as zero, it will represent infinity.

If use_charge is set to True, the following formula is evaluated instead:

\[V_{\mathrm{RF}}(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.

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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 (Boolean Operations and, or, not, optional) - evaluate pair potential using particle charges (default: False)

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

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)
class hoomd.md.pair.SLJ(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Shifted Lennard-Jones pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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.

SLJ specifies that a shifted Lennard-Jones type pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{SLJ}}(r) = & 4 \varepsilon \left[ \left( \frac{\sigma}{r - \Delta} \right)^{12} - \left( \frac{\sigma}{r - \Delta} \right)^{6} \right]; & r < (r_{\mathrm{cut}} + \Delta) \\ = & 0; & r \ge (r_{\mathrm{cut}} + \Delta) \\ \end{eqnarray*}

where \(\Delta = (d_i + d_j)/2 - 1\) and \(d_i\) is the diameter of particle \(i\).

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

Attention

Due to the way that SLJ modifies the cutoff criteria, a smoothing mode of xplor is not supported.

Set the max_diameter property of the neighbor list object to the largest particle diameter in the system (where diameter is a per-particle property of the same name in hoomd.State).

Warning

Failure to set max_diameter will result in missing pair interactions.

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]

Example:

nl = nlist.Cell()
nl.max_diameter = 2.0
slj = pair.SLJ(default_r_cut=3.0, nlist=nl)
slj.params[('A', 'B')] = dict(epsilon=2.0, r_cut=3.0)
slj.r_cut[('B', 'B')] = 2**(1.0/6.0)
class hoomd.md.pair.TWF(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Pair potential model for globular proteins.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list.

  • default_r_cut (float) – Default cutoff radius (in distance units).

  • default_r_on (float) – Default turn-on radius (in distance units).

  • mode (str) – Energy shifting/smoothing mode.

This potential was introduced by Ten-wolde and Daan Frenkel in 1997 for studying globular protein crystallization. The potential has the following form:

\begin{equation} V_{\mathrm{TWF}}(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]} \end{equation}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

The potential comes from Pieter Rein ten Wolde and Daan Frenkel 1997.

params

The LJ potential parameters. The dictionary has the following keys:

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

  • sigma (float, required) - particle size \(\sigma\) \([length]\)

  • alpha (float, required) - controls well-width \(\alpha\) \([dimensionless]\)

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

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
class hoomd.md.pair.Table(nlist, default_r_cut=None, default_r_on=0.0)

Tabulated pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list

  • default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).

  • default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).

Table specifies that a tabulated pair potential should be applied between every non-excluded particle pair in the simulation in the range \([r_{\mathrm{min}}, r_{\mathrm{cut}})\)

Note

For potentials that diverge near r=0, to set r_min to a non-zero value.

The force \(\vec{F}\) is:

\begin{eqnarray*} \vec{F}(\vec{r}) = & 0; & r < r_{\mathrm{min}} \\ = & F(r)\hat{r}; & r_{\mathrm{min}} \le r < r_{\mathrm{max}} \\ = & 0; & r \ge r_{\mathrm{max}} \\ \end{eqnarray*}

and the potential \(V(r)\) is:

\begin{eqnarray*} V(r) = & 0; & r < r_{\mathrm{min}} \\ = & V(r); & r_{\mathrm{min}} \le r < r_{\mathrm{max}} \\ = & 0; & r \ge r_{\mathrm{max}} \\ \end{eqnarray*}

where \(\vec{r}\) is the vector pointing from one particle to the other in the pair.

Provide \(F(r)\) and \(V(r)\) on an evenly space set of 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 V}{\partial r}\).

Table does not support energy shifting or smoothing modes.

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

  • V ((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 V.

Type

TypeParameter [ tuple [particle_type, particle_type], dict]

Note

The implicitly defined \(r\) values are those that would be returned by numpy.linspace(r_min, r_cut, len(V), endpoint=False).

class hoomd.md.pair.Yukawa(nlist, default_r_cut=None, default_r_on=0.0, mode='none')

Yukawa pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – 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 specifies that a Yukawa pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{yukawa}}(r) = & \varepsilon \frac{ \exp \left( -\kappa r \right) }{r}; & r < r_{\mathrm{cut}} \\ = & 0; & r \ge r_{\mathrm{cut}} \\ \end{eqnarray*}

See Pair for details on how forces are calculated and the available energy shifting and smoothing modes.

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]

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
class hoomd.md.pair.ZBL(nlist, default_r_cut=None, default_r_on=0.0)

ZBL pair potential.

Parameters
  • nlist (hoomd.md.nlist.NList) – Neighbor list.

  • default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).

  • default_r_on (float) – Default turn-on radius \([\mathrm{length}]\).

ZBL specifies that a Ziegler-Biersack-Littmark pair potential should be applied between every non-excluded particle pair in the simulation.

\begin{eqnarray*} V_{\mathrm{ZBL}}(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) \\ + 0.5099 \exp \left( -0.9423 \frac{r_{ij}}{a_F} \right) \\ + 0.2802 \exp \left( -0.4029 \frac{r_{ij}}{a_F} \right) \\ + 0.02817 \exp \left( -0.2016 \frac{r_{ij}}{a_F} \right) \right]; & r < r_{\mathrm{cut}} \\ = & 0; & r > r_{\mathrm{cut}} \\ \end{eqnarray*}

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

See Pair for details on how forces are calculated. Note ZBL does not support energy shifting or smoothing.

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]

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)

Modules