Merge branch 'develop' of github.com:lammps/lammps into kk_update_4.0

This commit is contained in:
Stan Gerald Moore
2023-07-03 16:35:57 -06:00
279 changed files with 15910 additions and 8704 deletions

View File

@ -52,7 +52,7 @@ can be translated to different output format using the `Sphinx
incorporates programmer documentation extracted from the LAMMPS C++
sources through the `Doxygen <https://doxygen.nl/>`_ program. Currently
the translation to HTML, PDF (via LaTeX), ePUB (for many e-book readers)
and MOBI (for Amazon Kindle(tm) readers) are supported. For that to work a
and MOBI (for Amazon Kindle readers) are supported. For that to work a
Python interpreter version 3.8 or later, the ``doxygen`` tools and
internet access to download additional files and tools are required.
This download is usually only required once or after the documentation

View File

@ -153,11 +153,11 @@ KOKKOS, o = OPENMP, t = OPT.
* :doc:`sph/t/atom <compute_sph_t_atom>`
* :doc:`spin <compute_spin>`
* :doc:`stress/atom <compute_stress_atom>`
* :doc:`stress/cartesian <compute_stress_profile>`
* :doc:`stress/cylinder <compute_stress_profile>`
* :doc:`stress/cartesian <compute_stress_cartesian>`
* :doc:`stress/cylinder <compute_stress_curvilinear>`
* :doc:`stress/mop <compute_stress_mop>`
* :doc:`stress/mop/profile <compute_stress_mop>`
* :doc:`stress/spherical <compute_stress_profile>`
* :doc:`stress/spherical <compute_stress_curvilinear>`
* :doc:`stress/tally <compute_tally>`
* :doc:`tdpd/cc/atom <compute_tdpd_cc_atom>`
* :doc:`temp (k) <compute_temp>`

View File

@ -171,6 +171,7 @@ OPT.
* :doc:`pafi <fix_pafi>`
* :doc:`pair <fix_pair>`
* :doc:`phonon <fix_phonon>`
* :doc:`pimd/langevin <fix_pimd>`
* :doc:`pimd/nvt <fix_pimd>`
* :doc:`planeforce <fix_planeforce>`
* :doc:`plumed <fix_plumed>`

View File

@ -37,6 +37,7 @@ OPT.
*
* :doc:`adp (ko) <pair_adp>`
* :doc:`agni (o) <pair_agni>`
* :doc:`aip/water/2dm (t) <pair_aip_water_2dm>`
* :doc:`airebo (io) <pair_airebo>`
* :doc:`airebo/morse (io) <pair_airebo>`
* :doc:`amoeba (g) <pair_amoeba>`

View File

@ -41,7 +41,7 @@ warning.
LATTE package
-------------
.. deprecated:: TBD
.. deprecated:: 15Jun2023
The LATTE package with the fix latte command was removed from LAMMPS.
This functionality has been superseded by :doc:`fix mdi/qm <fix_mdi_qm>`

View File

@ -56,17 +56,6 @@ C++ in the ``examples/COUPLE/simple`` folder of the LAMMPS distribution.
and Ubuntu 18.04 LTS and not compatible. Either newer compilers
need to be installed or the Linux updated.
.. versionchanged:: 8Feb2023
.. note::
A contributed Fortran interface is available in the
``examples/COUPLE/fortran2`` folder. However, since the completion
of the :f:mod:`LIBLAMMPS` module, this interface is now deprecated,
no longer actively maintained and will likely be removed in the
future. Please see the ``README`` file in that folder for more
information about it and how to contact its author and maintainer.
----------
Creating or deleting a LAMMPS object
@ -203,40 +192,62 @@ Below is an example demonstrating some of the possible uses.
.. code-block:: fortran
PROGRAM testprop
USE LIBLAMMPS
USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int64_t
USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : OUTPUT_UNIT
TYPE(lammps) :: lmp
INTEGER(KIND=c_int64_t), POINTER :: natoms
REAL(KIND=c_double), POINTER :: dt
INTEGER(KIND=c_int64_t), POINTER :: ntimestep
REAL(KIND=c_double) :: pe, ke
PROGRAM testprop
USE LIBLAMMPS
USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int64_t, c_int
USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : OUTPUT_UNIT
TYPE(lammps) :: lmp
INTEGER(KIND=c_int64_t), POINTER :: natoms, ntimestep, bval
REAL(KIND=c_double), POINTER :: dt, dval
INTEGER(KIND=c_int), POINTER :: nfield, typ, ival
INTEGER(KIND=c_int) :: i
CHARACTER(LEN=11) :: key
REAL(KIND=c_double) :: pe, ke
lmp = lammps()
CALL lmp%file('in.sysinit')
natoms = lmp%extract_global('natoms')
WRITE(OUTPUT_UNIT,'(A,I0,A)') 'Running a simulation with ', natoms, ' atoms'
WRITE(OUTPUT_UNIT,'(I0,A,I0,A,I0,A)') lmp%extract_setting('nlocal'), &
' local and ', lmp%extract_setting('nghost'), ' ghost atoms. ', &
lmp%extract_setting('ntypes'), ' atom types'
lmp = lammps()
CALL lmp%file('in.sysinit')
natoms = lmp%extract_global('natoms')
WRITE(OUTPUT_UNIT,'(A,I0,A)') 'Running a simulation with ', natoms, ' atoms'
WRITE(OUTPUT_UNIT,'(I0,A,I0,A,I0,A)') lmp%extract_setting('nlocal'), &
' local and ', lmp%extract_setting('nghost'), ' ghost atoms. ', &
lmp%extract_setting('ntypes'), ' atom types'
CALL lmp%command('run 2 post no')
dt = lmp%extract_global('dt')
ntimestep = lmp%extract_global('ntimestep')
WRITE(OUTPUT_UNIT,'(A,I0,A,F4.1,A)') 'At step: ', ntimestep, &
' Changing timestep from', dt, ' to 0.5'
dt = 0.5_c_double
CALL lmp%command('run 2 post no')
CALL lmp%command('run 2 post no')
WRITE(OUTPUT_UNIT,'(A,I0)') 'At step: ', ntimestep
pe = lmp%get_thermo('pe')
ke = lmp%get_thermo('ke')
PRINT*, 'PE = ', pe
PRINT*, 'KE = ', ke
ntimestep = lmp%last_thermo('step', 0)
nfield = lmp%last_thermo('num', 0)
WRITE(OUTPUT_UNIT,'(A,I0,A,I0)') 'Last thermo output on step: ', ntimestep, &
', number of fields: ', nfield
DO i=1, nfield
key = lmp%last_thermo('keyword',i)
typ = lmp%last_thermo('type',i)
IF (typ == lmp%dtype%i32) THEN
ival = lmp%last_thermo('data',i)
WRITE(OUTPUT_UNIT,*) key, ':', ival
ELSE IF (typ == lmp%dtype%i64) THEN
bval = lmp%last_thermo('data',i)
WRITE(OUTPUT_UNIT,*) key, ':', bval
ELSE IF (typ == lmp%dtype%r64) THEN
dval = lmp%last_thermo('data',i)
WRITE(OUTPUT_UNIT,*) key, ':', dval
END IF
END DO
CALL lmp%close(.TRUE.)
END PROGRAM testprop
dt = lmp%extract_global('dt')
ntimestep = lmp%extract_global('ntimestep')
WRITE(OUTPUT_UNIT,'(A,I0,A,F4.1,A)') 'At step: ', ntimestep, &
' Changing timestep from', dt, ' to 0.5'
dt = 0.5_c_double
CALL lmp%command('run 2 post no')
WRITE(OUTPUT_UNIT,'(A,I0)') 'At step: ', ntimestep
pe = lmp%get_thermo('pe')
ke = lmp%get_thermo('ke')
WRITE(OUTPUT_UNIT,*) 'PE = ', pe
WRITE(OUTPUT_UNIT,*) 'KE = ', ke
CALL lmp%close(.TRUE.)
END PROGRAM testprop
---------------
@ -262,6 +273,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS.
:ftype style: type(lammps_style)
:f type: derived type to access lammps type constants
:ftype type: type(lammps_type)
:f dtype: derived type to access lammps data type constants
:ftype dtype: type(lammps_dtype)
:f close: :f:subr:`close`
:ftype close: subroutine
:f subroutine error: :f:subr:`error`
@ -278,6 +291,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS.
:ftype get_natoms: function
:f get_thermo: :f:func:`get_thermo`
:ftype get_thermo: function
:f last_thermo: :f:func:`last_thermo`
:ftype last_thermo: function
:f extract_box: :f:subr:`extract_box`
:ftype extract_box: subroutine
:f reset_box: :f:subr:`reset_box`
@ -587,6 +602,96 @@ Procedures Bound to the :f:type:`lammps` Derived Type
--------
.. f:function:: last_thermo(what, index)
This function will call :cpp:func:`lammps_last_thermo` and returns
either a string or a pointer to a cached copy of LAMMPS last thermodynamic
output, depending on the data requested through *what*. Note that *index*
uses 1-based indexing to access thermo output columns.
.. versionadded:: 15Jun2023
Note that this function actually does not return a value, but rather
associates the pointer on the left side of the assignment to point to
internal LAMMPS data (with the exception of string data, which are
copied and returned as ordinary Fortran strings). Pointers must be
of the correct data type to point to said data (typically
``INTEGER(c_int)``, ``INTEGER(c_int64_t)``, or ``REAL(c_double)``).
The pointer being associated with LAMMPS data is type-checked at
run-time via an overloaded assignment operator. The pointers
returned by this function point to temporary, read-only data that may
be overwritten at any time, so their target values need to be copied
to local storage if they are supposed to persist.
For example,
.. code-block:: fortran
PROGRAM thermo
USE LIBLAMMPS
USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int64_t, c_int
TYPE(lammps) :: lmp
INTEGER(KIND=c_int64_t), POINTER :: ntimestep, bval
REAL(KIND=c_double), POINTER :: dval
INTEGER(KIND=c_int), POINTER :: nfield, typ, ival
INTEGER(KIND=c_int) :: i
CHARACTER(LEN=11) :: key
lmp = lammps()
CALL lmp%file('in.sysinit')
ntimestep = lmp%last_thermo('step', 0)
nfield = lmp%last_thermo('num', 0)
PRINT*, 'Last thermo output on step: ', ntimestep, ' Number of fields: ', nfield
DO i=1, nfield
key = lmp%last_thermo('keyword',i)
typ = lmp%last_thermo('type',i)
IF (typ == lmp%dtype%i32) THEN
ival = lmp%last_thermo('data',i)
PRINT*, key, ':', ival
ELSE IF (typ == lmp%dtype%i64) THEN
bval = lmp%last_thermo('data',i)
PRINT*, key, ':', bval
ELSE IF (typ == lmp%dtype%r64) THEN
dval = lmp%last_thermo('data',i)
PRINT*, key, ':', dval
END IF
END DO
CALL lmp%close(.TRUE.)
END PROGRAM thermo
would extract the last timestep where thermo output was done and the number
of columns it printed. Then it loops over the columns to print out column
header keywords and the corresponding data.
.. note::
If :f:func:`last_thermo` returns a string, the string must have a length
greater than or equal to the length of the string (not including the
terminal ``NULL`` character) that LAMMPS returns. If the variable's
length is too short, the string will be truncated. As usual in Fortran,
strings are padded with spaces at the end. If you use an allocatable
string, the string **must be allocated** prior to calling this function.
:p character(len=\*) what: string with the name of the thermo keyword
:p integer(c_int) index: 1-based column index
:to: :cpp:func:`lammps_last_thermo`
:r pointer [polymorphic]: pointer to LAMMPS data. The left-hand side of the
assignment should be either a string (if expecting string data) or a
C-compatible pointer (e.g., ``INTEGER(c_int), POINTER :: nlocal``) to the
extracted property.
.. warning::
Modifying the data in the location pointed to by the returned pointer
may lead to inconsistent internal data and thus may cause failures,
crashes, or bogus simulations. In general, it is much better
to use a LAMMPS input command that sets or changes these parameters.
Using an input command will take care of all side effects and necessary
updates of settings derived from such settings.
--------
.. f:subroutine:: extract_box([boxlo][, boxhi][, xy][, yz][, xz][, pflags][, boxflag])
This subroutine will call :cpp:func:`lammps_extract_box`. All
@ -764,13 +869,14 @@ Procedures Bound to the :f:type:`lammps` Derived Type
.. note::
If :f:func:`extract_global` returns a string, the string must have length
greater than or equal to the length of the string (not including the
terminal ``NULL`` character) that LAMMPS returns. If the variable's
length is too short, the string will be truncated. As usual in Fortran,
strings are padded with spaces at the end. If you use an allocatable
string, the string **must be allocated** prior to calling this function,
but you can automatically reallocate it to the correct length after the
If :f:func:`extract_global` returns a string, the string must have
a length greater than or equal to the length of the string (not
including the terminal ``NULL`` character) that LAMMPS returns. If
the variable's length is too short, the string will be
truncated. As usual in Fortran, strings are padded with spaces at
the end. If you use an allocatable string, the string **must be
allocated** prior to calling this function, but you can
automatically reallocate it to the correct length after the
function returns, viz.,
.. code-block :: fortran

View File

@ -69,15 +69,13 @@ SPC/E with rigid bonds.
timestep 1.0
fix rigid all shake 0.0001 10 10000 b 1 a 1
minimize 0.0 0.0 1000 10000
run 0 post no
reset_timestep 0
velocity all create 300.0 5463576
fix integrate all nvt temp 300.0 300.0 1.0
fix integrate all nvt temp 300.0 300.0 100.0
thermo_style custom step temp press etotal density pe ke
thermo 1000
run 20000 upto
write_data tip4p.data nocoeff
write_data spce.data nocoeff
.. _spce_molecule:
.. code-block::

View File

@ -128,11 +128,11 @@ TIP3P with rigid bonds.
fix rigid all shake 0.001 10 10000 b 1 a 1
minimize 0.0 0.0 1000 10000
run 0 post no
reset_timestep 0
timestep 1.0
velocity all create 300.0 5463576
fix integrate all nvt temp 300 300 1.0
fix integrate all nvt temp 300 300 100.0
thermo_style custom step temp press etotal pe

View File

@ -180,17 +180,17 @@ file changed):
fix rigid all shake 0.001 10 10000 b 1 a 1
minimize 0.0 0.0 1000 10000
run 0 post no
reset_timestep 0
timestep 1.0
velocity all create 300.0 5463576
fix integrate all nvt temp 300 300 1.0
fix integrate all nvt temp 300 300 100.0
thermo_style custom step temp press etotal pe
thermo 1000
run 20000
write_data tip3p.data nocoeff
write_data tip4p-implicit.data nocoeff
Below is the code for a LAMMPS input file using the explicit method and
a TIP4P molecule file. Because of using :doc:`fix rigid/nvt/small
@ -203,6 +203,7 @@ rigid/nvt/small can identify rigid bodies by their molecule ID:
units real
atom_style charge
atom_modify map array
region box block -5 5 -5 5 -5 5
create_box 3 box
@ -219,14 +220,14 @@ rigid/nvt/small can identify rigid bodies by their molecule ID:
molecule water tip4p.mol
create_atoms 0 random 33 34564 NULL mol water 25367 overlap 1.33
timestep 0.1
fix integrate all rigid/nvt/small molecule temp 300.0 300.0 1.0
timestep 0.5
fix integrate all rigid/nvt/small molecule temp 300.0 300.0 100.0
velocity all create 300.0 5463576
thermo_style custom step temp press etotal density pe ke
thermo 1000
run 20000
write_data tip4p.data nocoeff
write_data tip4p-explicit.data nocoeff
.. _tip4p_molecule:
.. code-block::

View File

@ -91,6 +91,7 @@ ID:
units real
atom_style charge
atom_modify map array
region box block -5 5 -5 5 -5 5
create_box 3 box
@ -107,8 +108,8 @@ ID:
molecule water tip5p.mol
create_atoms 0 random 33 34564 NULL mol water 25367 overlap 1.33
timestep 0.20
fix integrate all rigid/nvt/small molecule temp 300.0 300.0 1.0
timestep 0.5
fix integrate all rigid/nvt/small molecule temp 300.0 300.0 100.0
reset_timestep 0
velocity all create 300.0 5463576

View File

@ -5,6 +5,7 @@ This section documents the following functions:
- :cpp:func:`lammps_get_natoms`
- :cpp:func:`lammps_get_thermo`
- :cpp:func:`lammps_last_thermo`
- :cpp:func:`lammps_extract_box`
- :cpp:func:`lammps_reset_box`
- :cpp:func:`lammps_memory_usage`
@ -81,6 +82,11 @@ subdomains and processors.
-----------------------
.. doxygenfunction:: lammps_last_thermo
:project: progguide
-----------------------
.. doxygenfunction:: lammps_extract_box
:project: progguide

View File

@ -53,6 +53,7 @@ against invalid accesses.
* :py:meth:`version() <lammps.lammps.version()>`: return the numerical version id, e.g. LAMMPS 2 Sep 2015 -> 20150902
* :py:meth:`get_thermo() <lammps.lammps.get_thermo()>`: return current value of a thermo keyword
* :py:meth:`last_thermo() <lammps.lammps.last_thermo()>`: return a dictionary of the last thermodynamic output
* :py:meth:`get_natoms() <lammps.lammps.get_natoms()>`: total # of atoms as int
* :py:meth:`reset_box() <lammps.lammps.reset_box()>`: reset the simulation box size
* :py:meth:`extract_setting() <lammps.lammps.extract_setting()>`: return a global setting
@ -60,6 +61,10 @@ against invalid accesses.
* :py:meth:`extract_box() <lammps.lammps.extract_box()>`: extract box info
* :py:meth:`create_atoms() <lammps.lammps.create_atoms()>`: create N atoms with IDs, types, x, v, and image flags
**Properties**:
* :py:attr:`last_thermo_step <lammps.lammps.last_thermo_step>`: the last timestep thermodynamic output was computed
.. tab:: PyLammps/IPyLammps API
In addition to the functions provided by :py:class:`lammps <lammps.lammps>`, :py:class:`PyLammps <lammps.PyLammps>` objects

View File

@ -24,14 +24,17 @@ Syntax
*x, y, z* = the center of mass position of the 2 atoms when the bond broke (distance units)
*x/ref, y/ref, z/ref* = the initial center of mass position of the 2 atoms (distance units)
*overlay/pair* value = none
*overlay/pair* value = *yes* or *no*
bonded particles will still interact with pair forces
*smooth* value = *yes* or *no*
smooths bond forces near the breaking point
*break/no*
indicates that bonds should not break during a run
*normalize* value = *yes* or *no*
normalizes normal and shear forces by the reference length
*break* value = *yes* or *no*
indicates whether bonds break during a run
Examples
""""""""
@ -136,16 +139,19 @@ or :doc:`read_restart <read_restart>` commands:
* :math:`\gamma_r` (force*distance/velocity units)
* :math:`\gamma_t` (force*distance/velocity units)
However, the *normalize* option will normalize the radial and shear forces
by :math:`r_0` such that :math:`k_r` and :math:`k_s` are unit less.
By default, pair forces are not calculated between bonded particles.
Pair forces can alternatively be overlaid on top of bond forces using
the *overlay/pair* keyword. These settings require specific
the *overlay/pair* option. These settings require specific
:doc:`special_bonds <special_bonds>` settings described in the
restrictions. Further details can be found in the `:doc: how to
<Howto_BPM>` page on BPMs.
.. versionadded:: 28Mar2023
If the *break/no* keyword is used, then LAMMPS assumes bonds should not break
If the *break* option is used, then LAMMPS assumes bonds should not break
during a simulation run. This will prevent some unnecessary calculation.
However, if a bond does break, it will trigger an error.
@ -251,7 +257,7 @@ Related commands
Default
"""""""
The option defaults are *smooth* = *yes*
The option defaults are *overlay/pair* = *no*, *smooth* = *yes*, *normalize* = *no*, and *break* = *yes*
----------

View File

@ -24,14 +24,17 @@ Syntax
*x, y, z* = the center of mass position of the 2 atoms when the bond broke (distance units)
*x/ref, y/ref, z/ref* = the initial center of mass position of the 2 atoms (distance units)
*overlay/pair* value = none
*overlay/pair* value = *yes* or *no*
bonded particles will still interact with pair forces
*smooth* value = *yes* or *no*
smooths bond forces near the breaking point
*break/no*
indicates that bonds should not break during a run
*normalize* value = *yes* or *no*
normalizes bond forces by the reference length
*break* value = *yes* or *no*
indicates whether bonds break during a run
Examples
""""""""
@ -66,7 +69,7 @@ particles based on a model described by Clemmer and Robbins
F = k (r - r_0) w
where :math:`k_r` is a stiffness, :math:`r` is the current distance
where :math:`k` is a stiffness, :math:`r` is the current distance
and :math:`r_0` is the initial distance between the two particles, and
:math:`w` is an optional smoothing factor discussed below. Bonds will
break at a strain of :math:`\epsilon_c`. This is done by setting by
@ -102,16 +105,19 @@ the data file or restart files read by the :doc:`read_data
* :math:`\epsilon_c` (unit less)
* :math:`\gamma` (force/velocity units)
However, the *normalize* option will normalize the elastic bond force by
:math:`r_0` such that :math:`k` is unit less.
By default, pair forces are not calculated between bonded particles.
Pair forces can alternatively be overlaid on top of bond forces using
the *overlay/pair* keyword. These settings require specific
the *overlay/pair* option. These settings require specific
:doc:`special_bonds <special_bonds>` settings described in the
restrictions. Further details can be found in the `:doc: how to
<Howto_BPM>` page on BPMs.
.. versionadded:: 28Mar2023
If the *break/no* keyword is used, then LAMMPS assumes bonds should not break
If the *break* option is used, then LAMMPS assumes bonds should not break
during a simulation run. This will prevent some unnecessary calculation.
However, if a bond does break, it will trigger an error.
@ -206,7 +212,7 @@ Related commands
Default
"""""""
The option defaults are *smooth* = *yes*
The option defaults are *overlay/pair* = *no*, *smooth* = *yes*, *normalize* = *no*, and *break* = *yes*
----------

View File

@ -307,11 +307,11 @@ The individual style names on the :doc:`Commands compute <Commands_compute>` pag
* :doc:`sph/t/atom <compute_sph_t_atom>` - per-atom internal temperature of Smooth-Particle Hydrodynamics atoms
* :doc:`spin <compute_spin>` - magnetic quantities for a system of atoms having spins
* :doc:`stress/atom <compute_stress_atom>` - stress tensor for each atom
* :doc:`stress/cartesian <compute_stress_profile>` - stress tensor in cartesian coordinates
* :doc:`stress/cylinder <compute_stress_profile>` - stress tensor in cylindrical coordinates
* :doc:`stress/cartesian <compute_stress_cartesian>` - stress tensor in cartesian coordinates
* :doc:`stress/cylinder <compute_stress_curvilinear>` - stress tensor in cylindrical coordinates
* :doc:`stress/mop <compute_stress_mop>` - normal components of the local stress tensor using the method of planes
* :doc:`stress/mop/profile <compute_stress_mop>` - profile of the normal components of the local stress tensor using the method of planes
* :doc:`stress/spherical <compute_stress_profile>` - stress tensor in spherical coordinates
* :doc:`stress/spherical <compute_stress_curvilinear>` - stress tensor in spherical coordinates
* :doc:`stress/tally <compute_tally>` - stress between two groups of atoms via the tally callback mechanism
* :doc:`tdpd/cc/atom <compute_tdpd_cc_atom>` - per-atom chemical concentration of a specified species for each tDPD particle
* :doc:`temp <compute_temp>` - temperature of group of atoms

View File

@ -76,7 +76,10 @@ The value *force* is the magnitude of the force acting between the
pair of atoms in the bond.
The values *fx*, *fy*, and *fz* are the xyz components of
*force* between the pair of atoms in the bond.
*force* between the pair of atoms in the bond. For bond styles that apply
non-central forces, such as :doc:`bond_style bpm/rotational
<bond_bpm_rotational>`, these values only include the :math:`(x,y,z)`
components of the normal force component.
The remaining properties are all computed for motion of the two atoms
relative to the center of mass (COM) velocity of the 2 atoms in the

View File

@ -25,7 +25,7 @@ Examples
Description
"""""""""""
.. versionadded:: TBD
.. versionadded:: 15Jun2023
Define a computation that counts the current number of atoms for each
atom type. Or the number of bonds (angles, dihedrals, impropers) for

View File

@ -146,13 +146,13 @@ m to :math:`M` (inclusive). A middle asterisk means all types from m to n
Output info
"""""""""""
This compute calculates a local vector of doubles and a scalar. The vector
stores the unique components of the first requested tensor in the order
:math:`xx`, :math:`yy`, :math:`zz`, :math:`xy`, :math:`xz`, :math:`yz`
followed by the same components for all subsequent tensors.
This compute calculates a global vector of doubles and a global scalar. The
vector stores the unique components of the first requested tensor in the
order :math:`xx`, :math:`yy`, :math:`zz`, :math:`xy`, :math:`xz`,
:math:`yz` followed by the same components for all subsequent tensors.
The length of the vector is therefore six times the number of requested
tensors. The scalar output is the number of pairwise interactions included in
the calculation of the fabric tensor.
tensors. The scalar output is the number of pairwise interactions included
in the calculation of the fabric tensor.
Restrictions
""""""""""""

View File

@ -66,7 +66,9 @@ The value *eng* is the interaction energy for the pair of atoms.
The value *force* is the force acting between the pair of atoms, which
is positive for a repulsive force and negative for an attractive
force. The values *fx*, *fy*, and *fz* are the :math:`(x,y,z)` components of
*force* on atom I.
*force* on atom I. For pair styles that apply non-central forces,
such as :doc:`granular pair styles <pair_gran>`, these values only include
the :math:`(x,y,z)` components of the normal force component.
A pair style may define additional pairwise quantities which can be
accessed as *p1* to *pN*, where :math:`N` is defined by the pair style.

View File

@ -0,0 +1,123 @@
.. index:: compute stress/cartesian
compute stress/cartesian command
==================================
Syntax
""""""
.. code-block:: LAMMPS
compute ID group-ID stress/cartesian args
* ID, group-ID are documented in :doc:`compute <compute>` command
* args = argument specific to the compute style
.. parsed-literal::
*stress/cartesian* args = dim1 bin_width1 dim2 bin_width2 keyword
dim1 = *x* or *y* or *z*
bin_width1 = width of the bin
dim2 = *x* or *y* or *z* or *NULL*
bin_width2 = width of the bin
keyword = *ke* or *pair* or *bond*
Examples
""""""""
.. code-block:: LAMMPS
compute 1 all stress/cartesian x 0.1 NULL 0
compute 1 all stress/cartesian y 0.1 z 0.1
compute 1 all stress/cartesian x 0.1 NULL 0 ke pair
Description
"""""""""""
Compute *stress/cartesian* defines computations that calculate profiles of the
diagonal components of the local stress tensor over one or two Cartesian
dimensions, as described in :ref:`(Ikeshoji)<Ikeshoji2>`. The stress tensor is
split into a kinetic contribution :math:`P^k` and a virial contribution
:math:`P^v`. The sum gives the total stress tensor :math:`P = P^k+P^v`.
This compute obeys momentum balance through fluid interfaces. They use the
Irving--Kirkwood contour, which is the straight line between particle pairs.
.. versionadded:: 15Jun2023
Added support for bond styles
This compute only supports pair and bond (no angle, dihedral, improper,
or kspace) forces. By default, if no extra keywords are specified, all
supported contributions to the stress are included (ke, pair, bond). If any
keywords are specified, then only those components are summed.
Output info
"""""""""""
The output columns for *stress/cartesian* are the position of the
center of the local volume in the first and second dimensions, number
density, :math:`P^k_{xx}`, :math:`P^k_{yy}`, :math:`P^k_{zz}`,
:math:`P^v_{xx}`, :math:`P^v_{yy}`, and :math:`P^v_{zz}`. There are 8
columns when one dimension is specified and 9 columns when two
dimensions are specified. The number of bins (rows) is
:math:`(L_1/b_1)(L_2/b_2)`, where :math:`L_1` and :math:`L_2` are the lengths
of the simulation box in the specified dimensions and :math:`b_1` and
:math:`b_2` are the specified bin widths. When only one dimension is
specified, the number of bins (rows) is :math:`L_1/b_1`.
This array can be output with :doc:`fix ave/time <fix_ave_time>`,
.. code-block:: LAMMPS
compute p all stress/cartesian x 0.1
fix 2 all ave/time 100 1 100 c_p[*] file dump_p.out mode vector
The values calculated by this compute are "intensive". The stress
values will be in pressure :doc:`units <units>`. The number density
values are in inverse volume :doc:`units <units>`.
NOTE 1: The local stress does not include any Lennard-Jones tail
corrections to the stress added by the
:doc:`pair_modify tail yes <pair_modify>`
command, since those are contributions to the global system pressure.
NOTE 2: The local stress profiles generated by these computes are
similar to those obtained by the
:doc:`method-of-planes (MOP) <compute_stress_mop>`.
A key difference is that compute
:doc:`stress/mop/profile <compute_stress_mop>`
considers particles crossing a set of planes, while
*stress/cartesian* computes averages for a set of small volumes.
Moreover, this compute computes the diagonal components of the stress
tensor :math:`P_{xx}`, :math:`P_{yy}`, and :math:`P_{zz}`, while
*stress/mop/profile* computes the components
:math:`P_{ix}`, :math:`P_{iy}`, and :math:`P_{iz}`, where :math:`i` is the
direction normal to the plane.
More information on the similarities and differences can be found in
:ref:`(Ikeshoji)<Ikeshoji2>`.
Restrictions
""""""""""""
These computes calculate the stress tensor contributions for pair and bond
forces only (no angle, dihedral, improper, or kspace force).
It requires pairwise force calculations not available for most
many-body pair styles.
These computes are part of the EXTRA-COMPUTE package. They are only
enabled if LAMMPS was built with that package. See the :doc:`Build
package <Build_package>` doc page for more info.
Related commands
""""""""""""""""
:doc:`compute stress/atom <compute_stress_atom>`, :doc:`compute pressure <compute_pressure>`,
:doc:`compute stress/mop/profile <compute_stress_mop>`, :doc:`compute stress/spherical <compute_stress_curvilinear>`,
:doc:`compute stress/cylinder <compute_stress_curvilinear>`
----------
.. _Ikeshoji2:
**(Ikeshoji)** Ikeshoji, Hafskjold, Furuholt, Mol Sim, 29, 101-109, (2003).

View File

@ -1,11 +1,6 @@
.. index:: compute stress/cartesian
.. index:: compute stress/cylinder
.. index:: compute stress/spherical
compute stress/cartesian command
==================================
compute stress/cylinder command
=================================
@ -20,15 +15,11 @@ Syntax
compute ID group-ID style args
* ID, group-ID are documented in :doc:`compute <compute>` command
* style = stress/cartesian or stress/spherical or stress/cylinder
* style = stress/spherical or stress/cylinder
* args = argument specific to the compute style
.. parsed-literal::
*stress/cartesian* args = dim bin_width
dim = *x* or *y* or *z*
bin_width = width of the bin
one or two dim/bin_width pairs may be appended
*stress/cylinder* args = zlo zh Rmax bin_width keyword
zlo = minimum z-boundary for cylinder
zhi = maximum z-boundary for cylinder
@ -46,8 +37,6 @@ Examples
.. code-block:: LAMMPS
compute 1 all stress/cartesian x 0.1
compute 1 all stress/cartesian y 0.25 z 0.1
compute 1 all stress/cylinder -10.0 10.0 15.0 0.25
compute 1 all stress/cylinder -10.0 10.0 15.0 0.25 ke no
compute 1 all stress/spherical 0 0 0 0.1 10
@ -55,41 +44,28 @@ Examples
Description
"""""""""""
Compute *stress/cartesian*, compute *stress/cylinder*, and compute
Compute *stress/cylinder*, and compute
*stress/spherical* define computations that calculate profiles of the
diagonal components of the local stress tensor in the specified
coordinate system. The stress tensor is split into a kinetic
contribution :math:`P^k` and a virial contribution :math:`P^v`. The sum
gives the total stress tensor :math:`P = P^k+P^v`. These computes can
for example be used to calculate the diagonal components of the local
stress tensor of interfaces with flat, cylindrical, or spherical
stress tensor of surfaces with cylindrical or spherical
symmetry. These computes obeys momentum balance through fluid
interfaces. They use the Irving--Kirkwood contour, which is the straight
line between particle pairs.
The *stress/cartesian* computes the stress profile along one or two
Cartesian coordinates, as described in :ref:`(Ikeshoji)<Ikeshoji2>`. The
compute *stress/cylinder* computes the stress profile along the
The compute *stress/cylinder* computes the stress profile along the
radial direction in cylindrical coordinates, as described in
:ref:`(Addington)<Addington1>`. The compute *stress/spherical*
computes the stress profile along the radial direction in spherical
coordinates, as described in :ref:`(Ikeshoji)<Ikeshoji2>`.
coordinates, as described in :ref:`(Ikeshoji)<Ikeshoji4>`.
Output info
"""""""""""
The output columns for *stress/cartesian* are the position of the
center of the local volume in the first and second dimensions, number
density, :math:`P^k_{xx}`, :math:`P^k_{yy}`, :math:`P^k_{zz}`,
:math:`P^v_{xx}`, :math:`P^v_{yy}`, and :math:`P^v_{zz}`. There are 8
columns when one dimension is specified and 9 columns when two
dimensions are specified. The number of bins (rows) is
:math:`(L_1/b_1)(L_2/b_2)`, where :math:`L_1` and :math:`L_2` are the lengths
of the simulation box in the specified dimensions and :math:`b_1` and
:math:`b_2` are the specified bin widths. When only one dimension is
specified, the number of bins (rows) is :math:`L_1/b_1`.
The default output columns for *stress/cylinder* are the radius to the
center of the cylindrical shell, number density, :math:`P^k_{rr}`,
:math:`P^k_{\phi\phi}`, :math:`P^k_{zz}`, :math:`P^v_{rr}`,
@ -111,7 +87,7 @@ This array can be output with :doc:`fix ave/time <fix_ave_time>`,
.. code-block:: LAMMPS
compute p all stress/cartesian x 0.1
compute 1 all stress/spherical 0 0 0 0.1 10
fix 2 all ave/time 100 1 100 c_p[*] file dump_p.out mode vector
The values calculated by this compute are "intensive". The stress
@ -123,16 +99,6 @@ corrections to the stress added by the
:doc:`pair_modify tail yes <pair_modify>`
command, since those are contributions to the global system pressure.
NOTE 2: The local stress profiles generated by these computes are
similar to those obtained by the
:doc:`method-of-planes (MOP) <compute_stress_mop>`.
A key difference
is that compute `stress/mop/profile <compute_stress_mop>`
considers particles crossing a set of planes, while
*stress/cartesian* computes averages for a set of small volumes.
More information on the similarities and differences can be found in
:ref:`(Ikeshoji)<Ikeshoji2>`.
Restrictions
""""""""""""
@ -150,7 +116,8 @@ package <Build_package>` doc page for more info.
Related commands
""""""""""""""""
:doc:`compute stress/atom <compute_stress_atom>`, :doc:`compute pressure <compute_pressure>`, :doc:`compute stress/mop/profile <compute_stress_mop>`
:doc:`compute stress/atom <compute_stress_atom>`, :doc:`compute pressure <compute_pressure>`,
:doc:`compute stress/mop/profile <compute_stress_mop>`, :doc:`compute stress/cartesian <compute_stress_cartesian>`
Default
"""""""
@ -159,7 +126,7 @@ The keyword default for ke in style *stress/cylinder* is yes.
----------
.. _Ikeshoji2:
.. _Ikeshoji4:
**(Ikeshoji)** Ikeshoji, Hafskjold, Furuholt, Mol Sim, 29, 101-109, (2003).

View File

@ -18,7 +18,7 @@ Syntax
* style = *stress/mop* or *stress/mop/profile*
* dir = *x* or *y* or *z* is the direction normal to the plane
* args = argument specific to the compute style
* keywords = *kin* or *conf* or *total* (one of more can be specified)
* keywords = *kin* or *conf* or *total* or *pair* or *bond* or *angle* (one or more can be specified)
.. parsed-literal::
@ -45,85 +45,112 @@ Examples
Description
"""""""""""
Compute *stress/mop* and compute *stress/mop/profile* define computations that
calculate components of the local stress tensor using the method of
planes :ref:`(Todd) <mop-todd>`. Specifically in compute *stress/mop* calculates 3
components are computed in directions *dir*,\ *x*\ ; *dir*,\ *y*\ ; and
*dir*,\ *z*\ ; where *dir* is the direction normal to the plane, while
in compute *stress/mop/profile* the profile of the stress is computed.
Compute *stress/mop* and compute *stress/mop/profile* define
computations that calculate components of the local stress tensor using
the method of planes :ref:`(Todd) <mop-todd>`. Specifically in compute
*stress/mop* calculates 3 components are computed in directions *dir*,\
*x*\ ; *dir*,\ *y*\ ; and *dir*,\ *z*\ ; where *dir* is the direction
normal to the plane, while in compute *stress/mop/profile* the profile
of the stress is computed.
Contrary to methods based on histograms of atomic stress (i.e., using
:doc:`compute stress/atom <compute_stress_atom>`), the method of planes is
compatible with mechanical balance in heterogeneous systems and at
:doc:`compute stress/atom <compute_stress_atom>`), the method of planes
is compatible with mechanical balance in heterogeneous systems and at
interfaces :ref:`(Todd) <mop-todd>`.
The stress tensor is the sum of a kinetic term and a configurational
term, which are given respectively by Eq. (21) and Eq. (16) in
:ref:`(Todd) <mop-todd>`. For the kinetic part, the algorithm considers that
atoms have crossed the plane if their positions at times :math:`t-\Delta t`
and :math:`t` are one on either side of the plane, and uses the velocity at
time :math:`t-\Delta t/2` given by the velocity Verlet algorithm.
:ref:`(Todd) <mop-todd>`. For the kinetic part, the algorithm considers
that atoms have crossed the plane if their positions at times
:math:`t-\Delta t` and :math:`t` are one on either side of the plane,
and uses the velocity at time :math:`t-\Delta t/2` given by the velocity
Verlet algorithm.
Between one and three keywords can be used to indicate which
contributions to the stress must be computed: kinetic stress (kin),
configurational stress (conf), and/or total stress (total).
.. versionadded:: 15Jun2023
NOTE 1: The configurational stress is computed considering all pairs of atoms where at least one atom belongs to group group-ID.
contributions from bond and angle potentials
Between one and six keywords can be used to indicate which contributions
to the stress must be computed: total stress (total), kinetic stress
(kin), configurational stress (conf), stress due to bond stretching
(bond), stress due to angle bending (angle) and/or due to pairwise
non-bonded interactions (pair). The angle keyword is currently
available only for the *stress/mop* command and **not** the
*stress/mop/profile* command.
NOTE 1: The configurational stress is computed considering all pairs of
atoms where at least one atom belongs to group group-ID.
NOTE 2: The local stress does not include any Lennard-Jones tail
corrections to the stress added by the :doc:`pair_modify tail yes <pair_modify>`
command, since those are contributions to the global system pressure.
corrections to the stress added by the :doc:`pair_modify tail yes
<pair_modify>` command, since those are contributions to the global
system pressure.
NOTE 3: The local stress profile generated by compute *stress/mop/profile*
is similar to that obtained by compute
:doc:`stress/cartesian <compute_stress_profile>`.
A key difference is that compute *stress/mop/profile* considers particles
crossing a set of planes, while compute *stress/cartesian* computes averages
for a set of small volumes. More information
on the similarities and differences can be found in
:ref:`(Ikeshoji)<Ikeshoji2>`.
NOTE 3: The local stress profile generated by compute
*stress/mop/profile* is similar to that obtained by compute
:doc:`stress/cartesian <compute_stress_cartesian>`.
A key difference is that compute *stress/mop/profile*
considers particles crossing a set of planes, while
*stress/cartesian* computes averages for a set
of small volumes.
Moreover, *stress/cartesian* compute computes the diagonal components of the stress
tensor :math:`P_{xx}`, :math:`P_{yy}`, and :math:`P_{zz}`, while
*stress/mop/profile* computes the components
:math:`P_{ix}`, :math:`P_{iy}`, and :math:`P_{iz}`, where :math:`i` is the
direction normal to the plane.
Output info
"""""""""""
Compute *stress/mop* calculates a global vector (indices starting at 1), with 3
values for each declared keyword (in the order the keywords have been
declared). For each keyword, the stress tensor components are ordered as
follows: stress_dir,x, stress_dir,y, and stress_dir,z.
Compute *stress/mop* calculates a global vector (indices starting at 1),
with 3 values for each declared keyword (in the order the keywords have
been declared). For each keyword, the stress tensor components are
ordered as follows: stress_dir,x, stress_dir,y, and stress_dir,z.
Compute *stress/mop/profile* instead calculates a global array, with 1 column
giving the position of the planes where the stress tensor was computed,
and with 3 columns of values for each declared keyword (in the order the
keywords have been declared). For each keyword, the profiles of stress
tensor components are ordered as follows: stress_dir,x; stress_dir,y;
and stress_dir,z.
Compute *stress/mop/profile* instead calculates a global array, with 1
column giving the position of the planes where the stress tensor was
computed, and with 3 columns of values for each declared keyword (in the
order the keywords have been declared). For each keyword, the profiles
of stress tensor components are ordered as follows: stress_dir,x;
stress_dir,y; and stress_dir,z.
The values are in pressure :doc:`units <units>`.
The values produced by this compute can be accessed by various :doc:`output commands <Howto_output>`.
For instance, the results can be written to a file using the
:doc:`fix ave/time <fix_ave_time>` command. Please see the example
in the examples/PACKAGES/mop folder.
The values produced by this compute can be accessed by various
:doc:`output commands <Howto_output>`. For instance, the results can be
written to a file using the :doc:`fix ave/time <fix_ave_time>`
command. Please see the example in the examples/PACKAGES/mop folder.
Restrictions
""""""""""""
These styles are part of the EXTRA-COMPUTE package. They are only enabled if
LAMMPS is built with that package. See the :doc:`Build package <Build_package>`
doc page on for more info.
These styles are part of the EXTRA-COMPUTE package. They are only
enabled if LAMMPS is built with that package. See the :doc:`Build
package <Build_package>` doc page on for more info.
The method is only implemented for 3d orthogonal simulation boxes whose
size does not change in time, and axis-aligned planes.
The method only works with two-body pair interactions, because it
requires the class method pair->single() to be implemented. In
particular, it does not work with more than two-body pair interactions,
intra-molecular interactions, and long range (kspace) interactions.
requires the class method ``Pair::single()`` to be implemented, which is
not possible for manybody potentials. In particular, compute
*stress/mop/profile* does not work with more than two-body pair
interactions, long range (kspace) interactions and
angle/dihedral/improper intramolecular interactions. Similarly, compute
*stress/mop* does not work with more than two-body pair interactions,
long range (kspace) interactions and dihedral/improper intramolecular
interactions but works with all bond interactions with the class method
single() implemented and all angle interactions with the class method
born_matrix() implemented.
Related commands
""""""""""""""""
:doc:`compute stress/atom <compute_stress_atom>`, :doc:`compute pressure <compute_pressure>`, :doc:`compute stress/cartesian <compute_stress_profile>`, :doc:`compute stress/cylinder <compute_stress_profile>`, :doc:`compute stress/spherical <compute_stress_profile>`
:doc:`compute stress/atom <compute_stress_atom>`,
:doc:`compute pressure <compute_pressure>`,
:doc:`compute stress/cartesian <compute_stress_cartesian>`,
:doc:`compute stress/cylinder <compute_stress_curvilinear>`,
:doc:`compute stress/spherical <compute_stress_curvilinear>`
Default
"""""""

View File

@ -753,9 +753,13 @@ run, this option is ignored since the output is already balanced.
----------
The *thermo* keyword only applies the dump styles *netcdf* and *yaml*.
It triggers writing of :doc:`thermo <thermo>` information to the dump file
alongside per-atom data. The values included in the dump file are
identical to the values specified by :doc:`thermo_style <thermo_style>`.
It triggers writing of :doc:`thermo <thermo>` information to the dump
file alongside per-atom data. The values included in the dump file are
cached values from the last thermo output and include the exact same the
values as specified by the :doc:`thermo_style <thermo_style>` command.
Because these are cached values, they are only up-to-date when dump
output is on a timestep that also has thermo output. Dump style *yaml*
will skip thermo output on incompatible steps.
----------

View File

@ -323,7 +323,8 @@ accelerated styles exist.
* :doc:`pafi <fix_pafi>` - constrained force averages on hyper-planes to compute free energies (PAFI)
* :doc:`pair <fix_pair>` - access per-atom info from pair styles
* :doc:`phonon <fix_phonon>` - calculate dynamical matrix from MD simulations
* :doc:`pimd/nvt <fix_pimd>` - Feynman path integral molecular dynamics with Nose-Hoover thermostat
* :doc:`pimd/langevin <fix_pimd>` - Feynman path-integral molecular dynamics with stochastic thermostat
* :doc:`pimd/nvt <fix_pimd>` - Feynman path-integral molecular dynamics with Nose-Hoover thermostat
* :doc:`planeforce <fix_planeforce>` - constrain atoms to move in a plane
* :doc:`plumed <fix_plumed>` - wrapper on PLUMED free energy library
* :doc:`poems <fix_poems>` - constrain clusters of atoms to move as coupled rigid bodies

View File

@ -61,7 +61,7 @@ heat conduction with a source term (see Fig.12 in :ref:`(Li2014) <Li2014b>`)
or diffusion with a source term (see Fig.1 in :ref:`(Li2015) <Li2015b>`), as
an analog of a periodic Poiseuille flow problem.
.. deprecated:: TBD
.. deprecated:: 15Jun2023
The *sphere* and *cuboid* keywords will be removed in a future version
of LAMMPS. The same functionality and more can be achieved with a region.

View File

@ -136,7 +136,7 @@ due to the electric field were a spring-like F = kx, then the energy
formula should be E = -0.5kx\^2. If you don't do this correctly, the
minimization will not converge properly.
.. versionadded:: TBD
.. versionadded:: 15Jun2023
The *potential* keyword can be used as an alternative to the *energy* keyword
to specify the name of an atom-style variable, which is used to compute the

View File

@ -85,7 +85,7 @@ columns 4-6 will store the "uinp" values.
.. code-block:: LAMMPS
pair_style amoeba
fix ex all pair amoeba 10 uind 0 uinp 0
fix ex all pair 10 amoeba uind 0 uinp 0
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""

View File

@ -1,5 +1,9 @@
.. index:: fix pimd/langevin
.. index:: fix pimd/nvt
fix pimd/langevin command
=========================
fix pimd/nvt command
====================
@ -8,72 +12,107 @@ Syntax
.. parsed-literal::
fix ID group-ID pimd/nvt keyword value ...
fix ID group-ID style keyword value ...
* ID, group-ID are documented in :doc:`fix <fix>` command
* pimd/nvt = style name of this fix command
* style = *pimd/langevin* or *pimd/nvt* = style name of this fix command
* zero or more keyword/value pairs may be appended
* keyword = *method* or *fmass* or *sp* or *temp* or *nhc*
* keywords for style *pimd/nvt*
.. parsed-literal::
*keywords* = *method* or *fmass* or *sp* or *temp* or *nhc*
*method* value = *pimd* or *nmpimd* or *cmd*
*fmass* value = scaling factor on mass
*sp* value = scaling factor on Planck constant
*temp* value = temperature (temperarate units)
*temp* value = temperature (temperature units)
*nhc* value = Nc = number of chains in Nose-Hoover thermostat
* keywords for style *pimd/langevin*
.. parsed-literal::
*keywords* = *method* or *integrator* or *ensemble* or *fmmode* or *fmass* or *scale* or *temp* or *thermostat* or *tau* or *iso* or *aniso* or *barostat* or *taup* or *fixcom* or *lj*
*method* value = *nmpimd*
*integrator* value = *obabo* or *baoab*
*fmmode* value = *physical* or *normal*
*fmass* value = scaling factor on mass
*temp* value = temperature (temperature unit)
temperature = target temperature of the thermostat
*thermostat* values = style seed
style value = *PILE_L*
seed = random number generator seed
*tau* value = thermostat damping parameter (time unit)
*scale* value = scaling factor of the damping times of non-centroid modes of PILE_L thermostat
*iso* or *aniso* values = pressure (pressure unit)
pressure = scalar external pressure of the barostat
*barostat* value = *BZP* or *MTTK*
*taup* value = barostat damping parameter (time unit)
*fixcom* value = *yes* or *no*
*lj* values = epsilon sigma mass planck mvv2e
epsilon = energy scale for reduced units (energy units)
sigma = length scale for reduced units (length units)
mass = mass scale for reduced units (mass units)
planck = Planck's constant for other unit style
mvv2e = mass * velocity^2 to energy conversion factor for other unit style
Examples
""""""""
.. code-block:: LAMMPS
fix 1 all pimd/nvt method nmpimd fmass 1.0 sp 2.0 temp 300.0 nhc 4
fix 1 all pimd/langevin ensemble npt integrator obabo temp 113.15 thermostat PILE_L 1234 tau 1.0 iso 1.0 barostat BZP taup 1.0
Description
"""""""""""
.. versionchanged:: 28Mar2023
Fix pimd was renamed to fix pimd/nvt.
Fix pimd was renamed to fix *pimd/nvt* and fix *pimd/langevin* was added.
This command performs quantum molecular dynamics simulations based on
the Feynman path integral to include effects of tunneling and
These fix commands perform quantum molecular dynamics simulations based
on the Feynman path-integral to include effects of tunneling and
zero-point motion. In this formalism, the isomorphism of a quantum
partition function for the original system to a classical partition
function for a ring-polymer system is exploited, to efficiently sample
configurations from the canonical ensemble :ref:`(Feynman) <Feynman>`.
The classical partition function and its components are given
by the following equations:
.. math::
Z = & \int d{\bf q} d{\bf p} \cdot \textrm{exp} [ -\beta H_{eff} ] \\
H_{eff} = & \bigg(\sum_{i=1}^P \frac{p_i^2}{2m_i}\bigg) + V_{eff} \\
H_{eff} = & \bigg(\sum_{i=1}^P \frac{p_i^2}{2M_i}\bigg) + V_{eff} \\
V_{eff} = & \sum_{i=1}^P \bigg[ \frac{mP}{2\beta^2 \hbar^2} (q_i - q_{i+1})^2 + \frac{1}{P} V(q_i)\bigg]
:math:`M_i` is the fictitious mass of the :math:`i`-th mode, and m is the actual mass of the atoms.
The interested user is referred to any of the numerous references on
this methodology, but briefly, each quantum particle in a path
integral simulation is represented by a ring-polymer of P quasi-beads,
labeled from 1 to P. During the simulation, each quasi-bead interacts
with beads on the other ring-polymers with the same imaginary time
index (the second term in the effective potential above). The
quasi-beads also interact with the two neighboring quasi-beads through
the spring potential in imaginary-time space (first term in effective
potential). To sample the canonical ensemble, a Nose-Hoover massive
chain thermostat is applied :ref:`(Tuckerman) <pimd-Tuckerman>`. With the
massive chain algorithm, a chain of NH thermostats is coupled to each
degree of freedom for each quasi-bead. The keyword *temp* sets the
target temperature for the system and the keyword *nhc* sets the
number *Nc* of thermostats in each chain. For example, for a
simulation of N particles with P beads in each ring-polymer, the total
number of NH thermostats would be 3 x N x P x Nc.
this methodology, but briefly, each quantum particle in a path integral
simulation is represented by a ring-polymer of P quasi-beads, labeled
from 1 to P. During the simulation, each quasi-bead interacts with
beads on the other ring-polymers with the same imaginary time index (the
second term in the effective potential above). The quasi-beads also
interact with the two neighboring quasi-beads through the spring
potential in imaginary-time space (first term in effective potential).
To sample the canonical ensemble, any thermostat can be applied.
Fix *pimd/nvt* applies a Nose-Hoover massive chain thermostat
:ref:`(Tuckerman) <pimd-Tuckerman>`. With the massive chain
algorithm, a chain of NH thermostats is coupled to each degree of
freedom for each quasi-bead. The keyword *temp* sets the target
temperature for the system and the keyword *nhc* sets the number *Nc* of
thermostats in each chain. For example, for a simulation of N particles
with P beads in each ring-polymer, the total number of NH thermostats
would be 3 x N x P x Nc.
Fix *pimd/langevin* implements a Langevin thermostat in the normal mode
representation, and also provides a barostat to sample the NPH/NPT ensembles.
.. note::
Fix pimd/nvt implements a complete velocity-verlet integrator
combined with NH massive chain thermostat, so no other time
integration fix should be used.
Both these *fix* styles implement a complete velocity-verlet integrator
combined with a thermostat, so no other time integration fix should be used.
The *method* keyword determines what style of PIMD is performed. A
value of *pimd* is standard PIMD. A value of *nmpimd* is for
@ -81,7 +120,7 @@ normal-mode PIMD. A value of *cmd* is for centroid molecular dynamics
(CMD). The difference between the styles is as follows.
In standard PIMD, the value used for a bead's fictitious mass is
arbitrary. A common choice is to use Mi = m/P, which results in the
arbitrary. A common choice is to use :math:`M_i = m/P`, which results in the
mass of the entire ring-polymer being equal to the real quantum
particle. But it can be difficult to efficiently integrate the
equations of motion for the stiff harmonic interactions in the ring
@ -97,6 +136,10 @@ normal-mode PIMD. A value of *cmd* is for centroid molecular dynamics
overall translation of the ring-polymer and is assigned the mass of
the real particle.
.. note::
Fix pimd/langevin only supports *method* value *nmpimd*. This should be enough
for most PIMD applications for quantum thermodynamics purpose.
Motion of the centroid can be effectively uncoupled from the other
normal modes by scaling the fictitious masses to achieve a partial
adiabatic separation. This is called a Centroid Molecular Dynamics
@ -108,17 +151,86 @@ normal-mode PIMD. A value of *cmd* is for centroid molecular dynamics
only the k > 0 modes are thermostatted, not the centroid degrees of
freedom.
The keyword *integrator* specifies the Trotter splitting method used by *fix pimd/langevin*.
See :ref:`(Liu) <Liu>` for a discussion on the OBABO and BAOAB splitting schemes. Typically
either of the two should work fine.
The keyword *fmass* sets a further scaling factor for the fictitious
masses of beads, which can be used for the Partial Adiabatic CMD
:ref:`(Hone) <Hone>`, or to be set as P, which results in the fictitious
masses to be equal to the real particle masses.
The keyword *fmmode* of *fix pimd/langevin* determines the mode of fictitious
mass preconditioning. There are two options: *physical* and *normal*. If *fmmode* is
*physical*, then the physical mass of the particles are used (and then multiplied by
*fmass*). If *fmmode* is *normal*, then the physical mass is first multiplied by the
eigenvalue of each normal mode, and then multiplied by *fmass*. More precisely, the
fictitious mass of *fix pimd/langevin* is determined by two factors: *fmmode* and *fmass*.
If *fmmode* is *physical*, then the fictitious mass is
.. math::
M_i = \mathrm{fmass} \times m
If *fmmode* is *normal*, then the fictitious mass is
.. math::
M_i = \mathrm{fmass} \times \lambda_i \times m
where :math:`\lambda_i` is the eigenvalue of the :math:`i`-th normal mode.
.. note::
Fictitious mass is only used in the momentum of the equation of motion
(:math:`\mathbf{p}_i=M_i\mathbf{v}_i`), and not used in the spring elastic energy
(:math:`\sum_{i=1}^P \frac{1}{2}m\omega_P^2(q_i - q_{i+1})^2`, :math:`m` is always the
actual mass of the particles).
The keyword *sp* is a scaling factor on Planck's constant, which can
be useful for debugging or other purposes. The default value of 1.0
is appropriate for most situations.
The keyword *ensemble* for fix style *pimd/langevin* determines which ensemble is it
going to sample. The value can be *nve* (microcanonical), *nvt* (canonical), *nph* (isoenthalpic),
and *npt* (isothermal-isobaric).
The keyword *temp* specifies temperature parameter for fix styles *pimd/nvt* and *pimd/langevin*. It should read
a positive floating-point number.
.. note::
For pimd simulations, a temperature values should be specified even for nve ensemble. Temperature will make a difference
for nve pimd, since the spring elastic frequency between the beads will be affected by the temperature.
The keyword *thermostat* reads *style* and *seed* of thermostat for fix style *pimd/langevin*. *style* can only
be *PILE_L* (path integral Langevin equation local thermostat, as described in :ref:`Ceriotti <Ceriotti2>`), and *seed* should a positive integer number, which serves as the seed of the pseudo random number generator.
.. note::
The fix style *pimd/langevin* uses the stochastic PILE_L thermostat to control temperature. This thermostat works on the normal modes
of the ring polymer. The *tau* parameter controls the centroid mode, and the *scale* parameter controls the non-centroid modes.
The keyword *tau* specifies the thermostat damping time parameter for fix style *pimd/langevin*. It is in time unit. It only works on the centroid mode.
The keyword *scale* specifies a scaling parameter for the damping times of the non-centroid modes for fix style *pimd/langevin*. The default
damping time of the non-centroid mode :math:`i` is :math:`\frac{P}{\beta\hbar}\sqrt{\lambda_i\times\mathrm{fmass}}` (*fmmode* is *physical*) or :math:`\frac{P}{\beta\hbar}\sqrt{\mathrm{fmass}}` (*fmmode* is *normal*). The damping times of all non-centroid modes are the default values divided by *scale*.
The barostat parameters for fix style *pimd/langevin* with *npt* or *nph* ensemble is specified using one of *iso* and *aniso*
keywords. A *pressure* value should be given with pressure unit. The keyword *iso* means couple all 3 diagonal components together when pressure is computed (hydrostatic pressure), and dilate/contract the dimensions together. The keyword *aniso* means x, y, and z dimensions are controlled independently using the Pxx, Pyy, and Pzz components of the stress tensor as the driving forces, and the specified scalar external pressure.
The keyword *barostat* reads *style* of barostat for fix style *pimd/langevin*. *style* can
be *BZP* (Bussi-Zykova-Parrinello, as described in :ref:`Bussi <Bussi>`) or *MTTK* (Martyna-Tuckerman-Tobias-Klein, as described in :ref:`Martyna1 <Martyna3>` and :ref:`Martyna2 <Martyna4>`).
The keyword *taup* specifies the barostat damping time parameter for fix style *pimd/langevin*. It is in time unit.
The keyword *fixcom* specifies whether the center-of-mass of the extended ring-polymer system is fixed during the pimd simulation.
Once *fixcom* is set to be *yes*, the center-of-mass velocity will be distracted from the centroid-mode velocities in each step.
The keyword *lj* should be used if :doc:`lj units <units>` is used for *fix pimd/langevin*. Typically one may want to use
reduced units to run the simulation, and then convert the results into some physical units (for example, :doc:`metal units <units>`). In this case, the 5 quantities in the physical mass units are needed: epsilon (energy scale), sigma (length scale), mass, Planck's constant, mvv2e (mass * velocity^2 to energy conversion factor). Planck's constant and mvv2e can be found in src/update.cpp. If there is no need to convert reduced units to physical units, set all these five value to 1.
The PIMD algorithm in LAMMPS is implemented as a hyper-parallel scheme
as described in :ref:`(Calhoun) <Calhoun>`. In LAMMPS this is done by using
as described in :ref:`Calhoun <Calhoun>`. In LAMMPS this is done by using
:doc:`multi-replica feature <Howto_replica>` in LAMMPS, where each
quasi-particle system is stored and simulated on a separate partition
of processors. The following diagram illustrates this approach. The
@ -152,22 +264,49 @@ related tasks for each of the partitions, e.g.
.. code-block:: LAMMPS
dump dcd all dcd 10 system_${ibead}.dcd
dump 1 all custom 100 ${ibead}.xyz id type x y z vx vy vz ix iy iz fx fy fz
restart 1000 system_${ibead}.restart1 system_${ibead}.restart2
read_restart system_${ibead}.restart2
.. note::
Fix *pimd/langevin* dumps the Cartesian coordinates, but dumps the velocities and
forces in the normal mode representation. If the Cartesian velocities and forces are
needed, it is easy to perform the transformation when doing post-processing.
It is recommended to dump the image flags (*ix iy iz*) for fix *pimd/langevin*. It
will be useful if you want to calculate some estimators during post-processing.
Major differences of *fix pimd/nvt* and *fix pimd/langevin* are:
#. *Fix pimd/nvt* includes Cartesian pimd, normal mode pimd, and centroid md. *Fix pimd/langevin* only intends to support normal mode pimd, as it is commonly enough for thermodynamic sampling.
#. *Fix pimd/nvt* uses Nose-Hoover chain thermostat. *Fix pimd/langevin* uses Langevin thermostat.
#. *Fix pimd/langevin* provides barostat, so the npt ensemble can be sampled. *Fix pimd/nvt* only support nvt ensemble.
#. *Fix pimd/langevin* provides several quantum estimators in output.
#. *Fix pimd/langevin* allows multiple processes for each bead. For *fix pimd/nvt*, there is a large chance that multi-process tasks for each bead may fail.
#. The dump of *fix pimd/nvt* are all Cartesian. *Fix pimd/langevin* dumps normal-mode velocities and forces, and Cartesian coordinates.
Initially, the inter-replica communication and normal mode transformation parts of *fix pimd/langevin* are written based on
those of *fix pimd/nvt*, but are significantly revised.
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
Fix pimd/nvt writes the state of the Nose/Hoover thermostat over all
Fix *pimd/nvt* writes the state of the Nose/Hoover thermostat over all
quasi-beads 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.
Fix *pimd/langevin* writes the state of the barostat overall beads to
:doc:`binary restart files <restart>`. Since it uses a stochastic thermostat,
the state of the thermostat is not written. However, the state of the system
can be restored by reading the restart file, except that it will re-initialize
the random number generator.
None of the :doc:`fix_modify <fix_modify>` options
are relevant to fix pimd/nvt.
Fix pimd/nvt computes a global 3-vector, which can be accessed by
Fix *pimd/nvt* computes a global 3-vector, which can be accessed by
various :doc:`output commands <Howto_output>`. The three quantities in
the global vector are:
@ -176,21 +315,80 @@ the global vector are:
#. the current value of the scalar virial estimator for the kinetic
energy of the quantum system :ref:`(Herman) <Herman>`.
The vector values calculated by fix pimd/nvt are "extensive", except for the
The vector values calculated by fix *pimd/nvt* are "extensive", except for the
temperature, which is "intensive".
No parameter of fix pimd/nvt can be used with the *start/stop* keywords
of the :doc:`run <run>` command. Fix pimd/nvt is not invoked during
Fix *pimd/langevin* computes a global vector of quantities, which
can be accessed by various :doc:`output commands <Howto_output>`. Note that
it outputs multiple log files, and different log files contain information
about different beads or modes (see detailed explanations below). If *ensemble*
is *nve* or *nvt*, the vector has 10 values:
#. kinetic energy of the normal mode
#. spring elastic energy of the normal mode
#. potential energy of the bead
#. total energy of all beads (conserved if *ensemble* is *nve*)
#. primitive kinetic energy estimator
#. virial energy estimator
#. centroid-virial energy estimator
#. primitive pressure estimator
#. thermodynamic pressure estimator
#. centroid-virial pressure estimator
The first 3 are different for different log files, and the others are the same for different log files.
If *ensemble* is *nph* or *npt*, the vector stores internal variables of the barostat. If *iso* is used,
the vector has 15 values:
#. kinetic energy of the normal mode
#. spring elastic energy of the normal mode
#. potential energy of the bead
#. total energy of all beads (conserved if *ensemble* is *nve*)
#. primitive kinetic energy estimator
#. virial energy estimator
#. centroid-virial energy estimator
#. primitive pressure estimator
#. thermodynamic pressure estimator
#. centroid-virial pressure estimator
#. barostat velocity
#. barostat kinetic energy
#. barostat potential energy
#. barostat cell Jacobian
#. enthalpy of the extended system (sum of 4, 12, 13, and 14; conserved if *ensemble* is *nph*)
If *aniso* or *x* or *y* or *z* is used for the barostat, the vector has 17 values:
#. kinetic energy of the normal mode
#. spring elastic energy of the normal mode
#. potential energy of the bead
#. total energy of all beads (conserved if *ensemble* is *nve*)
#. primitive kinetic energy estimator
#. virial energy estimator
#. centroid-virial energy estimator
#. primitive pressure estimator
#. thermodynamic pressure estimator
#. centroid-virial pressure estimator
#. x component of barostat velocity
#. y component of barostat velocity
#. z component of barostat velocity
#. barostat kinetic energy
#. barostat potential energy
#. barostat cell Jacobian
#. enthalpy of the extended system (sum of 4, 14, 15, and 16; conserved if *ensemble* is *nph*)
No parameter of fix *pimd/nvt* or *pimd/langevin* can be used with the *start/stop* keywords
of the :doc:`run <run>` command. Fix *pimd/nvt* or *pimd/langevin* is not invoked during
:doc:`energy minimization <minimize>`.
Restrictions
""""""""""""
This fix is part of the REPLICA package. It is only enabled if LAMMPS
was built with that package. See the :doc:`Build package
These fixes are part of the REPLICA package. They are only enabled if
LAMMPS was built with that package. See the :doc:`Build package
<Build_package>` page for more info.
Fix pimd/nvt cannot be used with :doc:`lj units <units>`.
Fix *pimd/nvt* cannot be used with :doc:`lj units <units>`.
Fix *pimd/langevin* can be used with :doc:`lj units <units>`. See the above part for how to use it.
A PIMD simulation can be initialized with a single data file read via
the :doc:`read_data <read_data>` command. However, this means all
@ -207,7 +405,7 @@ variable, e.g.
Default
"""""""
The keyword defaults for fix pimd/nvt are method = pimd, fmass = 1.0, sp
The keyword defaults for fix *pimd/nvt* are method = pimd, fmass = 1.0, sp
= 1.0, temp = 300.0, and nhc = 2.
----------
@ -243,3 +441,22 @@ Path Integrals, McGraw-Hill, New York (1965).
**(Herman)** M. F. Herman, E. J. Bruskin, B. J. Berne, J Chem Phys, 76, 5150 (1982).
.. _Bussi:
**(Bussi)** G. Bussi, T. Zykova-Timan, M. Parrinello, J Chem Phys, 130, 074101 (2009).
.. _Ceriotti3:
**(Ceriotti)** M. Ceriotti, M. Parrinello, T. Markland, D. Manolopoulos, J. Chem. Phys. 133, 124104 (2010).
.. _Martyna3:
**(Martyna1)** G. Martyna, D. Tobias, M. Klein, J. Chem. Phys. 101, 4177 (1994).
.. _Martyna4:
**(Martyna2)** G. Martyna, A. Hughes, M. Tuckerman, J. Chem. Phys. 110, 3275 (1999).
.. _Liujian:
**(Liu)** J. Liu, D. Li, X. Liu, J. Chem. Phys. 145, 024103 (2016).

View File

@ -44,7 +44,7 @@ one word. If it contains variables it must be enclosed in double
quotes to ensure they are not evaluated when the input script line is
read, but will instead be evaluated each time the string is printed.
.. versionadded:: TBD
.. versionadded:: 15Jun2023
support for vector style variables

View File

@ -25,7 +25,7 @@ Syntax
.. parsed-literal::
*cutoff* value = I J Cutoff
I, J = atom types
I, J = atom types (see asterisk form below)
Cutoff = Bond-order cutoff value for this pair of atom types
*element* value = Element1, Element2, ...
*position* value = posfreq filepos
@ -49,7 +49,7 @@ Examples
.. code-block:: LAMMPS
fix 1 all reaxff/species 10 10 100 species.out
fix 1 all reaxff/species 1 2 20 species.out cutoff 1 1 0.40 cutoff 1 2 0.55
fix 1 all reaxff/species 1 2 20 species.out cutoff 1 1 0.40 cutoff 1 2*3 0.55
fix 1 all reaxff/species 1 100 100 species.out element Au O H position 1000 AuOH.pos
fix 1 all reaxff/species 1 100 100 species.out delete species.del masslimit 0 50
@ -88,13 +88,24 @@ If the filename ends with ".gz", the output file is written in gzipped
format. A gzipped dump file will be about 3x smaller than the text version,
but will also take longer to write.
.. versionadded:: 15Jun2023
Support for wildcards added
Optional keyword *cutoff* can be assigned to change the minimum
bond-order values used in identifying chemical bonds between pairs of
atoms. Bond-order cutoffs should be carefully chosen, as bond-order
cutoffs that are too small may include too many bonds (which will
result in an error), while cutoffs that are too large will result in
fragmented molecules. The default cutoff of 0.3 usually gives good
results.
cutoffs that are too small may include too many bonds (which will result
in an error), while cutoffs that are too large will result in fragmented
molecules. The default cutoff of 0.3 usually gives good results. A
wildcard asterisk can be used in place of or in conjunction with the I,J
arguments to set the bond-order cutoff for multiple pairs of atom types.
This takes the form "\*" or "\*n" or "n\*" or "m\*n". If :math:`N` is
the number of atom types, then an asterisk with no numeric values means
all types from 1 to :math:`N`. A leading asterisk means all types from
1 to n (inclusive). A trailing asterisk means all types from n to
:math:`N` (inclusive). A middle asterisk means all types from m to n
(inclusive).
The optional keyword *element* can be used to specify the chemical
symbol printed for each LAMMPS atom type. The number of symbols must

View File

@ -49,7 +49,7 @@ Syntax
*intel* args = NPhi keyword value ...
Nphi = # of co-processors per node
zero or more keyword/value pairs may be appended
keywords = *mode* or *omp* or *lrt* or *balance* or *ghost* or *tpc* or *tptask* or *no_affinity*
keywords = *mode* or *omp* or *lrt* or *balance* or *ghost* or *tpc* or *tptask* or *pppm_table* or *no_affinity*
*mode* value = *single* or *mixed* or *double*
single = perform force calculations in single precision
mixed = perform force calculations in mixed precision
@ -68,6 +68,9 @@ Syntax
Ntpc = max number of co-processor threads per co-processor core (default = 4)
*tptask* value = Ntptask
Ntptask = max number of co-processor threads per MPI task (default = 240)
*pppm_table* value = *yes* or *no*
*yes* = Precompute pppm values in table (doesn't change accuracy)
*no* = Compute pppm values on the fly
*no_affinity* values = none
*kokkos* args = keyword value ...
zero or more keyword/value pairs may be appended
@ -392,13 +395,13 @@ is identical to the default *verlet* style aside from supporting the
LRT feature. This feature requires setting the pre-processor flag
-DLMP_INTEL_USELRT in the makefile when compiling LAMMPS.
The *balance* keyword sets the fraction of :doc:`pair style <pair_style>` work offloaded to the co-processor for split
values between 0.0 and 1.0 inclusive. While this fraction of work is
running on the co-processor, other calculations will run on the host,
including neighbor and pair calculations that are not offloaded, as
well as angle, bond, dihedral, kspace, and some MPI communications.
If *split* is set to -1, the fraction of work is dynamically adjusted
automatically throughout the run. This typically give performance
The *balance* keyword sets the fraction of :doc:`pair style <pair_style>` work
offloaded to the co-processor for split values between 0.0 and 1.0 inclusive.
While this fraction of work is running on the co-processor, other calculations
will run on the host, including neighbor and pair calculations that are not
offloaded, as well as angle, bond, dihedral, kspace, and some MPI
communications. If *split* is set to -1, the fraction of work is dynamically
adjusted automatically throughout the run. This typically give performance
within 5 to 10 percent of the optimal fixed fraction.
The *ghost* keyword determines whether or not ghost atoms, i.e. atoms
@ -430,6 +433,13 @@ with 16 threads, for a total of 128.
Note that the default settings for *tpc* and *tptask* are fine for
most problems, regardless of how many MPI tasks you assign to a Phi.
.. versionadded:: 15Jun2023
The *pppm_table* keyword with the argument yes allows to use a
pre-computed table to efficiently spread the charge to the PPPM grid.
This feature is enabled by default but can be turned off using the
keyword with the argument *no*.
The *no_affinity* keyword will turn off automatic setting of core
affinity for MPI tasks and OpenMP threads on the host when using
offload to a co-processor. Affinity settings are used when possible
@ -699,39 +709,37 @@ Related commands
Default
"""""""
For the GPU package, the default is Ngpu = 0 and the option defaults
are neigh = yes, newton = off, binsize = 0.0, split = 1.0, gpuID = 0
to Ngpu-1, tpa = 1, omp = 0, and platform=-1. These settings are made
automatically if the "-sf gpu" :doc:`command-line switch <Run_options>`
is used. If it is not used, you must invoke the package gpu command
in your input script or via the "-pk gpu" :doc:`command-line switch <Run_options>`.
For the GPU package, the default is Ngpu = 0 and the option defaults are neigh
= yes, newton = off, binsize = 0.0, split = 1.0, gpuID = 0 to Ngpu-1, tpa = 1,
omp = 0, and platform=-1. These settings are made automatically if the "-sf
gpu" :doc:`command-line switch <Run_options>` is used. If it is not used, you
must invoke the package gpu command in your input script or via the "-pk gpu"
:doc:`command-line switch <Run_options>`.
For the INTEL package, the default is Nphi = 1 and the option
defaults are omp = 0, mode = mixed, lrt = no, balance = -1, tpc = 4,
tptask = 240. The default ghost option is determined by the pair
style being used. This value is output to the screen in the offload
report at the end of each run. Note that all of these settings,
except "omp" and "mode", are ignored if LAMMPS was not built with Xeon
Phi co-processor support. These settings are made automatically if the
"-sf intel" :doc:`command-line switch <Run_options>` is used. If it is
not used, you must invoke the package intel command in your input
script or via the "-pk intel" :doc:`command-line switch <Run_options>`.
For the INTEL package, the default is Nphi = 1 and the option defaults are omp
= 0, mode = mixed, lrt = no, balance = -1, tpc = 4, tptask = 240, pppm_table =
yes. The default ghost option is determined by the pair style being used.
This value is output to the screen in the offload report at the end of each
run. Note that all of these settings, except "omp" and "mode", are ignored if
LAMMPS was not built with Xeon Phi co-processor support. These settings are
made automatically if the "-sf intel" :doc:`command-line switch <Run_options>`
is used. If it is not used, you must invoke the package intel command in your
input script or via the "-pk intel" :doc:`command-line switch <Run_options>`.
For the KOKKOS package, the option defaults for GPUs are neigh = full,
neigh/qeq = full, newton = off, binsize for GPUs = 2x LAMMPS default
value, comm = device, sort = device, neigh/transpose = off, gpu/aware =
on. When LAMMPS can safely detect that GPU-aware MPI is not available,
the default value of gpu/aware becomes "off". For CPUs or Xeon Phis, the
option defaults are neigh = half, neigh/qeq = half, newton = on, binsize
= 0.0, comm = no, and sort = no. The option neigh/thread = on when
there are 16K atoms or less on an MPI rank, otherwise it is "off". These
settings are made automatically by the required "-k on"
:doc:`command-line switch <Run_options>`. You can change them by using
the package kokkos command in your input script or via the :doc:`-pk
neigh/qeq = full, newton = off, binsize for GPUs = 2x LAMMPS default value,
comm = device, sort = device, neigh/transpose = off, gpu/aware = on. When
LAMMPS can safely detect that GPU-aware MPI is not available, the default value
of gpu/aware becomes "off". For CPUs or Xeon Phis, the option defaults are
neigh = half, neigh/qeq = half, newton = on, binsize = 0.0, comm = no, and sort
= no. The option neigh/thread = on when there are 16K atoms or less on an MPI
rank, otherwise it is "off". These settings are made automatically by the
required "-k on" :doc:`command-line switch <Run_options>`. You can change them
by using the package kokkos command in your input script or via the :doc:`-pk
kokkos command-line switch <Run_options>`.
For the OMP package, the default is Nthreads = 0 and the option
defaults are neigh = yes. These settings are made automatically if
the "-sf omp" :doc:`command-line switch <Run_options>` is used. If it
is not used, you must invoke the package omp command in your input
script or via the "-pk omp" :doc:`command-line switch <Run_options>`.
For the OMP package, the default is Nthreads = 0 and the option defaults are
neigh = yes. These settings are made automatically if the "-sf omp"
:doc:`command-line switch <Run_options>` is used. If it is not used, you must
invoke the package omp command in your input script or via the "-pk omp"
:doc:`command-line switch <Run_options>`.

View File

@ -0,0 +1,167 @@
.. index:: pair_style aip/water/2dm
.. index:: pair_style aip/water/2dm/opt
pair_style aip/water/2dm command
===================================
Accelerator Variant: *aip/water/2dm/opt*
Syntax
""""""
.. code-block:: LAMMPS
pair_style [hybrid/overlay ...] aip/water/2dm cutoff tap_flag
* cutoff = global cutoff (distance units)
* tap_flag = 0/1 to turn off/on the taper function
Examples
""""""""
.. code-block:: LAMMPS
pair_style hybrid/overlay aip/water/2dm 16.0 1
pair_coeff * * aip/water/2dm COH.aip.water.2dm C Ow Hw
pair_style hybrid/overlay aip/water/2dm 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm COH.aip.water.2dm C Ow Hw
Description
"""""""""""
.. versionadded:: 15Jun2023
The *aip/water/2dm* style computes the anisotropic interfacial potential
(AIP) potential for interfaces of water with two-dimensional (2D)
materials as described in :ref:`(Feng) <Feng>`.
.. math::
E = & \frac{1}{2} \sum_i \sum_{j \neq i} V_{ij} \\
V_{ij} = & {\rm Tap}(r_{ij})\left \{ e^{-\alpha (r_{ij}/\beta -1)}
\left [ \epsilon + f(\rho_{ij}) + f(\rho_{ji})\right ] -
\frac{1}{1+e^{-d\left [ \left ( r_{ij}/\left (s_R \cdot r^{eff} \right ) \right )-1 \right ]}}
\cdot \frac{C_6}{r^6_{ij}} \right \}\\
\rho_{ij}^2 = & r_{ij}^2 - ({\bf r}_{ij} \cdot {\bf n}_i)^2 \\
\rho_{ji}^2 = & r_{ij}^2 - ({\bf r}_{ij} \cdot {\bf n}_j)^2 \\
f(\rho) = & C e^{ -( \rho / \delta )^2 } \\
{\rm Tap}(r_{ij}) = & 20\left ( \frac{r_{ij}}{R_{cut}} \right )^7 -
70\left ( \frac{r_{ij}}{R_{cut}} \right )^6 +
84\left ( \frac{r_{ij}}{R_{cut}} \right )^5 -
35\left ( \frac{r_{ij}}{R_{cut}} \right )^4 + 1
Where :math:`\mathrm{Tap}(r_{ij})` is the taper function which provides
a continuous cutoff (up to third derivative) for interatomic separations
larger than :math:`r_c` :doc:`pair_style ilp_graphene_hbn
<pair_ilp_graphene_hbn>`.
.. note::
This pair style uses the atomic normal vector definition from
:ref:`(Feng) <Feng>`), where the atomic normal vectors of the
hydrogen atoms are assumed to lie along the corresponding
oxygen-hydrogen bonds and the normal vector of the central oxygen
atom is defined as their average.
The provided parameter file, ``COH.aip.water.2dm``, is intended for use
with *metal* :doc:`units <units>`, with energies in meV. Two additional
parameters, *S*, and *rcut* are included in the parameter file. *S* is
designed to facilitate scaling of energies; *rcut* is the cutoff for an
internal, short distance neighbor list that is generated for speeding up
the calculation of the normals for all atom pairs.
.. note::
The parameters presented in the provided parameter file,
``COH.aip.water.2dm``, are fitted with the taper function enabled by
setting the cutoff equal to 16.0 Angstrom. Using a different cutoff
or taper function setting should be carefully checked as they can
lead to significant errors. These parameters provide a good
description in both short- and long-range interaction regimes. This
is essential for simulations in high pressure regime (i.e., the
interlayer distance is smaller than the equilibrium distance).
This potential must be used in combination with hybrid/overlay. Other
interactions can be set to zero using :doc:`pair_coeff settings
<pair_coeff>` with the pair style set to *none*\ .
This pair style tallies a breakdown of the total interlayer potential
energy into sub-categories, which can be accessed via the :doc:`compute
pair <compute_pair>` command as a vector of values of length 2. The 2
values correspond to the following sub-categories:
1. *E_vdW* = vdW (attractive) energy
2. *E_Rep* = Repulsive energy
To print these quantities to the log file (with descriptive column
headings) the following commands could be included in an input script:
.. code-block:: LAMMPS
compute 0 all pair aip/water/2dm
variable Evdw equal c_0[1]
variable Erep equal c_0[2]
thermo_style custom step temp epair v_Erep v_Evdw
----------
.. include:: accel_styles.rst
----------
Mixing, shift, table, tail correction, restart, rRESPA info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
This pair style does not support the pair_modify mix, shift, table, and
tail options.
This pair style does not write their information to 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.
Restrictions
""""""""""""
This pair style is part of the INTERLAYER package. It is only enabled
if LAMMPS was built with that package. See the :doc:`Build package
<Build_package>` page for more info.
This pair style requires the newton setting to be *on* for pair
interactions.
The ``COH.aip.water.2dm`` potential file provided with LAMMPS is
parameterized for *metal* units. You can use this pair style with any
LAMMPS units, but you would need to create your own potential file with
parameters in the appropriate units, if your simulation does not use
*metal* units.
Related commands
""""""""""""""""
:doc:`pair_coeff <pair_coeff>`,
:doc:`pair_none <pair_none>`,
:doc:`pair_style hybrid/overlay <pair_hybrid>`,
:doc:`pair_style drip <pair_drip>`,
:doc:`pair_style ilp_tmd <pair_ilp_tmd>`,
:doc:`pair_style saip_metal <pair_saip_metal>`,
:doc:`pair_style ilp_graphene_hbn <pair_ilp_graphene_hbn>`,
:doc:`pair_style pair_kolmogorov_crespi_z <pair_kolmogorov_crespi_z>`,
:doc:`pair_style pair_kolmogorov_crespi_full <pair_kolmogorov_crespi_full>`,
:doc:`pair_style pair_lebedeva_z <pair_lebedeva_z>`,
:doc:`pair_style pair_coul_shield <pair_coul_shield>`.
Default
"""""""
tap_flag = 1
----------
.. _Feng:
**(Feng)** Z. Feng and W. Ouyang et al., J. Phys. Chem. C. 127, 8704-8713 (2023).

View File

@ -736,7 +736,7 @@ or
.. math::
E_{eff,ij} = \frac{E_{ij}}{2(1-\nu_{ij})}
E_{eff,ij} = \frac{E_{ij}}{2(1-\nu_{ij}^2)}
These pair styles write their information to :doc:`binary restart files <restart>`, so a pair_style command does not need to be
specified in an input script that reads a restart file.

View File

@ -155,8 +155,8 @@ interactions.
The BNCH.ILP potential file provided with LAMMPS (see the potentials
directory) are parameterized for *metal* units. You can use this
potential with any LAMMPS units, but you would need to create your
BNCH.ILP potential file with coefficients listed in the appropriate
potential with any LAMMPS units, but you would need to create your own
custom BNCH.ILP potential file with coefficients listed in the appropriate
units, if your simulation does not use *metal* units.
Related commands
@ -181,6 +181,14 @@ tap_flag = 1
----------
.. _Ouyang1:
**(Ouyang1)** W. Ouyang, D. Mandelli, M. Urbakh and O. Hod, Nano Lett. 18, 6009-6016 (2018).
.. _Ouyang2:
**(Ouyang2)** W. Ouyang et al., J. Chem. Theory Comput. 16(1), 666-676 (2020).
.. _Leven1:
**(Leven1)** I. Leven, I. Azuri, L. Kronik and O. Hod, J. Chem. Phys. 140, 104106 (2014).
@ -196,11 +204,3 @@ tap_flag = 1
.. _Kolmogorov2:
**(Kolmogorov)** A. N. Kolmogorov, V. H. Crespi, Phys. Rev. B 71, 235415 (2005).
.. _Ouyang1:
**(Ouyang1)** W. Ouyang, D. Mandelli, M. Urbakh and O. Hod, Nano Lett. 18, 6009-6016 (2018).
.. _Ouyang2:
**(Ouyang2)** W. Ouyang et al., J. Chem. Theory Comput. 16(1), 666-676 (2020).

View File

@ -135,8 +135,8 @@ interactions.
The TMD.ILP potential file provided with LAMMPS (see the potentials
directory) are parameterized for *metal* units. You can use this
potential with any LAMMPS units, but you would need to create your
BNCH.ILP potential file with coefficients listed in the appropriate
potential with any LAMMPS units, but you would need to create your own
custom TMD.ILP potential file with coefficients listed in the appropriate
units, if your simulation does not use *metal* units.
Related commands

View File

@ -129,7 +129,7 @@ interactions.
The CH.KC potential file provided with LAMMPS (see the potentials
folder) is parameterized for metal units. You can use this pair style
with any LAMMPS units, but you would need to create your own custom
CC.KC potential file with all coefficients converted to the appropriate
CH.KC potential file with all coefficients converted to the appropriate
units.
Related commands

View File

@ -61,9 +61,11 @@ Description
.. versionadded:: 8Feb2023
.. versionchanged:: TBD
added pair styles *lepton* and *lepton/coul*
added pair style lepton/sphere
.. versionchanged:: 15Jun2023
added pair style *lepton/sphere*
Pair styles *lepton*, *lepton/coul*, *lepton/sphere* compute pairwise
interactions between particles which depend on the distance and have a

View File

@ -33,7 +33,7 @@ Examples
Description
"""""""""""
.. versionadded:: TBD
.. versionadded:: 15Jun2023
The *lj/cut/sphere* style compute the standard 12/6 Lennard-Jones potential,
given by

View File

@ -33,7 +33,7 @@ Examples
Description
"""""""""""
.. versionadded:: TBD
.. versionadded:: 15Jun2023
The *lj/expand/sphere* style compute a 12/6 Lennard-Jones potential with
a distance shifted by :math:`\Delta = \frac{1}{2} (d_i + d_j)`, the

View File

@ -134,8 +134,8 @@ interactions.
The CHAu.ILP potential file provided with LAMMPS (see the potentials
directory) are parameterized for *metal* units. You can use this
potential with any LAMMPS units, but you would need to create your
BNCH.ILP potential file with coefficients listed in the appropriate
potential with any LAMMPS units, but you would need to create your own
custom CHAu.ILP potential file with coefficients listed in the appropriate
units, if your simulation does not use *metal* units.
Related commands

View File

@ -114,6 +114,7 @@ accelerated styles exist.
* :doc:`adp <pair_adp>` - angular dependent potential (ADP) of Mishin
* :doc:`agni <pair_agni>` - AGNI machine-learning potential
* :doc:`aip/water/2dm <pair_aip_water_2dm>` - anisotropic interfacial potential for water in 2d geometries
* :doc:`airebo <pair_airebo>` - AIREBO potential of Stuart
* :doc:`airebo/morse <pair_airebo>` - AIREBO with Morse instead of LJ
* :doc:`amoeba <pair_amoeba>` -

View File

@ -46,7 +46,7 @@ lines of output, the string can be enclosed in triple quotes, as in
the last example above. If the text string contains variables, they
will be evaluated and their current values printed.
.. versionadded:: TBD
.. versionadded:: 15Jun2023
support for vector style variables

View File

@ -1001,7 +1001,7 @@ map kind (atom, bond, angle, dihedral, or improper) and the second
argument is the label. The function returns the corresponding numeric
type or triggers an error if the queried label does not exist.
.. versionadded:: TBD
.. versionadded:: 15Jun2023
The is_typelabel(kind,label) function has the same arguments as
label2type(), but returns 1 if the type label has been assigned,
@ -1320,7 +1320,7 @@ Vectors" discussion above.
Vector Initialization
---------------------
.. versionadded:: TBD
.. versionadded:: 15Jun2023
*Vector*-style variables only can be initialized with a special
syntax, instead of using a formula. The syntax is a bracketed,
@ -1481,7 +1481,7 @@ commands
.. code-block:: LAMMPS
# delete_atoms random fraction 0.5 yes all NULL 49839
# run 0
# run 0 post no
variable t equal temp # this thermo keyword invokes a temperature compute
print "Temperature of system = $t"
run 1000