tx258/lammps-gran-kokkos
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

doc pages for AMOEBA/HIPPO

Browse Source
This commit is contained in:
Steve Plimpton
2022-04-18 16:38:41 -06:00
parent fccca3405a
commit 8932d9ffaa
18 changed files with 1015 additions and 268 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
doc/src/Commands_bond.rst
View File

@ -71,6 +71,7 @@ OPT.
*
*
*
* :doc:`amoeba <angle_amoeba>`
* :doc:`charmm (iko) <angle_charmm>`
* :doc:`class2 (ko) <angle_class2>`
* :doc:`class2/p6 <angle_class2>`
@ -149,6 +150,7 @@ OPT.
*
*
*
* :doc:`amoeba <improper_amoeba>`
* :doc:`class2 (ko) <improper_class2>`
* :doc:`cossq (o) <improper_cossq>`
* :doc:`cvff (io) <improper_cvff>`

2
doc/src/Commands_fix.rst
View File

@ -27,6 +27,8 @@ OPT.
* :doc:`adapt/fep <fix_adapt_fep>`
* :doc:`addforce <fix_addforce>`
* :doc:`addtorque <fix_addtorque>`
* :doc:`amoeba/bitorsion <fix_amoeba_bitorsion>`
* :doc:`amoeba/pitorsion <fix_amoeba_pitorsion>`
* :doc:`append/atoms <fix_append_atoms>`
* :doc:`atc <fix_atc>`
* :doc:`atom/swap <fix_atom_swap>`

2
doc/src/Commands_pair.rst
View File

@ -38,6 +38,7 @@ OPT.
* :doc:`agni (o) <pair_agni>`
* :doc:`airebo (io) <pair_airebo>`
* :doc:`airebo/morse (io) <pair_airebo>`
* :doc:`amoeba <pair_amoeba>`
* :doc:`atm <pair_atm>`
* :doc:`awpmd/cut <pair_awpmd>`
* :doc:`beck (go) <pair_beck>`
@ -122,6 +123,7 @@ OPT.
* :doc:`hbond/dreiding/lj (o) <pair_hbond_dreiding>`
* :doc:`hbond/dreiding/morse (o) <pair_hbond_dreiding>`
* :doc:`hdnnp <pair_hdnnp>`
* :doc:`hippo <pair_amoeba>`
* :doc:`ilp/graphene/hbn <pair_ilp_graphene_hbn>`
* :doc:`kolmogorov/crespi/full <pair_kolmogorov_crespi_full>`
* :doc:`kolmogorov/crespi/z <pair_kolmogorov_crespi_z>`

307
doc/src/Howto_amoeba.rst
View File

@ -2,36 +2,87 @@ AMOEBA and HIPPO force fields
=============================
The AMOEBA and HIPPO polarizable force fields were developed by Jay
Ponder's group at U Washington at St Louis. Their implementation in
LAMMPS was done using F90 code provided by the Ponder group from their
`Tinker MD code <https://dasher.wustl.edu/tinker/>`_
Ponder's group at the U Washington at St Louis. Their implementation
in LAMMPS was done using F90 code provided by the Ponder group from
their `Tinker MD code <https://dasher.wustl.edu/tinker/>`_.
NOTE: what version of AMOEBA and HIPPO does LAMMPS implement?
NOTE: What version of AMOEBA and HIPPO does LAMMPS implement?
These force fields can be used when polarization effects are desired
in simulations of water, organic molecules, and biomolecules including
proteins, provided that parameterizations (force field files) are
available for the systems you are interested in. Files in the LAMMPS
potentials directory with a "amoeba" or "hippo" suffix can be used.
The Tinker distribution and website may have other force field files.
proteins, provided that parameterizations (Tinker PRM force field
files) are available for the systems you are interested in. Files in
the LAMMPS potentials directory with a "amoeba" or "hippo" suffix can
be used. The Tinker distribution and website have additional force
field files as well.
NOTE: Can we include a few Tinker FF files in the LAMMPS distro,
e.g. ubiquitin.prm as potentials/ubiquitin.prm.amoeba ?
Note that currently, HIPPO can only be used for water systems, but
HIPPO files for a variety of small organic and biomolecules are in
preparation by the Ponder group. Those force field files will be
included in the LAMMPS distribution when available.
The :doc:`pair_style amoeba <pair_amoeba>` doc page gives a brief
description of the AMOEBA and HIPPO force fields. Further details for
AMOEBA are in these papers: :ref:`(Ren) <amoeba-Ren>`, :ref:`(Shi)
<amoeba-Shi>`. Further details for HIPPO are in this paper:
:ref:`(Rackers) <amoeba-Rackers>`.
----------
The AMOEBA and HIPPO force fields contain the following terms in their
energy (U) computation. Further details for AMOEBA are in these
papers: :ref:`(Ren) <amoeba-Ren>`, :ref:`(Shi) <amoeba-Shi>`. Further
details for HIPPO are in this paper: :ref:`(Rackers)
<amoeba-Rackers>`.
.. math::
U & = U_{intermolecular} + U_{intramolecular) \\
U_{intermolecular} & = U_{hal} + U_{repulsion} + U_{dispersion} + U_{multipole} + U_{polar} + U_{qxfer} \\
U_{intramolecular} & = U_{bond} + U_{angle} + U_{torsion} + U_{oop} + U_{b\theta} + U_{Urey-Bradley} + U_{pitorsion} + U_{bitorsion}
For intermolecular terms, the AMOEBA force field includes only the
:math:`U_{hal}`, :math:`U_{multipole}`, :math:`U_{polar}` terms. The
HIPPO force field includes all but the :math:`U_{hal}` term. In
LAMMPS, these are all computed by the :doc:`pair_style <pair_style>`
command.
For intramolecular terms, the :math:`U_{bond}`, :math:`U_{angle}`,
:math:`U_{torsion}`, :math:`U_{oop}` terms are computed by the
:doc:`bond_style class2 <bond_class2>` :doc:`angle_style amoeba
<angle_amoeba>`, :doc:`dihedral_style fourier <dihedral_fourier>`, and
:doc:`improper_style amoeba <improper_amoeba>` commands respectively.
The :doc:`angle_style amoeba <angle_amoeba>` command includes the
:math:`U_{b\theta}` bond-angle cross term, and the
:math:`U_{Urey-Bradley}` term for a bond contribution between the 1-3
atoms in the angle.
The :math:`U_{pitorsion}` term is computed by the :doc:`fix
amoeba/pitorsion <fix_pitorsion>` command. It computes 6-body
interaction between a pair of bonded atoms which each have 2
additional bond partners.
The :math:`U_{bitorsion}` term is computed by the :doc:`fix
amoeba/bitorsion <fix_bitorsion>` command. It computes 5-body
interaction between two 4-body torsions (dihedrals) which overlap,
having 3 atoms in common.
These command doc pages have additional details on the terms they
compute:
* :doc:`pair_style amoeba or hippo <pair_amoeba>`
* :doc:`bond_style class2 <bond_class2>`
* :doc:`angle_style amoeba <angle_amoeba>`
* :doc:`dihedral_style fourier <dihedral_fourier>`
* :doc:`improper_style amoeba <improper_amoeba>`
* :doc:`fix amoeba/pitorsion <fix_pitorsion>`
* :doc:`fix amoeba/bitorsion <fix_bitorsion>`
----------
To use the AMOEBA force field in LAMMPS you should use command like
these appropriately in your input script. The only change needed for
a HIPPO simulation is for the pair_style and pair_coeff commands. See
examples/amoeba for example input scripts for both AMOEBA and HIPPO.
To use the AMOEBA or HIPPO force fields in LAMMPS, use commands like
the following appropriately in your input script. The only change
needed for AMOEBA vs HIPPO simulation is for the :doc:`pair_style
<pair_style>` and :doc:`pair_coeff <pair_coeff>` commands, as shown
below. See examples/amoeba for example input scripts for both AMOEBA
and HIPPO.
.. code-block:: LAMMPS
@ -66,73 +117,183 @@ examples/amoeba for example input scripts for both AMOEBA and HIPPO.
special_bonds lj/coul 0.5 0.5 0.5 one/five yes # 1-5 neighbors
The data file read by the :doc:`read_data <read_data>` command should
be created by the tools/tinker/tinker2lmp.py conversion program
described below. It will create a section in the data file with the
header "Tinker Types". A :doc:`fix property/atom <fix_property_atom>`
command for the data must be specified befroe the read_data command.
In the example above the fix ID is amtype.
NOTE: that some options not needed for simpler systems
Similarly, if the system you are simulating defines AMOEBA/HIPPO
pitorsion or bitorsion interactions, there will be entries in the data
file for those interactions. They require a :doc:`fix
amoeba/pitortion <fix_amoeba_pitortion>` and :doc:`fix
amoeba/bitorsion <fix_amoeba_bitorsion>` command be defined. In the
example above, the IDs for these two fixes are "pit" and "bit".
NOTE: tinker2lmp.py tool
Of course, if the system being modeled does not have one or more of
the following -- bond, angle, dihedral, improper, pitorision,
bitorsion interactions -- then the corresponding style and fix
commands above do not need to be used. See the example scripts in
examples/amoeba for water systems as examples; they are simpler than
what is listed above.
NOTE: are some commands not needed if system is simple and not ubi2 ?
The two :doc:`fix property/atom <fix_property_atom>` commands with IDs
"extra" and "polaxe" are also needed to define internal per-atom
quantities used by the AMOEBA and HIPPO force fields.
Both AMOEBA and HIPPO use amoeba for bond/angle/dihedral styles,
assuming the molecular system has bonds, angles, or dihedrals. These
style commands must be be used before the data file is read, as the
data file defines the coefficients for the various bond/angle/dihedral
types. The pair_style command can come after the read_data command,
as no pair coefficients are defined in the data file.
The :doc:`pair_coeff <pair_coeff>` command used for either the AMOEBA
or HIPPO force field takes two arguments for Tinker force field files,
namely a PRM and KEY file. The keyfile can be specified as NULL and
default values for a various settings will be used. Note that these 2
files are meant to allow use of native Tinker files as-is. However
LAMMPS does not support all the options which can be included
in a Tinker PRM or KEY file. See specifis below.
The :doc:`fix property/atom <fix_property_atom>` commands are required
to create per-atom data used by the AMOEBA/HIPPO force fields, some of
which is defined in the LAMMPS data file. The LAMMPS data file should
be produced as a pre-processing step by the tinker2lmp.py Python
script in tools/amoeba with a command like one of the following.
The tools/amoeba/README file has more details on using this tool.
A :doc:`special_bonds <special_bonds>` command with the "one/five"
option is required, since the AMOEBA/HIPPO force fields define
weighting factors for not only 1-2, 1-3, 1-4 interactions, but also
1-5 interactions. This command will trigger a per-atom list of 1-5
neighbors to be generated. The AMOEBA and HIPPO force fields define
their own custom weighting factors for all the 1-2, 1-3, 1-4, 1-5
terms which in the Tinker PRM and KEY files; they can be different for
different terms in the force field.
.. code-block:: bash
In addition to the list above, these command doc pages have additional
details:
% python tinker2lmp.py -xyz tinker.xyz -amoeba protein.prm.amoeba -data lmp.data # AMOEBA for non-periodic systems
% python tinker2lmp.py -xyz tinker.xyz -amoeba protein.prm.amoeba -data lmp.data -pbc 10.0 20.0 15.0 # AMOEBA for periodic systems
% python tinker2lmp.py -xyz tinker.xyz -hippo water.prm.hippo -data lmp.data # HIPPO for non-periodic systems
% python tinker2lmp.py -xyz tinker.xyz -hippo water.prm.hippo -data lmp.data -pbc 10.0 20.0 15.0 # HIPPO for periodic systems
Note that two input files are needed and a LAMMPS data file (lmp.data)
is produced. The data file will have information on Tinker atom types
and AMOEBA/HIPPO force field parameters for bonds, angles, and
dihedrals.
The first input is an XYZ file listing all the atoms in the
system.
NOTE: is this a Tinker-augmented-XYZ format or standard? In either
case, how do we suggest LAMMPS users come up with these files?
The second input is an AMOEBA or HIPPO PRM (force field) file. The
format of these files is defined by Tinker. A few such files are
provided in the LAMMPS potentials directory. Others may be available
in the Tinker distribution or from the Ponder group.
The pair_coeff command should specify the same PRM file, and
optionally a Tinker-format KEY file. See the :doc:`pair_style amoeba
<pair_amoeba>` doc page for more information about Tinker PRM and KEY
files.
Finally, the :doc:`special_bonds <special_bonds>` command is used to
set all LJ and Coulombic 1-2, 1-3, 1-4 weighting factors to non-zero
and non-unity values, and to generate a per-atom list of 1-5 neighbors
as well. This is to insure all bond-topology neighbors are included
in the neighbor lists used by AMOEBA/HIPPO. These force fields apply
their own custom weighting factors to all these terms, including the
1-5 neighbors.
* :doc:`atom_style amoeba <atom_style>`
* :doc:`fix property/atom <fix_property_atom>`
* :doc:`special_bonds <special_bonds>`
----------
These command doc pages have additional details:
Tinker PRM and KEY files
* :doc:`pair_style amoeba or hippo <pair_ameoba>`
* :doc:`bond_style amoeba <bond_amoeba>`
* :doc:`angle_style amoeba <angle_charmm>`
* :doc:`dihedral_style amoeba <dihedral_amoeba>`
* :doc:`fix property/atom <fix_property_atom>`
* :doc:`special_bonds <special_bonds>`
A Tinker PRM file is composed of sections, each of which has multiple
lines. This is the list of sections LAMMPS knows how to parse and
use. Any other sections are skipped:
* Angle Bending Parameters
* Atom Type Definitions
* Atomic Multipole Parameters
* Bond Stretching Parameters
* Charge Penetration Parameters
* Charge Transfer Parameters
* Dipole Polarizability Parameters
* Dispersion Parameters
* Force Field Definition
* Literature References
* Out-of-Plane Bend Parameters
* Pauli Repulsion Parameters
* Pi-Torsion Parameters
* Stretch-Bend Parameters
* Torsion-Torsion Parameters
* Torsional Parameters
* Urey-Bradley Parameters
* Van der Waals Pair Parameters
* Van der Waals Parameters
A Tinker KEY file is composed of lines, each of which has a keyword,
which can be followed by zero or more parameters. This is the list of
keywords LAMMPS knows how to parse and use in the same manner Tinker
does. Any other keywords are skipped. The value following the equal
sign is the default value for the keyword if it is not specified, or
if the keyfile in the :doc:`pair_coeff <pair_coeff>` command is
specified as NULL:
* a-axis = 0.0
* b-axis = 0.0
* c-axis = 0.0
* ctrn-cutoff = 6.0
* ctrn-taper = 0.9 * ctrn-cutoff
* cutoff
* delta-halgren = 0.07
* dewald = no use of long-range dispersion unless specified
* dewald-alpha = 0.4
* dewald-cutoff = 7.0
* dispersion-cutoff = 9.0
* dispersion-taper = 9.0 * dispersion-cutoff
* dpme-grid
* dpme-order = 4
* ewald = no use of long-range electrostatics unless specified
* ewald-alpha = 0.4
* ewald-cutoff = 7.0
* gamma-halgren = 0.12
* mpole-cutoff = 9.0
* mpole-taper 0.65 * mpole-cutoff
* pcg-guess = enabled (default)
* pcg-noguess = disable pcg-guess
* pcg-noprecond = disable pcg-precond
* pcg-peek = 1.0
* pcg-precond = enabled (default)
* pewald-alpha = 0.4
* pme-grid
* pme-order = 5
* polar-eps = 1.0e-6
* polar-iter = 100
* polar-predict = no prediction operation unless specified
* ppme-order = 5
* repulsion-cutoff = 6.0
* repulsion-taper = 0.9 * repulsion-cutoff
* taper
* usolve-cutoff = 4.5
* usolve-diag = 2.0
* vdw-cutoff = 9.0
* vdw-taper = 0.9 * vdw-cutoff
----------
Tinker2lmp.py tool
This conversion tool is found in the tools/tinker directory.
As listed in examples/amoeba/README, these commands produce
the data files found in examples/amoeba, and also illustrate
all the options available to use with the tinker2lmp.py script:
.. code-block:: bash
% python tinker2lmp.py -xyz water_dimer.xyz -amoeba amoeba_water.prm -data data.water_dimer.amoeba # AMOEBA non-periodic system
% python tinker2lmp.py -xyz water_dimer.xyz -hippo hippo_water.prm -data data.water_dimer.hippo # HIPPO non-periodic system
% python tinker2lmp.py -xyz water_box.xyz -amoeba amoeba_water.prm -data data.water_box.amoeba -pbc 18.643 18.643 18.643 # AMOEBA periodic system
% python tinker2lmp.py -xyz water_box.xyz -hippo hippo_water.prm -data data.water_box.hippo -pbc 18.643 18.643 18.643 # HIPPO periodic system
% python tinker2lmp.py -xyz ubiquitin.xyz -amoeba amoeba_ubiquitin.prm -data data.ubiquitin.new -pbc 54.99 41.91 41.91 -bitorsion bitorsion.ubiquitin.data.new # system with bitorsions
Switches and their arguments may be specified in any order.
The -xyz switch is required and specifies an input XYZ file as an
argument. The format of this file is an extended XYZ format used by
Tinker for its input. Example *.xyz files are in the examples/amoeba
directory. The file lists the atoms in the system. Each atom has the
following information: Tinker species name (ignored by LAMMPS), xyz
coordinates, Tinker numeric type, and a list of atom IDs the atom is
bonded to.
NOTE: is this a Tinker-unique augmented XYZ format or standard? Where
can a LAMMPS user get or generate this file for a system they want
to simulate?
The -amoeba or -hippo switch is required. It specifies an input
AMOEBA or HIPPO PRM force field file as an argument. This should be
the same file used by the :doc:`pair_style <pair_style>` command in
the input script.
The -data switch is required. It specifies an output file name for
the LAMMPS data file that will be produced.
For periodic systems, the -pbc switch is required. It specifies the
periodic box size for each dimension (x,y,z). For a Tinker simulation
these are specified in the KEY file.
NOTE: What about a system with a free surface. What about a triclinic
box.
The -bitorsion switch is only needed if the system contains Tinker
bitorsion interactions. The data for each type of bitorsion
interaction will be written to the specified file, and read by the
:doc:`fix amoeba/bitorsion <fix_amoeba_bitorsion>` command. The data
includes 2d arrays of values to which splines are fit, and thus is not
compatible with the LAMMPS data file format.
----------

137
doc/src/angle_amoeba.rst Normal file
View File

@ -0,0 +1,137 @@
.. index:: angle_style amoeba
angle_style amoeba command
==========================
Syntax
""""""
.. code-block:: LAMMPS
angle_style amoeba
Examples
""""""""
.. code-block:: LAMMPS
angle_style amoeba
angle_coeff * 75.0 -25.0 1.0 0.3 0.02 0.003
angle_coeff * ba 3.6551 24.895 1.0119 1.5228
angle_coeff * ub -7.6 1.5537
Description
"""""""""""
The *amoeba* angle style uses the potential
.. math::
E & = E_a + E_{ba} + E_{ub} \\
E_a = K_2\left(\theta - \theta_0\right)^2 + K_3\left(\theta - \theta_0\right)^3 + K_4\left(\theta - \theta_0\right)^4 + K_5\left(\theta - \theta_0\right)^5 + K_6\left(\theta - \theta_0\right)^6
E_{ba} & = N_1 (r_{ij} - r_1) (\theta - \theta_0) + N_2(r_{jk} - r_2)(\theta - \theta_0)
E_{ub} & = K_{ub} (r_{ik} - r_{ub})^2)
where :math:`E_a` is the angle term, :math:`E_{ba}` is a bond-angle
term, ::math:`E_{ub}` is a Urey-Bradley bond term, :math:`\theta_0` is
the equilibrium angle, :math:`r_1` and :math:`r_2` are the equilibrium
bond lengths, and :math:`r_{ub}` is the equilibrium Urey-Bradley bond
length between atoms I and K in the angle IJK.
These formulas match how the Tinker MD code does its angle
calculations for the AMOEBA and HIPPO force fields. See the
doc:`Howto amoeba <howto_ameoba>` doc page for more information about
the implemention of AMOEBA and HIPPO in LAMMPS.
Note that the :math:`E_a` and :math:`E_{ba}` formulas are identical to
those used for the :doc:`angle_style class2/p6 <angle_class2_p6>`
command, however there is no bond-bond cross term formula for
:math:`E_{bb}`. Additionally, there is a :math:`E_{ub}` term for a
Urey-Bradley bond. It is effectively a harmonic bond between the I
and K atoms of angle IJK, even though that bond is not enumerated in
the "Bonds" section of the data file.
There are also two ways that Tinker computes the angle :math:`\theta`
in the :math:`E_a` formula. The first is the standard way of treating
IJK as an "in-plane" angle. The second is an "out-of-plane" method
which Tinker may use if the center atom J in the angle is bonded to
one additional atom in addition to I and K. In this case, all 4 atoms
are used to compute the :math:`E_a` formula, resulting in forces on
all 4 atoms. In the Tinker PRM file, these 2 options are denoted by
"angle" versus "anglep" entries in the "Angle Bending Parameters"
section of the PRM force field file. The *pflag* coefficient
described below selects between the 2 options.
----------
Coefficients for the :math:`E_a`, :math:`E_{bb}`, and :math:`E_{ub}`
formulas must be defined for each angle type via the :doc:`angle_coeff
<angle_coeff>` command as in the example above, or in the data file or
restart files read by the :doc:`read_data <read_data>` or
:doc:`read_restart <read_restart>` commands.
These are the 8 coefficients for the :math:`E_a` formula:
* pflag = 0 or 1
* ublag = 0 or 1
* :math:`\theta_0` (degrees)
* :math:`K_2` (energy)
* :math:`K_3` (energy)
* :math:`K_4` (energy)
* :math:`K_5` (energy)
* :math:`K_6` (energy)
A pflag value of 0 vs 1 selects between the "in-plane" and
"out-of-plane" options described above. Ubflag is 1 if there is a
Urey-Bradley term associated with this angle type, else it is 0.
:math:`\theta_0` is specified in degrees, but LAMMPS converts it to
radians internally; hence the various :math:`K` are effectively energy
per radian\^2 or radian\^3 or radian\^4 or radian\^5 or radian\^6.
For the :math:`E_{ba}` formula, each line in a :doc:`angle_coeff
<angle_coeff>` command in the input script lists 5 coefficients, the
first of which is "ba" to indicate they are BondAngle coefficients.
In a data file, these coefficients should be listed under a "BondAngle
Coeffs" heading and you must leave out the "ba", i.e. only list 4
coefficients after the angle type.
* ba
* :math:`N_1` (energy/distance\^2)
* :math:`N_2` (energy/distance\^2)
* :math:`r_1` (distance)
* :math:`r_2` (distance)
The :math:`\theta_0` value in the :math:`E_{ba}` formula is not specified,
since it is the same value from the :math:`E_a` formula.
For the :math:`E_{ub}` formula, each line in a :doc:`angle_coeff
<angle_coeff>` command in the input script lists 3 coefficients, the
first of which is "ub" to indicate they are UreyBradley coefficients.
In a data file, these coefficients should be listed under a
"UreyBradley Coeffs" heading and you must leave out the "ub",
i.e. only list 2 coefficients after the angle type.
* ub
* :math:`K_{ub}` (energy/distance\^2)
* :math:`r_{ub}` (distance)
----------
Restrictions
""""""""""""
This angle style can only be used if LAMMPS was built with the AMOEBA
package. See the :doc:`Build package <Build_package>` doc page for
more info.
Related commands
""""""""""""""""
:doc:`angle_coeff <angle_coeff>`
Default
"""""""
none

2
doc/src/angle_class2.rst
View File

@ -24,7 +24,7 @@ Examples
.. code-block:: LAMMPS
angle_style class2
angle_coeff * 75.0
angle_coeff * 75.0 25.0 0.3 0.002
angle_coeff 1 bb 10.5872 1.0119 1.5228
angle_coeff * ba 3.6551 24.895 1.0119 1.5228

1
doc/src/angle_style.rst
View File

@ -73,6 +73,7 @@ of (g,i,k,o,t) to indicate which accelerated styles exist.
* :doc:`zero <angle_zero>` - topology but no interactions
* :doc:`hybrid <angle_hybrid>` - define multiple styles of angle interactions
* :doc:`amoeba <angle_amoeba>` - AMOEBA angle
* :doc:`charmm <angle_charmm>` - CHARMM angle
* :doc:`class2 <angle_class2>` - COMPASS (class 2) angle
* :doc:`class2/p6 <angle_class2>` - COMPASS (class 2) angle expanded to 6th order

51
doc/src/atom_style.rst
View File

@ -10,7 +10,7 @@ Syntax
atom_style style args
* style = *angle* or *atomic* or *body* or *bond* or *charge* or *dipole* or *dpd* or *edpd* or *electron* or *ellipsoid* or *full* or *line* or *mdpd* or *molecular* or *oxdna* or *peri* or *smd* or *sph* or *sphere* or *spin* or *tdpd* or *tri* or *template* or *hybrid*
* style = *amoeba* or *angle* or *atomic* or *body* or *bond* or *charge* or *dipole* or *dpd* or *edpd* or *electron* or *ellipsoid* or *full* or *line* or *mdpd* or *molecular* or *oxdna* or *peri* or *smd* or *sph* or *sphere* or *spin* or *tdpd* or *tri* or *template* or *hybrid*
.. parsed-literal::
@ -77,6 +77,8 @@ coordinates, velocities, atom IDs and types. See the
:doc:`set <set>` commands for info on how to set these various
quantities.
+--------------+-----------------------------------------------------+--------------------------------------+
| *amoeba* | molecular + charge + 1/5 neighbors | AMOEBA/HIPPO polarized force fields |
+--------------+-----------------------------------------------------+--------------------------------------+
| *angle* | bonds and angles | bead-spring polymers with stiffness |
+--------------+-----------------------------------------------------+--------------------------------------+
@ -134,15 +136,18 @@ quantities.
.. note::
It is possible to add some attributes, such as a molecule ID, to
atom styles that do not have them via the :doc:`fix property/atom <fix_property_atom>` command. This command also
allows new custom attributes consisting of extra integer or
floating-point values to be added to atoms. See the :doc:`fix property/atom <fix_property_atom>` page for examples of cases
where this is useful and details on how to initialize, access, and
output the custom values.
atom styles that do not have them via the :doc:`fix property/atom
<fix_property_atom>` command. This command also allows new custom
attributes consisting of extra integer or floating-point values to
be added to atoms. See the :doc:`fix property/atom
<fix_property_atom>` page for examples of cases where this is
useful and details on how to initialize, access, and output the
custom values.
All of the above styles define point particles, except the *sphere*,
*ellipsoid*, *electron*, *peri*, *wavepacket*, *line*, *tri*, and
*body* styles, which define finite-size particles. See the :doc:`Howto spherical <Howto_spherical>` page for an overview of using
*body* styles, which define finite-size particles. See the
:doc:`Howto spherical <Howto_spherical>` page for an overview of using
finite-size particle models with LAMMPS.
All of the point-particle styles assign mass to particles on a
@ -266,16 +271,17 @@ showing the use of the *template* atom style versus *molecular*.
.. note::
When using the *template* style with a :doc:`molecule template <molecule>` that contains multiple molecules, you should
insure the atom types, bond types, angle_types, etc in all the
molecules are consistent. E.g. if one molecule represents H2O and
another CO2, then you probably do not want each molecule file to
define 2 atom types and a single bond type, because they will conflict
with each other when a mixture system of H2O and CO2 molecules is
defined, e.g. by the :doc:`read_data <read_data>` command. Rather the
H2O molecule should define atom types 1 and 2, and bond type 1. And
the CO2 molecule should define atom types 3 and 4 (or atom types 3 and
2 if a single oxygen type is desired), and bond type 2.
When using the *template* style with a :doc:`molecule template
<molecule>` that contains multiple molecules, you should insure the
atom types, bond types, angle_types, etc in all the molecules are
consistent. E.g. if one molecule represents H2O and another CO2,
then you probably do not want each molecule file to define 2 atom
types and a single bond type, because they will conflict with each
other when a mixture system of H2O and CO2 molecules is defined,
e.g. by the :doc:`read_data <read_data>` command. Rather the H2O
molecule should define atom types 1 and 2, and bond type 1. And
the CO2 molecule should define atom types 3 and 4 (or atom types 3
and 2 if a single oxygen type is desired), and bond type 2.
For the *body* style, the particles are arbitrary bodies with internal
attributes defined by the "style" of the bodies, which is specified by
@ -352,6 +358,8 @@ Many of the styles listed above are only enabled if LAMMPS was built
with a specific package, as listed below. See the :doc:`Build package
<Build_package>` page for more info.
The *amoeba* style is part of the AMOEBA package.
The *angle*, *bond*, *full*, *molecular*, and *template* styles are
part of the MOLECULE package.
@ -363,9 +371,11 @@ The *dipole* style is part of the DIPOLE package.
The *peri* style is part of the PERI package for Peridynamics.
The *oxdna* style is part of the CG-DNA package for coarse-grained simulation of DNA and RNA.
The *oxdna* style is part of the CG-DNA package for coarse-grained
simulation of DNA and RNA.
The *electron* style is part of the EFF package for :doc:`electronic force fields <pair_eff>`.
The *electron* style is part of the EFF package for :doc:`electronic
force fields <pair_eff>`.
The *dpd* style is part of the DPD-REACT package for dissipative
particle dynamics (DPD).
@ -376,7 +386,8 @@ dissipative particle dynamics (mDPD), and transport dissipative particle
dynamics (tDPD), respectively.
The *sph* style is part of the SPH package for smoothed particle
hydrodynamics (SPH). See `this PDF guide <PDF/SPH_LAMMPS_userguide.pdf>`_ to using SPH in LAMMPS.
hydrodynamics (SPH). See `this PDF guide
<PDF/SPH_LAMMPS_userguide.pdf>`_ to using SPH in LAMMPS.
The *mesont* style is part of the MESONT package.

4
doc/src/fix.rst
View File

@ -170,6 +170,8 @@ accelerated styles exist.
* :doc:`adapt/fep <fix_adapt_fep>` - enhanced version of fix adapt
* :doc:`addforce <fix_addforce>` - add a force to each atom
* :doc:`addtorque <fix_addtorque>` - add a torque to a group of atoms
* :doc:`amoeba/bitorsion <fix_amoeba_bitorsion>` - torsion/torsion terms in AMOEBA force field
* :doc:`amoeba/pitorsion <fix_amoeba_pitorsion>` - 6-body terms in AMOEBA force field
* :doc:`append/atoms <fix_append_atoms>` - append atoms to a running simulation
* :doc:`atc <fix_atc>` - initiates a coupled MD/FE simulation
* :doc:`atom/swap <fix_atom_swap>` - Monte Carlo atom type swapping
@ -194,7 +196,7 @@ accelerated styles exist.
* :doc:`box/relax <fix_box_relax>` - relax box size during energy minimization
* :doc:`charge/regulation <fix_charge_regulation>` - Monte Carlo sampling of charge regulation
* :doc:`client/md <fix_client_md>` - MD client for client/server simulations
* :doc:`cmap <fix_cmap>` - enables CMAP cross-terms of the CHARMM force field
* :doc:`cmap <fix_cmap>` - CMAP torsion/torsion terms in CHARMM force field
* :doc:`colvars <fix_colvars>` - interface to the collective variables "Colvars" library
* :doc:`controller <fix_controller>` - apply control loop feedback mechanism
* :doc:`deform <fix_deform>` - change the simulation box size/shape

167
doc/src/fix_amoeba_bitorsion.rst Normal file
View File

@ -0,0 +1,167 @@
.. index:: fix amoeba/bitorsion
fix amoeba/bitorsion command
================
Syntax
""""""
.. parsed-literal::
fix ID group-ID ameoba/bitorsion filename
* ID, group-ID are documented in :doc:`fix <fix>` command
* amoeba/bitorsion = style name of this fix command
* filename = force-field file with AMOEBA bitorsion coefficients
Examples
""""""""
.. code-block:: LAMMPS
fix bit all amoeba/bitorsion bitorsion.ubiquitin.data
read_data proteinX.data fix bit bitorsions BiTorsions
fix_modify bit energy yes
Description
"""""""""""
This command enables 5-body torsion/torsion interactions to be added
to simulations which use the AMOEBA and HIPPO force fields. It
matches how the Tinker MD code computes its torsion/torsion
interactions for the AMOEBA and HIPPO force fields. See the
doc:`Howto amoeba <howto_ameoba>` doc page for more information about
the implemention of AMOEBA and HIPPO in LAMMPS.
Bitorsion interactions add additional potential energy contributions
to pairs of overlapping phi-psi dihedrals of amino-acids, which are
important to properly represent their conformational behavior.
The examples/amoeba directory has a sample input script and data file
for ubiquitin, which illustrates use of the fix amoeba/bitorsion
command.
As in the example above, this fix should be used before reading a data
file that contains a listing of bitorsion interactions. The
*filename* specified should contain the bitorsion parameters for the
AMOEBA or HIPPO force field.
The data file read by the "read_data" must contain the topology of all
the bitorsion interactions, similar to the topology data for bonds,
angles, dihedrals, etc. Specially it should have a line like this in
its header section:
.. parsed-literal::
N bitorsions
where N is the number of bitorsion 5-body interactions. It should
also have a section in the body of the data file like this with N
lines:
.. parsed-literal::
BiTorsions
1 1 8 10 12 18 20
2 5 18 20 22 25 27
[...]
N 3 314 315 317 318 330
The first column is an index from 1 to N to enumerate the bitorsion
5-atom tuples; it is ignored by LAMMPS. The second column is the
"type" of the interaction; it is an index into the bitorsion force
field file. The remaining 5 columns are the atom IDs of the atoms in
the two 4-atom dihedrals that overlap to create the bitorsion 5-body
interaction. Note that the "bitorsions" and "BiTorsions" keywords for
the header and body sections match those specified in the read_data
command following the data file name; see the :doc:`read_data
<read_data>` page for more details.
The data file should be generated by using the
tools/tinker/tinker2lmp.py conversion script which creates a LAMMPS
data file from Tinker input files, including its PRM file which
contains the parameters necessary for computing bitorsion
interactions. The script must be invoked with the optional
"-bitorsion" flag to do this; see the example for the ubiquitin system
in the tools/tinker/README file. The same conversion script also
creates the file of bitorsion coefficient data which is read by this
command.
The potential energy associated with bitorsion interactions can be
output as described below. It can also be included in the total
potential energy of the system, as output by the :doc:`thermo_style
<thermo_style>` command, if the :doc:`fix_modify energy <fix_modify>`
command is used, as in the example above. See the note below about
how to include the bitorsion energy when performing an :doc:`energy
minimization <minimize>`.
----------
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
This fix writes the list of bitorsion interactions to :doc:`binary
restart files <restart>`. See the :doc:`read_restart <read_restart>`
command for info on how to re-specify a fix in an input script that
reads a restart file, so that the operation of the fix continues in an
uninterrupted fashion.
The :doc:`fix_modify <fix_modify>` *energy* option is supported by
this fix to add the potential energy of the bitorsion interactions to
both the global potential energy and peratom potential energies of the
system as part of :doc:`thermodynamic output <thermo_style>` or output
by the :doc:`compute pe/atom <compute_pe_atom>` command. The default
setting for this fix is :doc:`fix_modify energy yes <fix_modify>`.
The :doc:`fix_modify <fix_modify>` *virial* option is supported by
this fix to add the contribution due to the bitorsion interactions to
both the global pressure and per-atom stress of the system via the
:doc:`compute pressure <compute_pressure>` and :doc:`compute
stress/atom <compute_stress_atom>` commands. The former can be
accessed by :doc:`thermodynamic output <thermo_style>`. The default
setting for this fix is :doc:`fix_modify virial yes <fix_modify>`.
This fix computes a global scalar which can be accessed by various
:doc:`output commands <Howto_output>`. The scalar is the potential
energy discussed above. The scalar value calculated by this fix is
"extensive".
No parameter of this fix can be used with the *start/stop* keywords of
the :doc:`run <run>` command.
The forces due to this fix are imposed during an energy minimization,
invoked by the :doc:`minimize <minimize>` command.
The :doc:`fix_modify <fix_modify>` *respa* option is supported by this
fix. This allows to set at which level of the :doc:`r-RESPA
<run_style>` integrator the fix is adding its forces. Default is the
outermost level.
.. note::
For energy minimization, if you want the potential energy
associated with the bitorsion terms forces to be included in the
total potential energy of the system (the quantity being
minimized), you MUST not disable the :doc:`fix_modify <fix_modify>`
*energy* option for this fix.
Restrictions
""""""""""""
To function as expected this fix command must be issued *before* a
:doc:`read_data <read_data>` command but *after* a :doc:`read_restart
<read_restart>` command.
This fix can only be used if LAMMPS was built with the AMOEBA package.
See the :doc:`Build package <Build_package>` page for more info.
Related commands
""""""""""""""""
:doc:`fix_modify <fix_modify>`, :doc:`read_data <read_data>`
Default
"""""""
none

186
doc/src/fix_amoeba_pitorsion.rst Normal file
View File

@ -0,0 +1,186 @@
.. index:: fix amoeba/pitorsion
fix amoeba/pitorsion command
================
Syntax
""""""
.. parsed-literal::
fix ID group-ID ameoba/pitorsion
* ID, group-ID are documented in :doc:`fix <fix>` command
* amoeba/pitorsion = style name of this fix command
Examples
""""""""
.. code-block:: LAMMPS
fix pit all amoeba/pitorsion
read_data proteinX.data fix pit "pitorsion types" "PiTorsion Coeffs" &
fix pit pitorsions PiTorsions
fix_modify pit energy yes
Description
"""""""""""
This command enables 6-body pitorsion interactions to be added to
simulations which use the AMOEBA and HIPPO force fields. It matches
how the Tinker MD code computes its pitorsion interactions for the
AMOEBA and HIPPO force fields. See the doc:`Howto amoeba
<howto_ameoba>` doc page for more information about the implemention
of AMOEBA and HIPPO in LAMMPS.
Pitorsion interactions add additional potential energy contributions
to 6-tuples of atoms IJKLMN which have a bond between atoms K and L,
where both K and L are additionally bonded to exactly two other atoms.
Namely K is also bonded to I and J. And L is also bonded to M and N.
The examples/amoeba directory has a sample input script and data file
for ubiquitin, which illustrates use of the fix amoeba/pitorsion
command.
As in the example above, this fix should be used before reading a data
file that contains a listing of pitorsion interactions. The
*filename* specified should contain the pitorsion parameters for the
AMOEBA or HIPPO force field.
The data file read by the "read_data" must contain the topology of all
the pitorsion interactions, similar to the topology data for bonds,
angles, dihedrals, etc. Specially it should have two lines like these
in its header section:
.. parsed-literal::
M pitorsion types
N pitorsions
where N is the number of pitorsion 5-body interactions and M is the
number of pitorsion types. It should also have two sections in the body
of the data file like these with M and N lines:
.. parsed-literal::
PiTorsion Coeffs
1 6.85
2 10.2
[...]
M 6.85
.. parsed-literal::
PiTorsions
1 1 8 10 12 18 20
2 5 18 20 22 25 27
[...]
N 3 314 315 317 318 330
For PiTorsion Coeffs, the first column is an index from 1 to M to
enumerate the pitorsion types. The second column is the single
coefficient needed for each type.
For PiTorsions, the first column is an index from 1 to N to enumerate
the pitorsion 5-atom tuples; it is ignored by LAMMPS. The second
column is the "type" of the interaction; it is an index into the
PiTorsion Coeffs. The remaining 5 columns are the atom IDs of the
atoms in the two 4-atom dihedrals that overlap to create the pitorsion
5-body interaction.
Note that the "pitorsion types" and "pitorsions" and "PiTorsion
Coeffs" and "PiTorsions" keywords for the header and body sections
match those specified in the read_data command following the data file
name; see the :doc:`read_data <read_data>` page for more details.
A data file containing pitorsion 6-body interactions can be generated
from using the tinker2lmp.py script in the tools/tinker directory of
the LAMMPS distribution. See the examples for its use for the
ubiquitin system in the tools/tinker/README file.
The data file should be generated by using the
tools/tinker/tinker2lmp.py conversion script which creates a LAMMPS
data file from Tinker input files, including its PRM file which
contains the parameters necessary for computing pitorsion
interactions. See the example for the ubiquitin system in the
tools/tinker/README file.
The potential energy associated with pitorsion interactions can be
output as described below. It can also be included in the total
potential energy of the system, as output by the :doc:`thermo_style
<thermo_style>` command, if the :doc:`fix_modify energy <fix_modify>`
command is used, as in the example above. See the note below about
how to include the pitorsion energy when performing an :doc:`energy
minimization <minimize>`.
----------
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
This fix writes the list of pitorsion interactions to :doc:`binary
restart files <restart>`. See the :doc:`read_restart <read_restart>`
command for info on how to re-specify a fix in an input script that
reads a restart file, so that the operation of the fix continues in an
uninterrupted fashion.
The :doc:`fix_modify <fix_modify>` *energy* option is supported by
this fix to add the potential energy of the pitorsion interactions to
both the global potential energy and peratom potential energies of the
system as part of :doc:`thermodynamic output <thermo_style>` or output
by the :doc:`compute pe/atom <compute_pe_atom>` command. The default
setting for this fix is :doc:`fix_modify energy yes <fix_modify>`.
The :doc:`fix_modify <fix_modify>` *virial* option is supported by
this fix to add the contribution due to the pitorsion interactions to
both the global pressure and per-atom stress of the system via the
:doc:`compute pressure <compute_pressure>` and :doc:`compute
stress/atom <compute_stress_atom>` commands. The former can be
accessed by :doc:`thermodynamic output <thermo_style>`. The default
setting for this fix is :doc:`fix_modify virial yes <fix_modify>`.
This fix computes a global scalar which can be accessed by various
:doc:`output commands <Howto_output>`. The scalar is the potential
energy discussed above. The scalar value calculated by this fix is
"extensive".
No parameter of this fix can be used with the *start/stop* keywords of
the :doc:`run <run>` command.
The forces due to this fix are imposed during an energy minimization,
invoked by the :doc:`minimize <minimize>` command.
The :doc:`fix_modify <fix_modify>` *respa* option is supported by this
fix. This allows to set at which level of the :doc:`r-RESPA
<run_style>` integrator the fix is adding its forces. Default is the
outermost level.
.. note::
For energy minimization, if you want the potential energy
associated with the pitorsion terms forces to be included in the
total potential energy of the system (the quantity being
minimized), you MUST not disable the :doc:`fix_modify <fix_modify>`
*energy* option for this fix.
Restrictions
""""""""""""
To function as expected this fix command must be issued *before* a
:doc:`read_data <read_data>` command but *after* a :doc:`read_restart
<read_restart>` command.
This fix can only be used if LAMMPS was built with the AMOEBA package.
See the :doc:`Build package <Build_package>` page for more info.
Related commands
""""""""""""""""
:doc:`fix_modify <fix_modify>`, :doc:`read_data <read_data>`
Default
"""""""
none

66
doc/src/fix_cmap.rst
View File

@ -26,13 +26,13 @@ Examples
Description
"""""""""""
This command enables CMAP cross-terms to be added to simulations which
use the CHARMM force field. These are relevant for any CHARMM model
of a peptide or protein sequences that is 3 or more amino-acid
residues long; see :ref:`(Buck) <Buck>` and :ref:`(Brooks) <Brooks2>`
for details, including the analytic energy expressions for CMAP
interactions. The CMAP cross-terms add additional potential energy
contributions to pairs of overlapping phi-psi dihedrals of
This command enables CMAP 5-body interactions to be added to
simulations which use the CHARMM force field. These are relevant for
any CHARMM model of a peptide or protein sequences that is 3 or more
amino-acid residues long; see :ref:`(Buck) <Buck>` and :ref:`(Brooks)
<Brooks2>` for details, including the analytic energy expressions for
CMAP interactions. The CMAP 5-body terms add additional potential
energy contributions to pairs of overlapping phi-psi dihedrals of
amino-acids, which are important to properly represent their
conformational behavior.
@ -47,15 +47,15 @@ lammps/potentials directory: charmm22.cmap and charmm36.cmap.
The data file read by the "read_data" must contain the topology of all
the CMAP interactions, similar to the topology data for bonds, angles,
dihedrals, etc. Specially it should have a line like this
in its header section:
dihedrals, etc. Specially it should have a line like this in its
header section:
.. parsed-literal::
N crossterms
where N is the number of CMAP cross-terms. It should also have a section
in the body of the data file like this with N lines:
where N is the number of CMAP 5-body interactions. It should also
have a section in the body of the data file like this with N lines:
.. parsed-literal::
@ -66,28 +66,29 @@ in the body of the data file like this with N lines:
[...]
N 3 314 315 317 318 330
The first column is an index from 1 to N to enumerate the CMAP terms;
it is ignored by LAMMPS. The second column is the "type" of the
interaction; it is an index into the CMAP force field file. The
The first column is an index from 1 to N to enumerate the CMAP 5-atom
tuples; it is ignored by LAMMPS. The second column is the "type" of
the interaction; it is an index into the CMAP force field file. The
remaining 5 columns are the atom IDs of the atoms in the two 4-atom
dihedrals that overlap to create the CMAP 5-body interaction. Note
that the "crossterm" and "CMAP" keywords for the header and body
sections match those specified in the read_data command following the
data file name; see the :doc:`read_data <read_data>` page for
more details.
dihedrals that overlap to create the CMAP interaction. Note that the
"crossterm" and "CMAP" keywords for the header and body sections match
those specified in the read_data command following the data file name;
see the :doc:`read_data <read_data>` page for more details.
A data file containing CMAP cross-terms can be generated from a PDB
file using the charmm2lammps.pl script in the tools/ch2lmp directory
of the LAMMPS distribution. The script must be invoked with the
optional "-cmap" flag to do this; see the tools/ch2lmp/README file for
more information.
A data file containing CMAP 5-body interactions can be generated from
a PDB file using the charmm2lammps.pl script in the tools/ch2lmp
directory of the LAMMPS distribution. The script must be invoked with
the optional "-cmap" flag to do this; see the tools/ch2lmp/README file
for more information. The same conversion script also creates the
file of CMAP coefficient data which is read by this command.
The potential energy associated with CMAP interactions can be output
as described below. It can also be included in the total potential
energy of the system, as output by the
:doc:`thermo_style <thermo_style>` command, if the :doc:`fix_modify energy <fix_modify>` command is used, as in the example above. See
the note below about how to include the CMAP energy when performing an
:doc:`energy minimization <minimize>`.
energy of the system, as output by the :doc:`thermo_style
<thermo_style>` command, if the :doc:`fix_modify energy <fix_modify>`
command is used, as in the example above. See the note below about
how to include the CMAP energy when performing an :doc:`energy
minimization <minimize>`.
----------
@ -134,10 +135,11 @@ outermost level.
.. note::
If you want the potential energy associated with the CMAP terms
forces to be included in the total potential energy of the system
(the quantity being minimized), you MUST not disable the
:doc:`fix_modify <fix_modify>` *energy* option for this fix.
For energy minimization, if you want the potential energy
associated with the CMAP terms forces to be included in the total
potential energy of the system (the quantity being minimized), you
MUST not disable the :doc:`fix_modify <fix_modify>` *energy* option
for this fix.
Restrictions
""""""""""""

77
doc/src/improper_amoeba.rst Normal file
View File

@ -0,0 +1,77 @@
.. index:: improper_style amoeba
improper_style harmonic command
===============================
Syntax
""""""
.. code-block:: LAMMPS
improper_style amoeba
Examples
""""""""
.. code-block:: LAMMPS
improper_style amoeba
improper_coeff 1 49.6
Description
"""""""""""
The *amoeba* improper style uses the potential
.. math::
E = K (\chi)^2
where :math:`\chi` is the improper angle and :math:`K` is a prefactor.
Note that the usual 1/2 factor is included in :math:`K`.
This formula seems like a simplified version of the formula for the
:doc:`improper_style harmonic <improper_harmonic>` command with
:math:`\chi_0` = 0.0. However the computation of the angle
:math:`\chi` is done differently to match how the Tinker MD code
computes its out-of-plane improper for the AMOEBA and HIPPO force
fields. See the doc:`Howto amoeba <howto_ameoba>` doc page for more
information about the implemention of AMOEBA and HIPPO in LAMMPS.
If the 4 atoms in an improper quadruplet (listed in the data file read
by the :doc:`read_data <read_data>` command) are ordered I,J,K,L then
atoms I,K,L form a plane and atom J is out-of-place. The angle
:math:`\chi_0` is computed as the Allinger angle which is defined as
the angle between the plane of I,K,L, and the vector from atom I to
atom J.
The following coefficient must be defined for each improper type via
the :doc:`improper_coeff <improper_coeff>` command as in the example
above, or in the data file or restart files read by the
:doc:`read_data <read_data>` or :doc:`read_restart <read_restart>`
commands:
* :math:`K` (energy)
Note that the angle :math:`\chi` is computed in radians; hence
:math:`K` is effectively energy per radian\^2.
----------
Restrictions
""""""""""""
This improper style can only be used if LAMMPS was built with the
AMOEBA package. See the :doc:`Build package <Build_package>` doc page
for more info.
Related commands
""""""""""""""""
:doc:`improper_coeff <improper_coeff>`, `improper_harmonic
:doc:<improper_harmonic>`
Default
"""""""
none

1
doc/src/improper_style.rst
View File

@ -77,6 +77,7 @@ more of (g,i,k,o,t) to indicate which accelerated styles exist.
* :doc:`zero <improper_zero>` - topology but no interactions
* :doc:`hybrid <improper_hybrid>` - define multiple styles of improper interactions
* :doc:`amoeba <improper_amoeba>` - AMOEBA out-of-plane improper
* :doc:`class2 <improper_class2>` - COMPASS (class 2) improper
* :doc:`cossq <improper_cossq>` - improper with a cosine squared term
* :doc:`cvff <improper_cvff>` - CVFF improper

163
doc/src/pair_amoeba.rst
View File

@ -33,100 +33,73 @@ Examples
Additional info
"""""""""""""""
doc:`Howto amoeba <howto_ameoba>`
doc:`bond_style amoeba <bond_amoeba>`
doc:`angle_style amoeba <angle_amoeba>`
doc:`dihedral_style amoeba <dihedral_amoeba>`
examples/amoeba
tools/amoeba
potentials/\*.prm.ameoba
potentials/\*.key.ameoba
potentials/\*.prm.hippo
potentials/\*.key.hippo
* doc:`Howto amoeba <howto_ameoba>`
* examples/amoeba
* tools/amoeba
* potentials/*.amoeba
* potentials/*.hippo
Description
"""""""""""
NOTE: edit this paragraph however you wish including adding or changing
citations (see bottom of this file)
The *amoeba* style computes the AMOEBA polarizeable field formulated
by Jay Ponder's group at U Washington at St Louis :ref:`(Ren)
<amoeba-Ren>`, :ref:`(Shi) <amoeba-Shi>`. The *hippo* style
computes the HIPPO polarizeable force field , an extension to AMOEBA,
formulated by Josh Rackers and collaborators in the Ponder group
:ref:`(Rackers) <amoeba-Rackers>`.
by Jay Ponder's group at the U Washington at St Louis :ref:`(Ren)
<amoeba-Ren>`, :ref:`(Shi) <amoeba-Shi>`. The *hippo* style computes
the HIPPO polarizeable force field, an extension to AMOEBA, formulated
by Josh Rackers and collaborators in the Ponder group :ref:`(Rackers)
<amoeba-Rackers>`.
These force fields can be used when polarization effects are desired
in simulations of water, organic molecules, and biomolecules including
proteins, provided that parameterizations (force field files) are
available for the systems you are interested in. Files in the LAMMPS
potentials with a "amoeba" or "hippo" suffix can be used. The Tinker
distribution and website may have other such files.
proteins, provided that parameterizations (Tinker PRM force field
files) are available for the systems you are interested in. Files in
the LAMMPS potentials directory with a "amoeba" or "hippo" suffix can
be used. The Tinker distribution and website have additional force
field files as well.
NOTE: replace these LaTeX formulas with a list of AMOEBA and HIPPO
terms in some simple notation. If desired, a more detailed
mathematical expression for each term could also be included below
the initial 2 formulas.
The AMOEBA force field can be written as a collection of terms
As discussed on the doc:`Howto amoeba <howto_ameoba>` doc page, the
intermolecular (non-bonded) portion of the AMOEBA force field contains
these terms:
.. math::
E = 4 \epsilon \left[ \left(\frac{\sigma}{r}\right)^{12} -
\left(\frac{\sigma}{r}\right)^6 \right]
\qquad r < r_c
U_{amoeba} = U_{hal} + U_{multipole} + U_{polar}
:math:`r_c` is the cutoff, blah is the dipole, etc
The HIPPO force field is similar but alters some of the terms
while the HIPPO force field contains these terms:
.. math::
E = 4 \epsilon \left[ \left(\frac{\sigma}{r}\right)^{12} -
\left(\frac{\sigma}{r}\right)^6 \right]
\qquad r < r_c
U_{hippo} = U_{repulsion} + U_{dispersion} + U_{multipole} + U_{polar} + U_{qxfer}
:math:`r_c` is the cutoff, blah is the dipole, etc
Conceptually, these terms compute the following interactions:
* :math:`U_{hal}` = buffered 14-7 van der Waals with offsets applied to hydrogen atoms
* :math:`U_{repulsion}` = Pauli repulsion due to rearrangement of electron density
* :math:`U_{dispersion}` = dispersion between correlated, instantaneous induced dipole moments
* :math:`U_{multipole}` = electrostatics between permanent point charges, dipoles, and quadrupoles
* :math:`U_{polar}` = electronic polarization bewteen induced point dipoles
* :math:`U_{qxfer}` = charge transfer effects
NOTE: Add a sentence for each term to explain what physical effects
the FFs are encoding. E.g. The polar term iteratively computes an
induced dipole for each atom, then calculates dipole-dipole
interactions ...
See the AMOEBA and HIPPO papers for further details.
Note that the AMOEBA versus HIPPO force fields typically compute the
same term differently using their own formulas. The references on
this doc page give full details for both force fields.
.. note::
The AMOEBA and HIPPO force fields compute long-range charge, dipole,
and quadrupole interactions (NOTE: also long-range dispersion?).
However, unlike other models with long-range interactions in LAMMPS,
this does not require use of a KSpace style via the
and quadrupole interactions as well as long-range dispersion
effects. However, unlike other models with long-range interactions
in LAMMPS, this does not require use of a KSpace style via the
:doc:`kspace_style <kspace_style>` command. That is because for
AMOEBA and HIPPO the long-range computations are intertwined with
the pairwise compuations. So these pair style include both short-
the pairwise computations. So these pair style include both short-
and long-range computations. This means the energy and virial
computed by the pair style as well as the "Pair" timing reported by
LAMMPS will include the long-range components.
.. note::
As explained on the :doc:`Howto amoeba <howto_ameoba>` doc page, use
of these pair styles to run a simulation with the AMOEBA or HIPPO
force fields requires your input script to use the :doc:`atom_style
hybrid full amoeba <atom_style>` atom style, AMOEBA versions of
bond/angle/dihedral styles, the :doc:`special_bonds one/five
<special_bonds>` option, and the :doc:`fix property/atom one/five
<fix_property>` command to define several additional per-atom
properties. The latter requires a "Tinker Types" section be
included in the LAMMPS data file. This can be auto-generated using
the tools/amoeba/tinker2lmp.py Python script. See the :doc:`Howto
amoeba <howto_ameoba>` doc page and tools/amoeba/README file for
details on using that tool.
LAMMPS will include the long-range calculations.
The implementation of the AMOEBA and HIPPO force fields in LAMMPS was
done using code provided by the Ponder group from their `Tinker MD
code <https://dasher.wustl.edu/tinker/>`_ written in F90.
done using F90 code provided by the Ponder group from their `Tinker MD
code <https://dasher.wustl.edu/tinker/>`_.
NOTE: what version of AMOEBA and HIPPO does LAMMPS implement?
@ -140,22 +113,22 @@ Only a single pair_coeff command is used with either the *amoeba* and
pair_coeff * * ../potentials/protein.prm.amoeba ../potentials/protein.key.amoeba
pair_coeff * * ../potentials/water.prm.hippo ../potentials/water.key.hippo
See examples of these files in the potentials directory.
Examples of the PRM files are in the potentials directory with an
*.amoeba or *.hippo suffix. The examples/amoeba directory has
examples of both PRM and KEY files.
The format of a PRM file is a collection of sections, each with
multiple lines. These are the sections which LAMMPS recognizes:
A Tinker PRM file composed of sections, each of which has multiple
lines. A Tinker KEY file is composed of lines, each of which has a keyword,
which can be followed by zero or more parameters.
NOTE: need to list these, possibly there are some LAMMPS skips
or which need to be removed for LAMMPS to use the PRM file?
The list of PRM sections and KEY keywords which LAMMPS recognizes are
listed on the doc:`Howto amoeba <howto_ameoba>` doc page. If not
recognized, the section or keyword is skipped.
The format of a KEY file is a series of lines, with one parameter and
its value per line. These are the parameters which LAMMPS recognizes:
NOTE: need to list these, possibly there are some keywords LAMMPS skips
or which need to be removed for LAMMPS to use the KEY file?
NOTE: Any other info about PRM and KEY files we should explain
for LAMMPS users here?
Note that if the KEY file is specified as NULL, then no file is
required; default values for various AMOEBA/HIPPO settings are used.
The doc:`Howto amoeba <howto_ameoba>` doc page also gives the default
settings.
----------
@ -184,14 +157,38 @@ enabled if LAMMPS was built with that package. See the :doc:`Build
package <Build_package>` doc page for more info.
The AMOEBA and HIPPO potential (PRM) and KEY files provided with
LAMMPS in the potentials directory are Tinker files parameterized for
Tinker units. Their numeric parameters are converted by LAMMPS to its
real units :doc:`units <units>`. You can only use these pair styles
with real units.
LAMMPS in the potentials and examples/amoeva directories are Tinker
files parameterized for Tinker units. Their numeric parameters are
converted by LAMMPS to its real units :doc:`units <units>`. Thus uou
can only use these pair styles with real units.
These potentials do not yet calculate per-atom energy or virial
contributions.
As explained on the :doc:`Howto amoeba <howto_ameoba>` doc page, use
of these pair styles to run a simulation with the AMOEBA or HIPPO
force fields requires several things.
The first is a data file generated by the tools/tinker/tinker2lmp.py
conversion script which uses Tinker file force field file input to
create a data file compatible with LAMMPS.
The second is use of these commands:
* :doc:`atom_style amoeba <atom_style>`
* :doc:`fix property/atom <fix_property_atom>`
* :doc:`special_bonds one/five <special_bonds>`
And third, depending on the model being simulated, these
commands for intramolecular interactions may also be required:
* :doc:`bond_style class2 <bond_class2>`
* :doc:`angle_style amoeba <angle_amoeba>`
* :doc:`dihedral_style fourier <dihedral_fourier>`
* :doc:`improper_style amoeba <improper_amoeba>`
* :doc:`fix amoeba/pitorsion <fix_pitorsion>`
* :doc:`fix amoeba/bitorsion <fix_bitorsion>`
----------
Related commands

104
doc/src/special_bonds.rst
View File

@ -11,7 +11,7 @@ Syntax
special_bonds keyword values ...
* one or more keyword/value pairs may be appended
* keyword = *amber* or *charmm* or *dreiding* or *fene* or *lj/coul* or *lj* or *coul* or *angle* or *dihedral*
* keyword = *amber* or *charmm* or *dreiding* or *fene* or *lj/coul* or *lj* or *coul* or *angle* or *dihedral* or *one/five*
.. parsed-literal::
@ -27,6 +27,7 @@ Syntax
w1,w2,w3 = weights (0.0 to 1.0) on pairwise Coulombic interactions
*angle* value = *yes* or *no*
*dihedral* value = *yes* or *no*
*one/five* value = *yes* or *no*
Examples
""""""""
@ -45,10 +46,10 @@ Description
Set weighting coefficients for pairwise energy and force contributions
between pairs of atoms that are also permanently bonded to each other,
either directly or via one or two intermediate bonds. These weighting
factors are used by nearly all :doc:`pair styles <pair_style>` in LAMMPS
that compute simple pairwise interactions. Permanent bonds between
atoms are specified by defining the bond topology in the data file
read by the :doc:`read_data <read_data>` command. Typically a
factors are used by nearly all :doc:`pair styles <pair_style>` in
LAMMPS that compute simple pairwise interactions. Permanent bonds
between atoms are specified by defining the bond topology in the data
file read by the :doc:`read_data <read_data>` command. Typically a
:doc:`bond_style <bond_style>` command is also used to define a bond
potential. The rationale for using these weighting factors is that
the interaction between a pair of bonded atoms is all (or mostly)
@ -58,31 +59,34 @@ atoms should be excluded (or reduced by a weighting factor).
.. note::
These weighting factors are NOT used by :doc:`pair styles <pair_style>` that compute many-body interactions, since the
"bonds" that result from such interactions are not permanent, but are
created and broken dynamically as atom conformations change. Examples
of pair styles in this category are EAM, MEAM, Stillinger-Weber,
Tersoff, COMB, AIREBO, and ReaxFF. In fact, it generally makes no
sense to define permanent bonds between atoms that interact via these
potentials, though such bonds may exist elsewhere in your system,
e.g. when using the :doc:`pair_style hybrid <pair_hybrid>` command.
Thus LAMMPS ignores special_bonds settings when many-body potentials
are calculated. Please note, that the existence of explicit bonds
for atoms that are described by a many-body potential will alter the
neighbor list and thus can render the computation of those interactions
invalid, since those pairs are not only used to determine direct
pairwise interactions but also neighbors of neighbors and more.
The recommended course of action is to remove such bonds, or - if
that is not possible - use a special bonds setting of 1.0 1.0 1.0.
These weighting factors are NOT used by :doc:`pair styles
<pair_style>` that compute many-body interactions, since the
"bonds" that result from such interactions are not permanent, but
are created and broken dynamically as atom conformations change.
Examples of pair styles in this category are EAM, MEAM,
Stillinger-Weber, Tersoff, COMB, AIREBO, and ReaxFF. In fact, it
generally makes no sense to define permanent bonds between atoms
that interact via these potentials, though such bonds may exist
elsewhere in your system, e.g. when using the :doc:`pair_style
hybrid <pair_hybrid>` command. Thus LAMMPS ignores special_bonds
settings when many-body potentials are calculated. Please note,
that the existence of explicit bonds for atoms that are described
by a many-body potential will alter the neighbor list and thus can
render the computation of those interactions invalid, since those
pairs are not only used to determine direct pairwise interactions
but also neighbors of neighbors and more. The recommended course
of action is to remove such bonds, or - if that is not possible -
use a special bonds setting of 1.0 1.0 1.0.
.. note::
Unlike some commands in LAMMPS, you cannot use this command
multiple times in an incremental fashion: e.g. to first set the LJ
settings and then the Coulombic ones. Each time you use this command
it sets all the coefficients to default values and only overrides the
one you specify, so you should set all the options you need each time
you use it. See more details at the bottom of this page.
settings and then the Coulombic ones. Each time you use this
command it sets all the coefficients to default values and only
overrides the one you specify, so you should set all the options
you need each time you use it. See more details at the bottom of
this page.
The Coulomb factors are applied to any Coulomb (charge interaction)
term that the potential calculates. The LJ factors are applied to the
@ -94,14 +98,14 @@ exclude it completely.
The first of the 3 coefficients (LJ or Coulombic) is the weighting
factor on 1-2 atom pairs, which are pairs of atoms directly bonded to
each other. The second coefficient is the weighting factor on 1-3 atom
pairs which are those separated by 2 bonds (e.g. the two H atoms in a
water molecule). The third coefficient is the weighting factor on 1-4
atom pairs which are those separated by 3 bonds (e.g. the first and fourth
atoms in a dihedral interaction). Thus if the 1-2 coefficient is set
to 0.0, then the pairwise interaction is effectively turned off for
all pairs of atoms bonded to each other. If it is set to 1.0, then
that interaction will be at full strength.
each other. The second coefficient is the weighting factor on 1-3
atom pairs which are those separated by 2 bonds (e.g. the two H atoms
in a water molecule). The third coefficient is the weighting factor
on 1-4 atom pairs which are those separated by 3 bonds (e.g. the first
and fourth atoms in a dihedral interaction). Thus if the 1-2
coefficient is set to 0.0, then the pairwise interaction is
effectively turned off for all pairs of atoms bonded to each other.
If it is set to 1.0, then that interaction will be at full strength.
.. note::
@ -184,21 +188,27 @@ interaction between atoms 2 and 5 will be unaffected (full weighting
of 1.0). If the *dihedral* keyword is specified as *no* which is the
default, then the 2,5 interaction will also be weighted by 0.5.
The *one/five* keyword enable calculation and storage of a list of 1-5
neighbors in the molecular topology for each atom. It is required by
some pair styles, such as :doc:`pair_style amoeba <pair_style>` and
:doc:`pair_style hippo <pair_style>`.
----------
.. note::
LAMMPS stores and maintains a data structure with a list of the
first, second, and third neighbors of each atom (within the bond topology of
the system). If new bonds are created (or molecules added containing
atoms with more special neighbors), the size of this list needs to
grow. Note that adding a single bond always adds a new first neighbor
but may also induce \*many\* new second and third neighbors, depending on the
molecular topology of your system. Using the *extra/special/per/atom*
keyword to either :doc:`read_data <read_data>` or :doc:`create_box <create_box>`
reserves empty space in the list for this N additional first, second, or third
neighbors to be added. If you do not do this, you may get an error
when bonds (or molecules) are added.
first, second, and third neighbors of each atom (within the bond
topology of the system). If new bonds are created (or molecules
added containing atoms with more special neighbors), the size of
this list needs to grow. Note that adding a single bond always
adds a new first neighbor but may also induce \*many\* new second
and third neighbors, depending on the molecular topology of your
system. Using the *extra/special/per/atom* keyword to either
:doc:`read_data <read_data>` or :doc:`create_box <create_box>`
reserves empty space in the list for this N additional first,
second, or third neighbors to be added. If you do not do this, you
may get an error when bonds (or molecules) are added.
----------
@ -226,16 +236,16 @@ In the first case you end up with (after the second command):
LJ: 0.0 0.0 0.0
Coul: 0.0 0.0 1.0
while only in the second case, you get the desired settings of:
while only in the second case do you get the desired settings of:
.. parsed-literal::
LJ: 0.0 1.0 1.0
Coul: 0.0 0.0 1.0
This happens because the LJ (and Coul) settings are reset to
their default values before modifying them, each time the
*special_bonds* command is issued.
This happens because the LJ (and Coul) settings are reset to their
default values before modifying them, each time the *special_bonds*
command is issued.
Restrictions
""""""""""""

3
src/AMOEBA/amoeba_file.cpp
View File

@ -493,9 +493,6 @@ void PairAmoeba::read_keyfile(char *filename)
} else if (strcmp(keyword,"pcg-peek") == 0) {
if (nwords != 2) error->all(FLERR,"AMOEBA keyfile line is invalid");
pcgpeek = utils::numeric(FLERR,words[1],true,lmp);
} else if (strcmp(keyword,"usolve-diag") == 0) {
if (nwords != 2) error->all(FLERR,"AMOEBA keyfile line is invalid");
udiag = utils::numeric(FLERR,words[1],true,lmp);
} else {}

8
src/AMOEBA/pair_amoeba.cpp
View File

@ -620,14 +620,6 @@ void PairAmoeba::coeff(int narg, char **arg)
if (narg == 3) read_keyfile(NULL);
else read_keyfile(arg[3]);
// required for now, until 1-5 settings are implemented in LAMMPS
//if (special_hal[4] != 1.0 || special_repel[4] != 1.0 ||
// special_disp[4] != 1.0 || special_mpole[4] != 1.0 ||
// special_polar_pscale[4] != 1.0 || special_polar_piscale[4] != 1.0 ||
// special_polar_wscale[4] != 1.0)
// error->all(FLERR,"AMOEBA 1-5 weights must be 1.0 for now");
// compute Vdwl mixing rules, only for AMOEBA
if (amoeba) {