From 0e114b10415ffd682937c91bb91c412ac7c3e9f6 Mon Sep 17 00:00:00 2001 From: Christoph Scherer Date: Thu, 26 May 2022 13:35:40 +0200 Subject: [PATCH] doc/src/pair_sw_3b_table.rst: documentation for added pair style sw/3b/table --- doc/src/pair_sw_3b_table.rst | 309 +++++++++++++++++++++++++++++++++++ 1 file changed, 309 insertions(+) create mode 100644 doc/src/pair_sw_3b_table.rst diff --git a/doc/src/pair_sw_3b_table.rst b/doc/src/pair_sw_3b_table.rst new file mode 100644 index 0000000000..6ee13eb095 --- /dev/null +++ b/doc/src/pair_sw_3b_table.rst @@ -0,0 +1,309 @@ +.. index:: pair_style sw_table + +pair_style sw/3b/table command +=========================== + +Syntax +"""""" + +.. code-block:: LAMMPS + + pair_style style + +* style = *sw/3b/table* + + +Examples +"""""""" + +.. code-block:: LAMMPS + + pair_style sw/3b/table + pair_coeff * * spce.sw type + pair_coeff * * GaN.sw Ga N Ga + + +Description +""""""""""" + +The *sw/3b/table* style is a modification of the original :doc:`pair_style sw `. It has been +developed for coarse-grained simulations (of water) (:ref:`Scherer1 `), +but can be employed for all kinds of systems. It computes a modified 3-body :ref:`Stillinger-Weber ` +potential for the energy E of a system of atoms as + +.. math:: + + E & = \sum_i \sum_{j > i} \phi_2 (r_{ij}) + + \sum_i \sum_{j \neq i} \sum_{k > j} + \phi_3 (r_{ij}, r_{ik}, \theta_{ijk}) \\ + \phi_2(r_{ij}) & = A_{ij} \epsilon_{ij} \left[ B_{ij} (\frac{\sigma_{ij}}{r_{ij}})^{p_{ij}} - + (\frac{\sigma_{ij}}{r_{ij}})^{q_{ij}} \right] + \exp \left( \frac{\sigma_{ij}}{r_{ij} - a_{ij} \sigma_{ij}} \right) \\ + \phi_3(r_{ij},r_{ik},\theta_{ijk}) & = f^{\textrm{3b}}\left(\theta_{ijk}\right) + \exp \left( \frac{\gamma_{ij} \sigma_{ij}}{r_{ij} - a_{ij} \sigma_{ij}} \right) + \exp \left( \frac{\gamma_{ik} \sigma_{ik}}{r_{ik} - a_{ik} \sigma_{ik}} \right) + +where :math:`\phi_2` is a two-body term and :math:`\phi_3` is a +three-body term. The summations in the formula are over all neighbors J +and K of atom I within a cutoff distance :math:`a \sigma`. +In contrast to the original *sw* style, *sw/3b/table* allows for a flexible +three-body term :math:`f^{\textrm{3b}}\left(\theta_{ijk}\right)` which is read in +as a tabulated interaction. It can be parametrized with the csg_fmatch app of VOTCA +as available at: . + +Only a single pair_coeff command is used with the *sw/3b/table* style +which specifies a modified Stillinger-Weber potential file with parameters for all +needed elements. These are mapped to LAMMPS atom types by specifying +N_el additional arguments after the ".sw" filename in the pair_coeff command, +where N_el is the number of LAMMPS atom types: + +* ".sw" filename +* N_el element names = mapping of SW elements to atom types + +See the :doc:`pair_coeff ` page for alternate ways +to specify the path for the potential file. + +As an example, imagine a file SiC.sw has Stillinger-Weber values for +Si and C. If your LAMMPS simulation has 4 atoms types and you want +the first 3 to be Si, and the fourth to be C, you would use the following +pair_coeff command: + +.. code-block:: LAMMPS + + pair_coeff * * SiC.sw Si Si Si C + +The first 2 arguments must be \* \* so as to span all LAMMPS atom types. +The first three Si arguments map LAMMPS atom types 1,2,3 to the Si +element in the SW file. The final C argument maps LAMMPS atom type 4 +to the C element in the SW file. If a mapping value is specified as +NULL, the mapping is not performed. This can be used when a *sw/3b/table* +potential is used as part of the *hybrid* pair style. The NULL values +are placeholders for atom types that will be used with other +potentials. + +The (modified) Stillinger-Weber files have a ".sw" suffix. Lines that are not blank or +comments (starting with #) define parameters for a triplet of +elements. The parameters in a single entry correspond to the two-body +and three-body coefficients in the formula above. Here, also the suffix +".sw" is used though the original Stillinger-Weber file format is supplemented +with four additional lines per parameter block to specify the tabulated +three-body interaction. A single entry then contains: + +* element 1 (the center atom in a 3-body interaction) +* element 2 +* element 3 +* :math:`\epsilon` (energy units) +* :math:`\sigma` (distance units) +* a +* :math:`\lambda` +* :math:`\gamma` +* :math:`\cos\theta_0` +* A +* B +* p +* q +* tol +* filename +* keyword +* style +* N + +The A, B, p, and q parameters are used only for two-body interactions. +The :math:`\lambda` and :math:`\cos\theta_0` parameters, only used +for three-body interactions in the original Stillinger-Weber style, are read +in but ignored in this modified pair style. The :math:`\epsilon` parameter is only used +for two-body interactions in this modified pair style and not for the three-body +terms. The :math:`\sigma` and *a* parameters are used for both two-body and three-body +interactions. :math:`\gamma` is used only in the three-body +interactions, but is defined for pairs of atoms. The non-annotated +parameters are unitless. + +LAMMPS introduces an additional performance-optimization parameter tol +that is used for both two-body and three-body interactions. In the +Stillinger-Weber potential, the interaction energies become negligibly +small at atomic separations substantially less than the theoretical +cutoff distances. LAMMPS therefore defines a virtual cutoff distance +based on a user defined tolerance tol. The use of the virtual cutoff +distance in constructing atom neighbor lists can significantly reduce +the neighbor list sizes and therefore the computational cost. LAMMPS +provides a *tol* value for each of the three-body entries so that they +can be separately controlled. If tol = 0.0, then the standard +Stillinger-Weber cutoff is used. + +The additional parameters *filename*, *keyword*, *style*, and *N* refer to +the tabulated angular potential :math:`f^{\textrm{3b}}\left(\theta_{ijk}\right)`. +The tabulated angular potential has to be of the format as used in the +:doc:`angle_style table ` command: + +An interpolation tables of length *N* is created. The +interpolation is done in one of 2 *styles*: *linear* or *spline*. +For the *linear* style, the angle is used to find 2 surrounding table +values from which an energy or its derivative is computed by linear +interpolation. For the *spline* style, a cubic spline coefficients are computed and +stored at each of the *N* values in the table. The angle is used to +find the appropriate set of coefficients which are used to evaluate a +cubic polynomial which computes the energy or derivative. + +The *filename* specifies the file containing the tabulated energy and +derivative values of :math:`f^{\textrm{3b}}\left(\theta_{ijk}\right)`. +The *keyword* then specifies a section of the file. The +format of this file is as follows (without the +parenthesized comments): + +.. parsed-literal:: + + # Angle potential for harmonic (one or more comment or blank lines) + + HAM (keyword is the first text on line) + N 181 FP 0 0 EQ 90.0 (N, FP, EQ parameters) + (blank line) + 1 0.0 200.5 2.5 (index, angle, energy, derivative) + 2 1.0 198.0 2.5 + ... + 181 180.0 0.0 0.0 + +A section begins with a non-blank line whose first character is not a +"#"; blank lines or lines starting with "#" can be used as comments +between sections. The first line begins with a keyword which +identifies the section. The next line lists (in any +order) one or more parameters for the table. Each parameter is a +keyword followed by one or more numeric values. + +The parameter "N" is required and its value is the number of table +entries that follow. Note that this may be different than the *N* +specified in the Stillinger-Weber potential file. Let +Nsw = *N* in the ".sw" file, and Nfile = "N" in the +tabulated angular file. What LAMMPS does is a preliminary interpolation by +creating splines using the Nfile tabulated values as nodal points. It +uses these to interpolate as needed to generate energy and derivative +values at Ntable different points. The resulting tables of length +Nsw are then used as described above, when computing energy and +force for individual angles and their atoms. This means that if you +want the interpolation tables of length Nsw to match exactly what +is in the tabulated file (with effectively no preliminary +interpolation), you should set Nsw = Nfile. + +The "FP" parameter is optional. If used, it is followed by two values +fplo and fphi, which are the second derivatives at the innermost and +outermost angle settings. These values are needed by the spline +construction routines. If not specified by the "FP" parameter, they +are estimated (less accurately) by the first two and last two +derivative values in the table. + +The "EQ" parameter is also optional. If used, it is followed by a the +equilibrium angle value, which is used, for example, by the :doc:`fix shake ` command. If not used, the equilibrium angle is +set to 180.0. + +Following a blank line, the next N lines of the angular table file list the tabulated values. +On each line, the first value is the index from 1 to N, the second value is +the angle value (in degrees), the third value is the energy (in energy +units), and the fourth is -dE/d(theta) (also in energy units). The third +term is the energy of the 3-atom configuration for the specified +angle. The last term is the derivative of the energy with respect to +the angle (in degrees, not radians). Thus the units of the last term +are still energy, not force. The angle values must increase from one +line to the next. The angle values must also begin with 0.0 and end +with 180.0, i.e. span the full range of possible angles. + +Note that one angular potential file can contain many sections, each with a tabulated +potential. LAMMPS reads the file section by section until it finds +one that matches the specified *keyword* of appropriate section of the ".sw" file. + +The Stillinger-Weber potential file must contain entries for all the +elements listed in the pair_coeff command. It can also contain +entries for additional elements not being used in a particular +simulation; LAMMPS ignores those entries. + +For a single-element simulation, only a single entry is required +(e.g. SiSiSi). For a two-element simulation, the file must contain 8 +entries (for SiSiSi, SiSiC, SiCSi, SiCC, CSiSi, CSiC, CCSi, CCC), that +specify SW parameters for all permutations of the two elements +interacting in three-body configurations. Thus for 3 elements, 27 +entries would be required, etc. + +As annotated above, the first element in the entry is the center atom +in a three-body interaction. Thus an entry for SiCC means a Si atom +with 2 C atoms as neighbors. The parameter values used for the +two-body interaction come from the entry where the second and third +elements are the same. Thus the two-body parameters for Si +interacting with C, comes from the SiCC entry. The three-body +angular potential :math:`f^{\textrm{3b}}\left(\theta_{ijk}\right)` +can in principle be specific to the three elements of the +configuration. However, the user must ensure that it makes physically sense. +Note also that the function :math:`\phi_3` contains two exponential +screening factors with parameter values from the ij pair and ik +pairs. So :math:`\phi_3` for a C atom bonded to a Si atom and a second C atom +will depend on the three-body parameters for the CSiC entry, and also +on the two-body parameters for the CCC and CSiSi entries. Since the +order of the two neighbors is arbitrary, the three-body parameters +and the tabulated angular potential for +entries CSiC and CCSi should be the same. Similarly, the two-body +parameters for entries SiCC and CSiSi should also be the same. The +parameters used only for two-body interactions (A, B, p, and q) in +entries whose second and third element are different (e.g. SiCSi) are not +used for anything and can be set to 0.0 if desired. +This is also true for the parameters in :math:`\phi_3` that are +taken from the ij and ik pairs (:math:`\sigma`, *a*, :math:`\gamma`) + +---------- + +.. include:: accel_styles.rst + +.. note:: + + When using the INTEL package with this style, there is an additional + 5 to 10 percent performance improvement when the Stillinger-Weber + parameters p and q are set to 4 and 0 respectively. These + parameters are common for modeling silicon and water. + +---------- + +Mixing, shift, table, tail correction, restart, rRESPA info +""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" + +For atom type pairs I,J and I != J, where types I and J correspond to +two different element types, mixing is performed by LAMMPS as +described above from values in the potential file, but not for the +tabulated angular potential file. + +This pair style does not support the :doc:`pair_modify ` +shift, table, and tail options. + +This pair style does not write its information to :doc:`binary restart files `, since it is stored in potential files. +Thus, you need to re-specify the pair_style and pair_coeff commands in an input +script that reads a restart file. + +This pair style can only be used via the *pair* keyword of the +:doc:`run_style respa ` command. It does not support the +*inner*, *middle*, *outer* keywords. + +---------- + +Restrictions +"""""""""""" + +This is a user pair style. For more information, see :ref:`Scherer1 `. It is only enabled +if LAMMPS was explicitly built with it. + +This pair style requires the :doc:`newton ` setting to be "on" +for pair interactions. + +For an example of an extended Stillinger-Weber potential file, have a look at the tutorial +in the tutorial folder. + +Related commands +"""""""""""""""" + +:doc:`pair_coeff ` + + +---------- + +.. _Stillinger2: + +**(Stillinger)** Stillinger and Weber, Phys Rev B, 31, 5262 (1985). + +.. _Scherer1: + +**(Scherer1)** C. Scherer and D. Andrienko, Phys. Chem. Chem. Phys. 20, 22387-22394 (2018). +