md.many_body¶
Overview
Base class triplet force. |
|
Reversible crosslinker three-body force. |
|
Soft force for simulating a van der Waals liquid. |
|
Tersoff force. |
Details
Implement many body potentials.
Triplet 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)\).
Triplet
applies a short range cutoff for performance and assumes that both
\(U(\vec{r}_{ij}, \vec{r}_{ik})\) and its derivatives are 0 when
\(r_{ij} > r_\mathrm{cut}\) or \(r_{ik} > r_\mathrm{cut}\).
Specifically, the force \(\vec{F}\) applied to each particle \(i\) is:
The per particle energy terms are:
- class hoomd.md.many_body.RevCross(nlist, default_r_cut=None)¶
Bases:
Triplet
Reversible crosslinker three-body force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
RevCross
computes the revcross three-body force on every particle in the simulation state. Despite the fact that the revcross potential accounts for the effects of third bodies, it is actually just a combination of two body potential terms. It can thus use type-pair parameters similar to those of the pair potentials.The RevCross potential has been described in detail in S. Ciarella and W.G. Ellenbroek 2019. It is based on a generalized Lennard-Jones pairwise attraction to form bonds between interacting particles:
\[\begin{split}U_{ij}(r) = \begin{cases} 4 \varepsilon \left[ \left( \dfrac{ \sigma}{r_{ij}} \right)^{2n} - \left( \dfrac{ \sigma}{r_{ij}} \right)^{n} \right] & r < r_\mathrm{cut} \\ 0 & r \ge r_\mathrm{cut} \end{cases}\end{split}\]Then an additional three-body repulsion is evaluated to compensate the bond energies imposing single bond per particle condition:
\[v^{\left( 3b \right)}_{ijk} = \lambda_3 \epsilon \hat{v}^{ \left( 2b \right)}_{ij} \left(\vec{r}_{ij}\right) \cdot \hat{v}^{ \left( 2b \right)}_{ik} \left(\vec{r}_{ik}\right)~,\]where the two body potential is rewritten as:
\[\begin{split}\hat{v}^{ \left( 2b \right)}_{ij}\left(\vec{r}_{ij}\right) = \begin{cases} 1 & r \le r_{min} \\ - \dfrac{v_{ij}\left(\vec{r}_{ij}\right)}{\epsilon} & r > r_{min} \\ \end{cases}\end{split}\]Attention
The RevCross potential models an asymmetric interaction between two different chemical moieties that can form a reversible bond. This requires the definition of (at least) two different types of particles. A reversible bond is only possible between two different species, otherwise \(v^{\left( 3b \right)}_{ijk}\), would prevent any bond. In our example we then set the interactions for types A and B with
rev_c.params[[('A','B'),('A','B')]]
to{"sigma": 0.0, "n": 0, "epsilon": 0, "lambda3": 0}
and the only non-zero energy only between the different types with settingrev_c.params[('A','B')]
to{"sigma":1, "n": 100, "epsilon": 100, "lambda3": 1}
. Notice that the number of the minority species corresponds to the maximum number of bonds.This three-body term also tunes the energy required for a bond swap through the unitless coefficient \(\lambda_3\) . In S. Ciarella and W.G. Ellenbroek 2019 is explained that setting \(\lambda_3=1\) corresponds to no energy requirement to initiate bond swap, while this energy barrier scales roughly as \(\beta \Delta E_\text{sw} =\beta \varepsilon(\lambda_3-1)\).
Note
Choosing \(\lambda<1\) pushes the system to cluster because the three-body term is not enough to compensate the energy of multiple bonds, so it may cause nonphysical situations.
- params¶
The revcross 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}]\)lambda3
(float
, required) - \(\lambda_3\) \([\mathrm{dimensionless}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
Example:
nl = md.nlist.Cell() bond_swap = md.many_body.RevCross(default_r_cut=1.3,nlist=nl) bond_swap.params[('A', 'A'), ('B', 'B')] = { "sigma":0,"n": 0, "epsilon": 0, "lambda3": 0} # a bond can be made only between A-B and not A-A or B-B bond_swap.params[('A','B')] = { "sigma": 1, "n": 100, "epsilon": 10, "lambda3": 1}
- class hoomd.md.many_body.SquareDensity(nlist, default_r_cut=None)¶
Bases:
Triplet
Soft force for simulating a van der Waals liquid.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
SquareDensity
that the square density three-body force should on every particle in the simulation state.The self energy per particle takes the form:
\[\Psi^{ex} = B (\rho - A)^2\]which gives a pair-wise additive, three-body force:
\[\vec{f}_{ij} = \left( B (n_i - A) + B (n_j - A) \right) w'_{ij} \vec{e}_{ij}\]Here, \(w_{ij}\) is a quadratic, normalized weighting function,
\[w(x) = \frac{15}{2 \pi r_{c,\mathrm{weight}}^3} (1-r/r_{c,\mathrm{weight}})^2\]The local density at the location of particle i is defined as
\[n_i = \sum\limits_{j\neq i} w_{ij} \left(\big| \vec r_i - \vec r_j \big|\right)\]- params¶
The SquareDensity potential parameters. The dictionary has the following keys:
A
(float
, required) - \(A\) - mean density (default:0) \([\mathrm{length}^{-2}]\) in 2D and \([\mathrm{length}^{-3}]\) in 3DB
(float
, required) - \(B\) - coefficient of the harmonic density term \([\mathrm{energy} \cdot \mathrm{length}^4]\) in 2D and \([\mathrm{energy} \cdot \mathrm{length}^6]\) in 3D
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
Example:
nl = nlist.Cell() sqd = md.many_body.SquareDensity(nl, default_r_cut=3.0) sqd.params[('A', 'B')] = dict(A=1.0, B=2.0) sqd.params[('B', 'B')] = dict(A=2.0, B=2.0, default_r_on=1.0)
For further details regarding this multibody potential, see
[1] P. B. Warren, “Vapor-liquid coexistence in many-body dissipative particle dynamics” Phys. Rev. E. Stat. Nonlin. Soft Matter Phys., vol. 68, no. 6 Pt 2, p. 066702, 2003.
- class hoomd.md.many_body.Tersoff(nlist, default_r_cut=None)¶
Bases:
Triplet
Tersoff force.
- Parameters:
nlist (hoomd.md.nlist.NeighborList) – Neighbor list
default_r_cut (float) – Default cutoff radius \([\mathrm{length}]\).
The Tersoff potential is a bond-order potential based on the Morse potential that accounts for the weakening of individual bonds with increasing coordination number. It does this by computing a modifier to the attractive term of the potential. The modifier contains the effects of third-bodies on the bond energies. The potential also includes a smoothing function around the cutoff. The implemented smoothing function is exponential in nature as opposed to the sinusoid used by J. Tersoff 1988.
Tersoff
computes the Tersoff three-body force on every particle in the simulation state. Despite the fact that the Tersoff potential accounts for the effects of third bodies, the species of the third body is irrelevant. It can thus use type-pair parameters similar to those of the pair potentials:\[U_{ij}(r) = \frac{1}{2} f_C(r_{ij}) \left[f_R(r_{ij}) + b_{ij}f_A(r_{ij})\right]\]where
\[ \begin{align}\begin{aligned}f_R(r) = A_1e^{\lambda_1(r_D-r)}\\f_A(r) = A_2e^{\lambda_2(r_D-r)}\end{aligned}\end{align} \]\[\begin{split}f_C(r) = \begin{cases} 1 & r < r_{\mathrm{cut}} - r_{CT} \\ \exp \left[-\alpha\frac{x(r)^3}{x(r)^3 - 1} \right] & r_{\mathrm{cut}} - r_{CT} < r < r_{\mathrm{cut}} \\ 0 & r > r_{\mathrm{cut}} \end{cases}\end{split}\]\[b_{ij} = (1 + \gamma^n\chi_{ij}^n)^{\frac{-1}{2n}}\]In the definition of \(f_C(r)\), there is a quantity \(x(r)\), which is defined as
\[x(r) = \frac{r - (r_{\mathrm{cut}} - r_{CT})}{r_{CT}}\]which ensures continuity between the different regions of the potential. In the definition of \(b_{ij}\), there is a quantity \(\chi_{ij}\) which is defined as
\[ \begin{align}\begin{aligned}\chi_{ij} = \sum_{k \neq i,j} f_C(r_{ik}) \cdot e^{\lambda_3^3 |r_{ij} - r_{ik}|^3} \cdot g(\theta_{ijk})\\g(\theta_{ijk}) = 1 + \frac{c^2}{d^2} - \frac{c^2}{d^2 + |m - \cos(\theta_{ijk})|^2}\end{aligned}\end{align} \]- params¶
The Tersoff potential parameters. The dictionary has the following keys:
magnitudes
(tuple[float
,float
]) - \((A_1, A_2)\) - Magnitudes of the repulsive and attractive terms (default: (1.0, 1.0)) \([\mathrm{energy}]\)exp_factors
(tuple[float
,float
]) - \((\lambda_1, \lambda_2)\) - exponential factors of the repulsive and attractive terms (default: 2.0) \([\mathrm{length}^{-1}]\)lambda3
(float
) - \(\lambda_3\) - exponential factor in \(\chi_{ij}\) (default: 0.0) \([\mathrm{length}^{-1}]\)dimer_r
(float
) - \(r_D\) - length shift in attractive and repulsive terms (default: 1.5) \([\mathrm{length}]\)cutoff_thickness
(float
) - \(r_{CT}\) - distance which defines the different regions of the potential (default: 0.2) \([\mathrm{length}]\)alpha
(float
) - \(\alpha\) - decay rate of the cutoff term \(f_C(r)\) (default: 3.0) \([\mathrm{dimensionless}]\)n
(float
) - \(n\) - power in \(b_{ij}\) (default: 0.0) \([\mathrm{dimensionless}]\)gamma
(float
) - \(\gamma\) - coefficient in \(b_{ij}\) (default: 0.0) \([\mathrm{dimensionless}]\)c
(float
) - \(c\) - coefficient in \(g(\theta)\) (default: 0.0) \([\mathrm{dimensionless}]\)d
(float
) - \(d\) - coefficient in \(g(\theta)\) (default: 1.0) \([\mathrm{dimensionless}]\)m
(float
) - \(m\) - coefficient in \(g(\theta)\) (default: 0.0) \([\mathrm{dimensionless}]\)
Type:
TypeParameter
[tuple
[particle_type
,particle_type
],dict
]
Example:
nl = md.nlist.Cell() tersoff = md.many_body.Tersoff(default_r_cut=1.3, nlist=nl) tersoff.params[('A', 'B')] = dict(magnitudes=(2.0, 1.0), lambda3=5.0)
- class hoomd.md.many_body.Triplet(nlist, default_r_cut=None)¶
Bases:
Force
Base class triplet force.
Triplet
is the base class for many-body triplet forces.Warning
This class should not be instantiated by users. The class can be used for
isinstance
orissubclass
checks.- 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
])
- nlist¶
Neighbor list used to compute the triplet potential.
Warning
Currently HOOMD-blue does not support reverse force communication between MPI domains on the GPU. Since reverse force communication is required for the calculation of three-body forces, attempting to use this potential on the GPU with MPI will result in an error.