1CPN Implementation in Lammps¶

Unfortunately, this page is still a bit rough. Its goal is to give an idea of the custom potentials specific to 1CPN that we implemented into LAMMPS and what the different parameters mean.

In this regard, this page is currently okay, but I’ll conceed that it could use a bit of work.

Linking 1CPN with LAMMPS¶

To link 1CPN to LAMMPS please check out Quick Start.

1CPN Performance across many processors¶

The key to good 1CPN performance across multiple processors is the LAMMPS fix balance command.

Custom 1CPN Potentials¶

While developing the 1CPN model, many custom potentials were necessary

Pair Zewdie¶

Give functional form of potential? And an explanation of the pair_coeffs arguments

Usage:

pair_style zewdie ${pe000} ${pecc2} ${pe220} ${pe222} ${pe224} ${ps000} ${pscc2} ${ps220} ${ps222} ${ps224}

pair_coeff   1 1 zewdie ${pe0} ${ps0} ${cutoff}

pe000 - See $\epsilon_{000}$ in (3)
pecc2 - See $\epsilon_{cc2}$ in (3)
pe220 - See $\epsilon_{220}$ in (3)
pe222 - See $\epsilon_{222}$ in (3)
pe224 - See $\epsilon_{224}$ in (3)
ps000 - See $\sigma_{000}$ in (2)
pscc2 - See $\sigma_{cc2}$ in (2)
ps220 - See $\sigma_{220}$ in (2)
ps222 - See $\sigma_{222}$ in (2)
ps224 - See $\sigma_{224}$ in (2)
pe0 - See $\epsilon_{0}$ in (3)
ps0 - See $\sigma_{0}$ in (2)

(1)¶\[{U}_{Zewdie} ({r}_{ij},{\hat{f}}_i , {\hat{f}}_j \,;\, \sigma_0, \epsilon_0) = 4 \epsilon \left[ \left(\frac{\sigma_0}{r_{ij}-\sigma+\sigma_0}\right)^{12}- \left(\frac{\sigma_0}{r_{ij}-\sigma+\sigma_0)}\right)^6 \right]\]

where

(2)¶\[\sigma = \sigma_{0} [\sigma_{000}S_{000} + \sigma_{cc2}(S_{022} + S_{202}) + \sigma_{220}S_{220} + \sigma_{222}S_{222} + \sigma_{224}S_{224}]\]

(3)¶\[\epsilon = \epsilon_{0} [\epsilon_{000}S_{000} + \epsilon_{cc2}(S_{022} + S_{202}) + \epsilon_{220}S_{220} + \epsilon_{222}S_{222} + \epsilon_{224}S_{224}].\]

and the S-functions are defined as

\[\begin{split}S_{000} &= 1,\\ S_{202} &= (3 a_1^2 - 1) / 2\sqrt{5},\\ S_{022} &= (3 a_2^2 - 1) / 2\sqrt{5},\\ S_{220} &= (3 a_0^2 - 1) / 2\sqrt{5},\\ S_{222} &= \dfrac{1}{\sqrt{70}}(2 - 3 a_1^2 - 3 a_2^2 - 3 a_0^2 + 9 a_0 a_1 a_2),\\ S_{224} &= \dfrac{1}{4\sqrt{70}}(1 + 2 a_0^2 - 5 a_1^2 - 5 a_2^2 - 20 a_0 a_1 a_2 + 35 a_1^2 a_2^2)\\\end{split}\]

with

\[\begin{split}a_0 &= {\hat{f}}_i \cdot {\hat{f}}_j, \\ a_1 &= {\hat{f}}_i \cdot {\hat r}_{ij}, \\ a_2 &= {\hat{f}}_j \cdot {\hat{r}}_{ij} \\\end{split}\]

Discuss how sphere-ellipse interactions are handled

Pair Gauss Aniso¶

Usage:

pair_style gauss/aniso ${gauss_rcut}

pair_coeff 2 3  sigma d0 r0 theta0 phi0 Ktheta Kphi

sigma - Gaussian width. Given by $\sigma$ in (4)
d0 - Gaussian depth. Given by $d_0$ in (4)
r0 - Gaussian Center Position. Given by $r_0$ in (4)
theta0 - Equilibrium position with with respect to $\theta$, where $\theta = \arccos(\hat{r}_{ij} \cdot \hat{u}_i)$. See $\theta_0$ in (5)
phi0 - Equilibrium position with with respect to $\phi$, where $\phi = \arccos(\hat{r}_{ij} \cdot \hat{f}_i)$. See $\phi_0$ in (5)
Ktheta - Width of modulating function with respect to $\theta$. See $K_\theta$ in (5)
Kphi - Width of modulating function with respect to $\phi$. See $K_\phi$ in (5)

Warning

In the implementation, of this potential, the atom_types of sites i and j cannot be the same.

This is because the ith particle is always chosen to be the atom with the lower type index.

Not symmetric.

For example

(4)¶\[\begin{split}U_{gauss,aniso} &= f(K_\theta,\Delta \theta) f(K_\phi, \Delta \phi) \, {U}_{gauss} \nonumber \\ &= f(K_\theta,\Delta \theta) f(K_\phi, \Delta \phi) \left(-d_0 e^{-(r-r_0)^2 / 2\sigma^2}\right)\end{split}\]

(5)¶\[\begin{split}f(K_\theta,\Delta \theta) = \begin{cases} 1 &-\frac{\pi}{2K_\theta} < \Delta \theta < \frac{\pi}{2K_\theta} \\ 1-\cos^2\left( K_\theta \Delta \theta \right) & \frac{-\pi}{K_\theta}<\Delta \theta < \frac{-\pi}{2K_\theta} \textrm{ or } \frac{\pi}{2K_\theta} < \Delta \theta < \frac{\pi}{K_\theta}\\ 0 &\Delta \theta < -\frac{\pi}{K_\theta} \textrm{ or } \Delta \theta > \frac{\pi}{K_\theta} \end{cases}\end{split}\]

Angle Orient¶

Usage::

angle_style orient
angle_coeff 1 angle_{f,v,u} ktheta1 ktheta2 kphi theta1 theta2 phi

angle_vector - Possible values angle_f, angle_v, or angle_u. Defines whether $\hat{w} = \{\hat{f},\hat{v},\hat{u}\}$.
ktheta1 - Spring constant for deformations in $\theta_1$. See (6).
ktheta2 - Spring constant for deformations in $\theta_2$. See (6).
kphi - Spring constant for deformations in $\phi$. See (6).
theta1 - Equilibrium value of $\theta_1$. See $\theta_{1,0}$ in (6).
theta2 - Equilibrium value of $\theta_2$. See $\theta_{2,0}$ in (6).
phi - Equilibrium value of $\phi$. See $\phi_{0}$ in (6).

(6)¶\[U = \frac{1}{2} \left( k_{\theta_1} (\theta_1 - \theta_{1,0})^2 + k_{\theta_2} (\theta_2 - \theta_{2,0})^2 + k_{\phi} (\phi - \phi_{0})^2 \right)\]

where $\theta_1, \theta_2, \phi$ are given by

\[\begin{split}\theta_1 &= \arccos (\hat{w}_i \cdot \hat{r}_{ij}) \\ \theta_2 &= \arccos (\hat{w}_j \cdot \hat{r}_{ij}) \\ \phi &= \arccos (\hat{w}_i \cdot \hat{w}_j)\end{split}\]

Angle Orient Cosine¶

Usage::

angle_style orient/cosine
angle_coeff 1 angle_{f,v,u} ktheta1 ktheta2 kphi theta1 theta2 phi

angle_vector - Possible values angle_f, angle_v, or angle_u. Defines whether $\hat{w} = \{\hat{f},\hat{v},\hat{u}\}$.
ktheta1 - Spring constant for deformations in $\theta_1$. See (7).
ktheta2 - Spring constant for deformations in $\theta_2$. See (7).
kphi - Spring constant for deformations in $\phi$. See (7).
theta1 - Equilibrium value of $\theta_1$. See $\theta_{1,0}$ in (7).
theta2 - Equilibrium value of $\theta_2$. See $\theta_{2,0}$ in (7).
phi - Equilibrium value of $\phi$. See $\phi_{0}$ in (7).

(7)¶\[U = \frac{1}{2} \left[ k_{\theta_1} \left(1 - \cos(\theta_1 - \theta_{1,0})\right) + k_{\theta_2} \left(1 - \cos(\theta_2 - \theta_{2,0})\right) + k_{\phi} \left(1 - \cos(\phi - \phi_{0})\right) \right]\]

where $\theta_1, \theta_2, \phi$ are given by

\[\begin{split}\theta_1 &= \arccos (\hat{w}_i \cdot \hat{r}_{ij}) \\ \theta_2 &= \arccos (\hat{w}_j \cdot \hat{r}_{ij}) \\ \phi &= \arccos (\hat{w}_i \cdot \hat{w}_j)\end{split}\]

Angle WLC Twist¶

Usage::

angle_style wlctwist
angle_coeff 2 wlctwist ${kalign} ${ktwist} ${omega0}

kalign - Alignment spring constant
ktwist - Twist spring constant
omega0 - Equilibrium Twist

This potential is based off of the implementation of Brackley et al.

(8)¶\[U = k_{\omega} \left(1-\cos(\omega_i - \omega_{0})\right) + k_{\psi} \left(1-\cos (\psi_i)\right)\]

Maintaing Compatability with LAMMPS¶

In order to maintain compatability of 1CPN with the most recent version of LAMMPS it is helpful to know which core LAMMPS potential the 1CPN potentials were derived from. By seeing what changed between the core LAMMPS potentials (i.e. diff old and new version), it is typically straightforward to make the necessary minor changes to the 1CPN potential to allow 1CPN-LAMMPS to compile.

When trying a new version of LAMMPS, be sure to run the integration tests in ${D_1CPN}/test/integ_tests, to make sure the model is behaving correctly

1CPN Potential Original Lammps Potential

pair_style zewdie pair_style gayberne

pair_style gauss/aniso pair_style gayberne, pair_style gauss

angle_style wlctwist angle_style wlctwist (Brackley2014)

angle_style orient angle_style wlctwist

angle_style orient/cosine angle_style wlctwist