Merge branch 'develop' into kokkospod

.github/CODEOWNERS

@@ -38,6 +38,7 @@ src/ML-HDNNP/*        @singraber
 src/ML-IAP/*          @athomps
 src/ML-PACE/*         @yury-lysogorskiy
 src/ML-POD/*          @exapde
+src/ML-UF3/*          @monk-04
 src/MOFFF/*           @hheenen
 src/MOLFILE/*         @akohlmey
 src/NETCDF/*          @pastewka
@@ -72,6 +73,8 @@ src/MC/fix_sgcmc.*                  @athomps
 src/REAXFF/compute_reaxff_atom.*    @rbberger
 src/KOKKOS/compute_reaxff_atom_kokkos.*    @rbberger
 src/REPLICA/fix_pimd_langevin.*     @Yi-FanLi
+src/DPD-BASIC/pair_dpd_coul_slater_long.* @Eddy-Barraud
+src/GPU/pair_dpd_coul_slater_long.* @Eddy-Barraud
 
 # core LAMMPS classes
 src/lammps.*              @sjplimp

cmake/CMakeLists.txt

@@ -256,6 +256,7 @@ set(STANDARD_PACKAGES
   DRUDE
   EFF
   ELECTRODE
+  EXTRA-COMMAND
   EXTRA-COMPUTE
   EXTRA-DUMP
   EXTRA-FIX
@@ -285,6 +286,7 @@ set(STANDARD_PACKAGES
   ML-RANN
   ML-SNAP
   ML-POD
+  ML-UF3
   MOFFF
   MOLECULE
   MOLFILE
@@ -689,7 +691,7 @@ endif()
 # packages which selectively include variants based on enabled styles
 # e.g. accelerator packages
 ######################################################################
-foreach(PKG_WITH_INCL CORESHELL DPD-SMOOTH MC MISC PHONON QEQ OPENMP KOKKOS OPT INTEL GPU)
+foreach(PKG_WITH_INCL CORESHELL DPD-BASIC DPD-SMOOTH MC MISC PHONON QEQ OPENMP KOKKOS OPT INTEL GPU)
   if(PKG_${PKG_WITH_INCL})
     include(Packages/${PKG_WITH_INCL})
   endif()

cmake/Modules/Packages/DPD-BASIC.cmake

@@ -0,0 +1,9 @@
+# pair style dpd/coul/slater/long may only be installed if also KSPACE is installed
+if(NOT PKG_KSPACE)
+  get_property(LAMMPS_PAIR_HEADERS GLOBAL PROPERTY PAIR)
+  list(REMOVE_ITEM LAMMPS_PAIR_HEADERS ${LAMMPS_SOURCE_DIR}/DPD-BASIC/pair_dpd_coul_slater_long.h)
+  set_property(GLOBAL PROPERTY PAIR "${LAMMPS_PAIR_HEADERS}")
+  get_target_property(LAMMPS_SOURCES lammps SOURCES)
+  list(REMOVE_ITEM LAMMPS_SOURCES ${LAMMPS_SOURCE_DIR}/DPD-BASIC/pair_dpd_coul_slater_long.cpp)
+  set_property(TARGET lammps PROPERTY SOURCES "${LAMMPS_SOURCES}")
+endif()

cmake/Modules/Packages/PLUMED.cmake

@@ -1,5 +1,9 @@
 # Plumed2 support for PLUMED package
 
+# for supporting multiple concurrent plumed2 installations for debugging and testing
+set(PLUMED_SUFFIX "" CACHE STRING "Suffix for Plumed2 library")
+mark_as_advanced(PLUMED_SUFFIX)
+
 if(BUILD_MPI)
   set(PLUMED_CONFIG_MPI "--enable-mpi")
   set(PLUMED_CONFIG_CC  ${CMAKE_MPI_C_COMPILER})
@@ -21,9 +25,11 @@ else()
   set(PLUMED_CONFIG_OMP "--disable-openmp")
 endif()
 
-set(PLUMED_URL "https://github.com/plumed/plumed2/releases/download/v2.8.3/plumed-src-2.8.3.tgz"
+# Note: must also adjust check for supported API versions in
+# fix_plumed.cpp when version changes from v2.n.x to v2.n+1.y
+set(PLUMED_URL "https://github.com/plumed/plumed2/releases/download/v2.9.1/plumed-src-2.9.1.tgz"
   CACHE STRING "URL for PLUMED tarball")
-set(PLUMED_MD5 "76d23cd394eba9e6530316ed1184e219" CACHE STRING "MD5 checksum of PLUMED tarball")
+set(PLUMED_MD5 "c3b2d31479c1e9ce211719d40e9efbd7" CACHE STRING "MD5 checksum of PLUMED tarball")
 
 mark_as_advanced(PLUMED_URL)
 mark_as_advanced(PLUMED_MD5)
@@ -151,15 +157,15 @@ else()
     file(MAKE_DIRECTORY ${INSTALL_DIR}/include)
   else()
     find_package(PkgConfig REQUIRED)
-    pkg_check_modules(PLUMED REQUIRED plumed)
+    pkg_check_modules(PLUMED REQUIRED plumed${PLUMED_SUFFIX})
     add_library(LAMMPS::PLUMED INTERFACE IMPORTED)
     if(PLUMED_MODE STREQUAL "STATIC")
-      include(${PLUMED_LIBDIR}/plumed/src/lib/Plumed.cmake.static)
+      include(${PLUMED_LIBDIR}/plumed${PLUMED_SUFFIX}/src/lib/Plumed.cmake.static)
     elseif(PLUMED_MODE STREQUAL "SHARED")
-      include(${PLUMED_LIBDIR}/plumed/src/lib/Plumed.cmake.shared)
+      include(${PLUMED_LIBDIR}/plumed${PLUMED_SUFFIX}/src/lib/Plumed.cmake.shared)
     elseif(PLUMED_MODE STREQUAL "RUNTIME")
-      set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_COMPILE_DEFINITIONS "__PLUMED_DEFAULT_KERNEL=${PLUMED_LIBDIR}/${CMAKE_SHARED_LIBRARY_PREFIX}plumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
+      set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_COMPILE_DEFINITIONS "__PLUMED_DEFAULT_KERNEL=${PLUMED_LIBDIR}/${CMAKE_SHARED_LIBRARY_PREFIX}plumed${PLUMED_SUFFIX}Kernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
-      include(${PLUMED_LIBDIR}/plumed/src/lib/Plumed.cmake.runtime)
+      include(${PLUMED_LIBDIR}/plumed${PLUMED_SUFFIX}/src/lib/Plumed.cmake.runtime)
     endif()
     set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_LINK_LIBRARIES "${PLUMED_LOAD}")
     set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_INCLUDE_DIRECTORIES "${PLUMED_INCLUDE_DIRS}")

cmake/packaging/build_linux_tgz.sh

@@ -59,12 +59,14 @@ done
 
 echo "Set up wrapper script"
 MYDIR=$(dirname "$0")
+cp ${MYDIR}/xdg-open ${DESTDIR}/bin
 cp ${MYDIR}/linux_wrapper.sh ${DESTDIR}/bin
 for s in ${DESTDIR}/bin/*
 do \
         EXE=$(basename $s)
         test ${EXE} = linux_wrapper.sh && continue
         test ${EXE} = qt.conf && continue
+        test ${EXE} = xdg-open && continue
         ln -s bin/linux_wrapper.sh ${DESTDIR}/${EXE}
 done
 

cmake/packaging/linux_wrapper.sh

@@ -4,15 +4,17 @@
 # reset locale to avoid problems with decimal numbers
 export LC_ALL=C
 
-BASEDIR=$(dirname "$0")
+BASEDIR="$(dirname "$0")"
-EXENAME=$(basename "$0")
+EXENAME="$(basename "$0")"
+
+PATH="${BASEDIR}/bin:${PATH}"
 
 # append to LD_LIBRARY_PATH to prefer local (newer) libs
-LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${BASEDIR}/lib
+LD_LIBRARY_PATH="${LD_LIBRARY_PATH}:${BASEDIR}/lib"
 
 # set some environment variables for LAMMPS etc.
-LAMMPS_POTENTIALS=${BASEDIR}/share/lammps/potentials
+LAMMPS_POTENTIALS="${BASEDIR}/share/lammps/potentials"
-MSI2LMP_LIBRARY=${BASEDIR}/share/lammps/frc_files
+MSI2LMP_LIBRARY="${BASEDIR}/share/lammps/frc_files"
-export LD_LIBRARY_PATH LAMMPS_POTENTIALS MSI2LMP_LIBRARY
+export LD_LIBRARY_PATH LAMMPS_POTENTIALS MSI2LMP_LIBRARY PATH
 
 exec "${BASEDIR}/bin/${EXENAME}" "$@"

cmake/presets/all_off.cmake

@@ -28,6 +28,7 @@ set(ALL_PACKAGES
   DRUDE
   ELECTRODE
   EFF
+  EXTRA-COMMAND
   EXTRA-COMPUTE
   EXTRA-DUMP
   EXTRA-FIX
@@ -60,6 +61,7 @@ set(ALL_PACKAGES
   ML-QUIP
   ML-RANN
   ML-SNAP
+  ML-UF3
   MOFFF
   MOLECULE
   MOLFILE

cmake/presets/all_on.cmake

@@ -30,6 +30,7 @@ set(ALL_PACKAGES
   DRUDE
   ELECTRODE
   EFF
+  EXTRA-COMMAND
   EXTRA-COMPUTE
   EXTRA-DUMP
   EXTRA-FIX
@@ -62,6 +63,7 @@ set(ALL_PACKAGES
   ML-QUIP
   ML-RANN
   ML-SNAP
+  ML-UF3
   MOFFF
   MOLECULE
   MOLFILE

cmake/presets/mingw-cross.cmake

@@ -24,6 +24,7 @@ set(WIN_PACKAGES
   DRUDE
   ELECTRODE
   EFF
+  EXTRA-COMMAND
   EXTRA-COMPUTE
   EXTRA-DUMP
   EXTRA-FIX
@@ -50,6 +51,7 @@ set(WIN_PACKAGES
   ML-POD
   ML-RANN
   ML-SNAP
+  ML-UF3
   MOFFF
   MOLECULE
   MOLFILE

cmake/presets/most.cmake

@@ -26,6 +26,7 @@ set(ALL_PACKAGES
   DRUDE
   EFF
   ELECTRODE
+  EXTRA-COMMAND
   EXTRA-COMPUTE
   EXTRA-DUMP
   EXTRA-FIX
@@ -45,6 +46,7 @@ set(ALL_PACKAGES
   ML-IAP
   ML-POD
   ML-SNAP
+  ML-UF3
   MOFFF
   MOLECULE
   OPENMP

cmake/presets/windows.cmake

@@ -22,6 +22,7 @@ set(WIN_PACKAGES
   DRUDE
   EFF
   ELECTRODE
+  EXTRA-COMMAND
   EXTRA-COMPUTE
   EXTRA-DUMP
   EXTRA-FIX
@@ -42,6 +43,7 @@ set(WIN_PACKAGES
   ML-IAP
   ML-POD
   ML-SNAP
+  ML-UF3
   MOFFF
   MOLECULE
   MOLFILE

doc/src/Commands_bond.rst

@@ -27,7 +27,7 @@ OPT.
 
    * :doc:`none <bond_none>`
    * :doc:`zero <bond_zero>`
-   * :doc:`hybrid <bond_hybrid>`
+   * :doc:`hybrid (k) <bond_hybrid>`
    *
    *
    *

doc/src/Commands_pair.rst

@@ -25,16 +25,16 @@ OPT.
 
    * :doc:`none <pair_none>`
    * :doc:`zero <pair_zero>`
-   * :doc:`hybrid (k) <pair_hybrid>`
+   * :doc:`hybrid (ko) <pair_hybrid>`
-   * :doc:`hybrid/overlay (k) <pair_hybrid>`
+   * :doc:`hybrid/molecular (o) <pair_hybrid>`
-   * :doc:`hybrid/scaled <pair_hybrid>`
+   * :doc:`hybrid/overlay (ko) <pair_hybrid>`
+   * :doc:`hybrid/scaled (o) <pair_hybrid>`
    * :doc:`kim <pair_kim>`
    * :doc:`list <pair_list>`
    * :doc:`tracker <pair_tracker>`
    *
    *
    *
-   *
    * :doc:`adp (ko) <pair_adp>`
    * :doc:`agni (o) <pair_agni>`
    * :doc:`aip/water/2dm (t) <pair_aip_water_2dm>`
@@ -94,9 +94,10 @@ OPT.
    * :doc:`coul/wolf (ko) <pair_coul>`
    * :doc:`coul/wolf/cs <pair_cs>`
    * :doc:`dpd (giko) <pair_dpd>`
-   * :doc:`dpd/fdt <pair_dpd_fdt>`
+   * :doc:`dpd/coul/slater/long (g) <pair_dpd_coul_slater_long>`
    * :doc:`dpd/ext (ko) <pair_dpd_ext>`
    * :doc:`dpd/ext/tstat (ko) <pair_dpd_ext>`
+   * :doc:`dpd/fdt <pair_dpd_fdt>`
    * :doc:`dpd/fdt/energy (k) <pair_dpd_fdt>`
    * :doc:`dpd/tstat (gko) <pair_dpd>`
    * :doc:`dsmc <pair_dsmc>`
@@ -269,7 +270,7 @@ OPT.
    * :doc:`smd/ulsph <pair_smd_ulsph>`
    * :doc:`smtbq <pair_smtbq>`
    * :doc:`snap (ik) <pair_snap>`
-   * :doc:`soft (go) <pair_soft>`
+   * :doc:`soft (gko) <pair_soft>`
    * :doc:`sph/heatconduction (g) <pair_sph_heatconduction>`
    * :doc:`sph/idealgas <pair_sph_idealgas>`
    * :doc:`sph/lj (g) <pair_sph_lj>`
@@ -303,6 +304,7 @@ OPT.
    * :doc:`tip4p/long/soft (o) <pair_fep_soft>`
    * :doc:`tri/lj <pair_tri_lj>`
    * :doc:`ufm (got) <pair_ufm>`
+   * :doc:`uf3 (k) <pair_uf3>`
    * :doc:`vashishta (gko) <pair_vashishta>`
    * :doc:`vashishta/table (o) <pair_vashishta>`
    * :doc:`wf/cut <pair_wf_cut>`

doc/src/Commands_removed.rst

@@ -148,6 +148,14 @@ performance characteristics on NVIDIA GPUs. Both, the KOKKOS
 and the :ref:`GPU package <PKG-GPU>` are maintained
 and allow running LAMMPS with GPU acceleration.
 
+i-PI tool
+---------
+
+.. versionchanged:: TBD
+
+The i-PI tool has been removed from the LAMMPS distribution.  Instead,
+instructions to install i-PI from PyPi via pip are provided.
+
 restart2data tool
 -----------------
 

doc/src/Fortran.rst

@@ -305,6 +305,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS.
    :ftype extract_setting: function
    :f extract_global: :f:func:`extract_global`
    :ftype extract_global: function
+   :f map_atom: :f:func:`map_atom`
+   :ftype map_atom: function
    :f extract_atom: :f:func:`extract_atom`
    :ftype extract_atom: function
    :f extract_compute: :f:func:`extract_compute`

doc/src/Library_properties.rst

@@ -13,6 +13,7 @@ This section documents the following functions:
 - :cpp:func:`lammps_extract_setting`
 - :cpp:func:`lammps_extract_global_datatype`
 - :cpp:func:`lammps_extract_global`
+- :cpp:func:`lammps_map_atom`
 
 --------------------
 
@@ -120,3 +121,8 @@ subdomains and processors.
 .. doxygenfunction:: lammps_extract_global
    :project: progguide
 
+-----------------------
+
+.. doxygenfunction:: lammps_map_atom
+   :project: progguide
+

doc/src/Packages_details.rst

@@ -52,6 +52,7 @@ page gives those details.
    * :ref:`DRUDE <PKG-DRUDE>`
    * :ref:`EFF <PKG-EFF>`
    * :ref:`ELECTRODE <PKG-ELECTRODE>`
+   * :ref:`EXTRA-COMMAND <PKG-EXTRA-COMMAND>`
    * :ref:`EXTRA-COMPUTE <PKG-EXTRA-COMPUTE>`
    * :ref:`EXTRA-DUMP <PKG-EXTRA-DUMP>`
    * :ref:`EXTRA-FIX <PKG-EXTRA-FIX>`
@@ -84,6 +85,7 @@ page gives those details.
    * :ref:`ML-QUIP <PKG-ML-QUIP>`
    * :ref:`ML-RANN <PKG-ML-RANN>`
    * :ref:`ML-SNAP <PKG-ML-SNAP>`
+   * :ref:`ML-UF3 <PKG-ML-UF3>`
    * :ref:`MOFFF <PKG-MOFFF>`
    * :ref:`MOLECULE <PKG-MOLECULE>`
    * :ref:`MOLFILE <PKG-MOLFILE>`
@@ -403,6 +405,7 @@ and :ref:`ASPHERE <PKG-ASPHERE>` packages are installed.
 * :doc:`bond_style oxdna2/\* <bond_oxdna>`
 * :doc:`bond_style oxrna2/\* <bond_oxdna>`
 * :doc:`fix nve/dotc/langevin <fix_nve_dotc_langevin>`
+* examples/PACKAGES/cgdna
 
 ----------
 
@@ -676,7 +679,12 @@ DPD-BASIC package
 Pair styles for the basic dissipative particle dynamics (DPD) method
 and DPD thermostatting.
 
-**Author:** Kurt Smith (U Pittsburgh), Martin Svoboda, Martin Lisal (ICPF and UJEP)
+Pair style :doc:`dpd/coul/slater/long <pair_dpd_coul_slater_long>` also
+includes smeared charges for coulomb interactions and thus requires the
+:ref:`KSPACE <PKG-KSPACE>` package to be installed to handle the long-range
+Coulomb part of the interactions.
+
+**Authors:** Kurt Smith (U Pittsburgh), Martin Svoboda, Martin Lisal (ICPF and UJEP), Eddy Barraud (IFPEN)
 
 **Supporting info:**
 
@@ -685,6 +693,7 @@ and DPD thermostatting.
 * :doc:`pair_style dpd/tstat <pair_dpd>`
 * :doc:`pair_style dpd/ext <pair_dpd_ext>`
 * :doc:`pair_style dpd/ext/tstat <pair_dpd_ext>`
+* :doc:`pair_style dpd/coul/slater/long <pair_dpd_coul_slater_long>`
 * examples/PACKAGES/dpd-basic
 
 ----------
@@ -886,6 +895,22 @@ This package has :ref:`specific installation instructions <electrode>` on the
 
 ----------
 
+.. _PKG-EXTRA-COMMAND:
+
+EXTRA-COMMAND package
+---------------------
+
+**Contents:**
+
+Additional command styles that are less commonly used.
+
+**Supporting info:**
+
+* src/EXTRA-COMMAND: filenames -> commands
+* :doc:`general commands <Commands_all>`
+
+----------
+
 .. _PKG-EXTRA-COMPUTE:
 
 EXTRA-COMPUTE package
@@ -1925,6 +1950,31 @@ computes which analyze attributes of the potential.
 
 ----------
 
+.. _PKG-ML-UF3:
+
+ML-UF3 package
+--------------
+
+**Contents:**
+
+A pair style for the ultra-fast force field potentials (UF3). UF3 is a
+methodology for deriving a highly accurate classical potential which is
+fast to evaluate and is fitted to a large archives of quantum mechanical
+(DFT) data.  The use of b-spline basis set in UF3 enables the rapid
+evaluation of 2-body and 3-body interactions.
+
+**Authors:** Ajinkya C Hire (University of Florida),
+Hendrik Krass (University of Constance),
+Matthias Rupp (Luxembourg Institute of Science and Technology),
+Richard Hennig (University of Florida)
+
+**Supporting info:**
+
+* src/ML-UF3: filenames -> commands
+* :doc:`pair_style uf3 <pair_uf3>`
+* examples/uf3
+* https://github.com/uf3/uf3
+
 .. _PKG-MOFFF:
 
 MOFFF package

doc/src/Packages_list.rst

@@ -158,6 +158,11 @@ whether an extra library is needed to build and use the package:
      - :doc:`fix electrode/conp <fix_electrode>`
      - PACKAGES/electrode
      - no
+   * - :ref:`EXTRA-COMMAND <PKG-EXTRA-COMMAND>`
+     - additional command styles
+     - :doc:`general commands <Commands_all>`
+     - n/a
+     - no
    * - :ref:`EXTRA-COMPUTE <PKG-EXTRA-COMPUTE>`
      - additional compute styles
      - :doc:`compute <compute>`
@@ -318,6 +323,11 @@ whether an extra library is needed to build and use the package:
      - :doc:`pair_style snap <pair_snap>`
      - snap
      - no
+   * - :ref:`ML-UF3 <PKG-ML-UF3>`
+     - quantum-fitted ultra fast potentials
+     - :doc:`pair_style uf3 <pair_uf3>`
+     - PACKAGES/uf3
+     - no
    * - :ref:`MOFFF <PKG-MOFFF>`
      - styles for `MOF-FF <MOFplus_>`_ force field
      - :doc:`pair_style buck6d/coul/gauss <pair_buck6d_coul_gauss>`

doc/src/Tools.rst

@@ -90,7 +90,7 @@ Miscellaneous tools
 
    * :ref:`LAMMPS coding standards <coding_standard>`
    * :ref:`emacs <emacs>`
-   * :ref:`i-pi <ipi>`
+   * :ref:`i-PI <ipi>`
    * :ref:`kate <kate>`
    * :ref:`LAMMPS shell <lammps_shell>`
    * :ref:`LAMMPS GUI <lammps_gui>`
@@ -376,21 +376,40 @@ See README file in the tools/fep directory.
 
 .. _ipi:
 
-i-pi tool
+i-PI tool
 -------------------
 
-The tools/i-pi directory contains a version of the i-PI package, with
+.. versionchanged:: TBD
-all the LAMMPS-unrelated files removed.  It is provided so that it can
+
-be used with the :doc:`fix ipi <fix_ipi>` command to perform
+The tools/i-pi directory used to contain a bundled version of the i-PI
-path-integral molecular dynamics (PIMD).
+software package for use with LAMMPS.  This version, however, was
+removed in 06/2024.
 
 The i-PI package was created and is maintained by Michele Ceriotti,
 michele.ceriotti at gmail.com, to interface to a variety of molecular
 dynamics codes.
 
-See the tools/i-pi/manual.pdf file for an overview of i-PI, and the
+i-PI is now available via PyPi using the pip package manager at:
-:doc:`fix ipi <fix_ipi>` page for further details on running PIMD
+https://pypi.org/project/ipi/
-calculations with LAMMPS.
+
+Here are the commands to set up a virtual environment and install
+i-PI into it with all its dependencies.
+
+.. code-block:: sh
+
+   python -m venv ipienv
+   source ipienv/bin/activate
+   pip install --upgrade pip
+   pip install ipi
+
+To install the development version from GitHub, please use:
+
+.. code-block:: sh
+
+   pip install git+https://github.com/i-pi/i-pi.git
+
+For further information, please consult the [i-PI home
+page](https://ipi-code.org).
 
 ----------
 
@@ -709,8 +728,8 @@ CMake is required.
 The LAMMPS GUI has been successfully compiled and tested on:
 
 - Ubuntu Linux 20.04LTS x86_64 using GCC 9, Qt version 5.12
-- Fedora Linux 38 x86\_64 using GCC 13 and Clang 16, Qt version 5.15LTS
+- Fedora Linux 40 x86\_64 using GCC 14 and Clang 17, Qt version 5.15LTS
-- Fedora Linux 38 x86\_64 using GCC 13, Qt version 6.5LTS
+- Fedora Linux 40 x86\_64 using GCC 14, Qt version 6.5LTS
 - Apple macOS 12 (Monterey) and macOS 13 (Ventura) with Xcode on arm64 and x86\_64, Qt version 5.15LTS
 - Windows 10 and 11 x86_64 with Visual Studio 2022 and Visual C++ 14.36, Qt version 5.15LTS
 - Windows 10 and 11 x86_64 with MinGW / GCC 10.0 cross-compiler on Fedora 38, Qt version 5.15LTS
@@ -752,22 +771,23 @@ if necessary.  When both, Qt5 and Qt6 are available, Qt6 will be preferred
 unless ``-D LAMMPS_GUI_USE_QT5=yes`` is set.
 
 It should be possible to build the LAMMPS GUI as a standalone
-compilation (e.g. when LAMMPS has been compiled with traditional make),
+compilation (e.g. when LAMMPS has been compiled with traditional make).
-then the CMake configuration needs to be told where to find the LAMMPS
+Then the CMake configuration needs to be told where to find the LAMMPS
 headers and the LAMMPS library, via ``-D
 LAMMPS_SOURCE_DIR=/path/to/lammps/src``.  CMake will try to guess a
 build folder with the LAMMPS library from that path, but it can also be
 set with ``-D LAMMPS_LIB_DIR=/path/to/lammps/lib``.
 
 Rather than linking to the LAMMPS library during compilation, it is also
-possible to compile the GUI with a plugin loader library that will load
+possible to compile the GUI with a plugin loader that will load
 the LAMMPS library dynamically at runtime during the start of the GUI
 from a shared library; e.g. ``liblammps.so`` or ``liblammps.dylib`` or
 ``liblammps.dll`` (depending on the operating system).  This has the
-advantage that the LAMMPS library can be updated LAMMPS without having
+advantage that the LAMMPS library can be built from updated or modified
-to recompile the GUI.  The ABI of the LAMMPS C-library interface is very
+LAMMPS source without having to recompile the GUI.  The ABI of the
-stable and generally backward compatible.  This feature is enabled by
+LAMMPS C-library interface is very stable and generally backward
-setting ``-D LAMMPS_GUI_USE_PLUGIN=on`` and then ``-D
+compatible.  This feature is enabled by setting
+``-D LAMMPS_GUI_USE_PLUGIN=on`` and then ``-D
 LAMMPS_PLUGINLIB_DIR=/path/to/lammps/plugin/loader``. Typically, this
 would be the ``examples/COUPLE/plugin`` folder of the LAMMPS
 distribution.
@@ -779,8 +799,8 @@ macOS
 """""
 
 When building on macOS, the build procedure will try to manufacture a
-drag-n-drop installer, LAMMPS-macOS-multiarch.dmg, when using the 'dmg'
+drag-n-drop installer, ``LAMMPS-macOS-multiarch.dmg``, when using the
-target (i.e. ``cmake --build <build dir> --target dmg`` or ``make dmg``.
+'dmg' target (i.e. ``cmake --build <build dir> --target dmg`` or ``make dmg``.
 
 To build multi-arch executables that will run on both, arm64 and x86_64
 architectures natively, it is necessary to set the CMake variable ``-D
@@ -820,11 +840,11 @@ and LAMMPS GUI can be launched from anywhere from the command line.
 The standard CMake build procedure can be applied and the
 ``mingw-cross.cmake`` preset used. By using ``mingw64-cmake`` the CMake
 command will automatically include a suitable CMake toolset file (the
-regular cmake command can be used after that).  After building the
+regular cmake command can be used after that to modify the configuration
-libraries and executables, you can build the target 'zip'
+settings, if needed).  After building the libraries and executables,
-(i.e. ``cmake --build <build dir> --target zip`` or ``make zip``
+you can build the target 'zip' (i.e. ``cmake --build <build dir> --target zip``
-to stage all installed files into a LAMMPS_GUI folder and then
+or ``make zip`` to stage all installed files into a LAMMPS_GUI folder
-run a script to copy all required dependencies, some other files,
+and then run a script to copy all required dependencies, some other files,
 and create a zip file from it.
 
 Linux

doc/src/bond_hybrid.rst

@@ -1,8 +1,11 @@
 .. index:: bond_style hybrid
+.. index:: bond_style hybrid/kk
 
 bond_style hybrid command
 =========================
 
+Accelerator Variants: *hybrid/kk*
+
 Syntax
 """"""
 
@@ -15,7 +18,7 @@ Syntax
 Examples
 """"""""
 
-.. code-block: LAMMPS
+.. code-block:: LAMMPS
 
    bond_style hybrid harmonic fene
    bond_coeff 1 harmonic 80.0 1.2
@@ -60,6 +63,10 @@ bond types.
 
 ----------
 
+.. include:: accel_styles.rst
+
+----------
+
 Restrictions
 """"""""""""
 

doc/src/bond_oxdna.rst

@@ -27,6 +27,7 @@ Examples
 
 .. code-block:: LAMMPS
 
+   # LJ units
    bond_style oxdna/fene
    bond_coeff * 2.0 0.25 0.7525
 
@@ -36,6 +37,32 @@ Examples
    bond_style oxrna2/fene
    bond_coeff * 2.0 0.25 0.76107
 
+   bond_style oxdna/fene
+   bond_coeff * oxdna_lj.cgdna
+
+   # Real units
+   bond_style oxdna/fene
+   bond_coeff * 11.92337812042065 2.1295 6.409795
+
+   bond_style oxdna2/fene
+   bond_coeff * 11.92337812042065 2.1295 6.4430152
+
+   bond_style oxrna2/fene
+   bond_coeff * 11.92337812042065 2.1295 6.482800913
+
+   bond_style oxrna2/fene
+   bond_coeff * oxrna2_real.cgdna
+
+.. note::
+
+   The coefficients in the above examples have to be kept fixed and
+   cannot be changed without reparameterizing the entire model. They are
+   provided in forms compatible with both *units lj* and *units real*
+   (see documentation of :doc:`units <units>`).  These can also be read
+   from a potential file with correct unit style by specifying the name
+   of the file. Several potential files for each unit style are included
+   in the ``potentials`` directory of the LAMMPS distribution.
+
 Description
 """""""""""
 
@@ -46,15 +73,14 @@ The *oxdna/fene*, *oxdna2/fene*, and *oxrna2/fene* bond styles use the potential
    E = - \frac{\epsilon}{2} \ln \left[ 1 - \left(\frac{r-r_0}{\Delta}\right)^2\right]
 
 to define a modified finite extensible nonlinear elastic (FENE)
-potential :ref:`(Ouldridge) <Ouldridge0>` to model the connectivity of the
+potential :ref:`(Ouldridge) <Ouldridge0>` to model the connectivity of
-phosphate backbone in the oxDNA/oxRNA force field for coarse-grained
+the phosphate backbone in the oxDNA/oxRNA force field for coarse-grained
 modelling of DNA/RNA.
 
 The following coefficients must be defined for the bond type via the
 :doc:`bond_coeff <bond_coeff>` command as given in the above example, or
-in the data file or restart files read by the
+in the data file or restart files read by the :doc:`read_data
-:doc:`read_data <read_data>` or :doc:`read_restart <read_restart>`
+<read_data>` or :doc:`read_restart <read_restart>` commands:
-commands:
 
 * :math:`\epsilon` (energy)
 * :math:`\Delta` (distance)
@@ -62,41 +88,40 @@ commands:
 
 .. note::
 
-   The oxDNA bond style has to be used together with the
+   The oxDNA bond style has to be used together with the corresponding
-   corresponding oxDNA pair styles for excluded volume interaction
+   oxDNA pair styles for excluded volume interaction *oxdna/excv* ,
-   *oxdna/excv* , stacking *oxdna/stk* , cross-stacking *oxdna/xstk* and
+   stacking *oxdna/stk* , cross-stacking *oxdna/xstk* and coaxial
-   coaxial stacking interaction *oxdna/coaxstk* as well as
+   stacking interaction *oxdna/coaxstk* as well as hydrogen-bonding
-   hydrogen-bonding interaction *oxdna/hbond* (see also documentation of
+   interaction *oxdna/hbond* (see also documentation of :doc:`pair_style
-   :doc:`pair_style oxdna/excv <pair_oxdna>`). For the oxDNA2
+   oxdna/excv <pair_oxdna>`). For the oxDNA2 :ref:`(Snodin) <Snodin0>`
-   :ref:`(Snodin) <Snodin0>` bond style the analogous pair styles
+   bond style the analogous pair styles *oxdna2/excv* , *oxdna2/stk* ,
-   *oxdna2/excv* , *oxdna2/stk* , *oxdna2/xstk* , *oxdna2/coaxstk* ,
+   *oxdna2/xstk* , *oxdna2/coaxstk* , *oxdna2/hbond* and an additional
-   *oxdna2/hbond* and an additional Debye-Hueckel pair style
+   Debye-Hueckel pair style *oxdna2/dh* have to be defined. The same
-   *oxdna2/dh* have to be defined. The same applies to the oxRNA2
+   applies to the oxRNA2 :ref:`(Sulc1) <Sulc01>` styles.
-   :ref:`(Sulc1) <Sulc01>` styles.
-   The coefficients in the above example have to be kept fixed and cannot
-   be changed without reparameterizing the entire model.
 
 .. note::
 
-   This bond style has to be used with the *atom_style hybrid bond ellipsoid oxdna*
+   This bond style has to be used with the *atom_style hybrid bond
-   (see documentation of :doc:`atom_style <atom_style>`). The *atom_style oxdna*
+   ellipsoid oxdna* (see documentation of :doc:`atom_style
-   stores the 3'-to-5' polarity of the nucleotide strand, which is set through
+   <atom_style>`). The *atom_style oxdna* stores the 3'-to-5' polarity
-   the bond topology in the data file. The first (second) atom in a bond definition
+   of the nucleotide strand, which is set through the bond topology in
-   is understood to point towards the 3'-end (5'-end) of the strand.
+   the data file. The first (second) atom in a bond definition is
+   understood to point towards the 3'-end (5'-end) of the strand.
 
 .. warning::
 
-   If data files are produced with :doc:`write_data <write_data>`, then the
+   If data files are produced with :doc:`write_data <write_data>`, then
-   :doc:`newton <newton>` command should be set to *newton on* or *newton off on*.
+   the :doc:`newton <newton>` command should be set to *newton on* or
-   Otherwise the data files will not have the same 3'-to-5' polarity as the
+   *newton off on*.  Otherwise the data files will not have the same
-   initial data file. This limitation does not apply to binary restart files
+   3'-to-5' polarity as the initial data file. This limitation does not
-   produced with :doc:`write_restart <write_restart>`.
+   apply to binary restart files produced with :doc:`write_restart
+   <write_restart>`.
 
 Example input and data files for DNA and RNA duplexes can be found in
-examples/PACKAGES/cgdna/examples/oxDNA/ , /oxDNA2/ and /oxRNA2/.  A simple python
+``examples/PACKAGES/cgdna/examples/oxDNA/`, `.../oxDNA2/`` and
-setup tool which creates single straight or helical DNA strands, DNA/RNA
+``.../oxRNA2/``.  A simple python setup tool which creates single
-duplexes or arrays of DNA/RNA duplexes can be found in
+straight or helical DNA strands, DNA/RNA duplexes or arrays of DNA/RNA
-examples/PACKAGES/cgdna/util/.
+duplexes can be found in ``examples/PACKAGES/cgdna/util/``.
 
 Please cite :ref:`(Henrich) <Henrich0>` in any publication that uses
 this implementation. An updated documentation that contains general information
@@ -113,6 +138,39 @@ and for sequence-specific hydrogen-bonding and stacking interactions
 
 ----------
 
+Potential file reading
+""""""""""""""""""""""
+
+For each style oxdna, oxdna2 and oxrna2, the first parameter argument
+can be a filename, and if it is, no further arguments should be
+supplied. Therefore the following command:
+
+.. code-block:: LAMMPS
+
+   bond_style oxdna/fene
+   bond_coeff * oxdna_lj.cgdna
+
+will be interpreted as a request to read the (FENE) potential
+:ref:`(Ouldridge) <Ouldridge0>` parameters from the file with the given
+name.  The file can define multiple potential parameters for both bonded
+and pair interactions, but for the above bonded interactions there must
+exist in the file a line of the form:
+
+.. code-block:: LAMMPS
+
+   *   fene    epsilon delta r0
+
+There are sample potential files for each unit style in the
+``potentials`` directory of the LAMMPS distribution. The potential file
+unit system must align with the units defined via the :doc:`units
+<units>` command. For conversion between different *LJ* and *real* unit
+systems for oxDNA, the python tool *lj2real.py* located in the
+``examples/PACKAGES/cgdna/util/`` directory can be used. This tool
+assumes similar file structure to the examples found in
+``examples/PACKAGES/cgdna/examples/``.
+
+----------
+
 Restrictions
 """"""""""""
 

doc/src/compute_rdf.rst

@@ -13,8 +13,8 @@ Syntax
 * ID, group-ID are documented in :doc:`compute <compute>` command
 * rdf = style name of this compute command
 * Nbin = number of RDF bins
-* itypeN = central atom type for Nth RDF histogram (see asterisk form below)
+* itypeN = central atom type for Nth RDF histogram (integer, type label, or asterisk form)
-* jtypeN = distribution atom type for Nth RDF histogram (see asterisk form below)
+* jtypeN = distribution atom type for Nth RDF histogram (integer, type label, or asterisk form)
 * zero or more keyword/value pairs may be appended
 * keyword = *cutoff*
 
@@ -96,14 +96,16 @@ is computed for :math:`g(r)` between all atom types.  If one or more pairs are
 listed, then a separate histogram is generated for each
 *itype*,\ *jtype* pair.
 
-The *itypeN* and *jtypeN* settings can be specified in one of two
+The *itypeN* and *jtypeN* settings can be specified in one of three
-ways.  An explicit numeric value can be used, as in the fourth example
+ways.  One or both of the types in the I,J pair can be a
-above.  Or a wild-card asterisk can be used to specify a range of atom
+:doc:`type label <Howto_type_labels>`.  Or an explicit numeric value can be
-types.  This takes the form "\*" or "\*n" or "m\*" or "m\*n".  If
+used, as in the fourth example above. Or a wild-card asterisk can be used
-:math:`N` is the number of atom types, then an asterisk with no numeric values
+to specify a range of atom types.  This takes the form "\*" or "\*n" or
-means all types from 1 to :math:`N`. A leading asterisk means all types from 1
+"m\*" or "m\*n". If :math:`N` is the number of atom types, then an asterisk
-to n (inclusive).  A trailing asterisk means all types from m to :math:`N`
+with no numeric values means all types from 1 to :math:`N`. A leading
-(inclusive).  A middle asterisk means all types from m to n (inclusive).
+asterisk means all types from 1 to n (inclusive).  A trailing asterisk
+means all types from m to :math:`N` (inclusive).  A middle asterisk means
+all types from m to n (inclusive).
 
 If both *itypeN* and *jtypeN* are single values, as in the fourth example
 above, this means that a :math:`g(r)` is computed where atoms of type *itypeN*

doc/src/create_atoms.rst

@@ -10,7 +10,7 @@ Syntax
 
    create_atoms type style args keyword values ...
 
-* type = atom type (1-Ntypes) of atoms to create (offset for molecule creation)
+* type = atom type (1-Ntypes or type label) of atoms to create (offset for molecule creation)
 * style = *box* or *region* or *single* or *mesh* or *random*
 
   .. parsed-literal::
@@ -37,7 +37,7 @@ Syntax
          seed = random # seed (positive integer)
        *basis* values = M itype
          M = which basis atom
-         itype = atom type (1-N) to assign to this basis atom
+         itype = atom type (1-Ntypes or type label) to assign to this basis atom
        *ratio* values = frac seed
          frac = fraction of lattice sites (0 to 1) to populate randomly
          seed = random # seed (positive integer)
@@ -74,6 +74,13 @@ Examples
 .. code-block:: LAMMPS
 
    create_atoms 1 box
+
+   labelmap atom 1 Pt
+   create_atoms Pt box
+
+   labelmap atom 1 C 2 Si
+   create_atoms C region regsphere basis Si C
+
    create_atoms 3 region regsphere basis 2 3
    create_atoms 3 region regsphere basis 2 3 ratio 0.5 74637
    create_atoms 3 single 0 0 5

doc/src/delete_bonds.rst

@@ -43,6 +43,9 @@ Examples
    delete_bonds all bond 0*3 special
    delete_bonds all stats
 
+   labelmap atom 4 hc
+   delete_bonds all atom hc special
+
 Description
 """""""""""
 
@@ -59,19 +62,20 @@ For all styles, by default, an interaction is only turned off (or on)
 if all the atoms involved are in the specified group.  See the *any*
 keyword to change the behavior.
 
-Several of the styles (\ *atom*, *bond*, *angle*, *dihedral*,
+Several of the styles (\ *atom*, *bond*, *angle*, *dihedral*, *improper*\ )
-*improper*\ ) take a *type* as an argument.  The specified *type* should
+take a *type* as an argument.  The specified *type* can be a
-be an integer from 0 to :math:`N`, where :math:`N` is the number of relevant
+:doc:`type label <Howto_type_labels>`.  Otherwise, the type should be an
+integer from 0 to :math:`N`, where :math:`N` is the number of relevant
 types (atom types, bond types, etc.).  A value of 0 is only relevant for
-style *bond*\ ; see details below.  In all cases, a wildcard asterisk
+style *bond*\ ; see details below.  For numeric types, a wildcard asterisk
 can be used in place of or in conjunction with the *type* argument to
 specify a range of types.  This takes the form "\*" or "\*n" or "m\*" or
-"m\*n".  If :math:`N` is the number of types, then an asterisk with no numeric
+"m\*n".  If :math:`N` is the number of types, then an asterisk with no
-values means all types from 0 to :math:`N`.  A leading asterisk means all
+numeric values means all types from 0 to :math:`N`.  A leading asterisk
-types from 0 to n (inclusive).  A trailing asterisk means all types
+means all types from 0 to n (inclusive).  A trailing asterisk means all
-from m to N (inclusive).  A middle asterisk means all types from m to
+types from m to N (inclusive).  A middle asterisk means all types from m to
-n (inclusive).  Note that it is fine to include a type of 0 for
+n (inclusive).  Note that it is fine to include a type of 0 for non-bond
-non-bond styles; it will simply be ignored.
+styles; it will simply be ignored.
 
 For style *multi* all bond, angle, dihedral, and improper interactions
 of any type, involving atoms in the group, are turned off.

doc/src/dump.rst

@@ -114,6 +114,7 @@ Syntax
            proc = ID of processor that owns atom
            procp1 = ID+1 of processor that owns atom
            type = atom type
+           typelabel = atom :doc:`type label <Howto_type_labels>`
            element = name of atom element, as defined by :doc:`dump_modify <dump_modify>` command
            mass = atom mass
            x,y,z = unscaled atom coordinates
@@ -470,8 +471,9 @@ followed by one line per atom with the atom type and the :math:`x`-,
 :math:`y`-, and :math:`z`-coordinate of that atom.  You can use the
 :doc:`dump_modify element <dump_modify>` option to change the output
 from using the (numerical) atom type to an element name (or some other
-label). This will help many visualization programs to guess bonds and
+label). This option will help many visualization programs to guess bonds
-colors.
+and colors. You can use the :doc:`dump_modify types labels <dump_modify>`
+option to replace numeric atom types with :doc:`type labels <Howto_type_labels>`.
 
 .. versionadded:: 22Dec2022
 
@@ -774,21 +776,21 @@ command creates a per-atom array with six columns:
 
 Per-atom attributes used as arguments to the *custom* and *cfg* styles:
 
-The *id*, *mol*, *proc*, *procp1*, *type*, *element*, *mass*, *vx*,
+The *id*, *mol*, *proc*, *procp1*, *type*, *typelabel*, *element*, *mass*,
-*vy*, *vz*, *fx*, *fy*, *fz*, *q* attributes are self-explanatory.
+*vx*, *vy*, *vz*, *fx*, *fy*, *fz*, *q* attributes are self-explanatory.
 
-*Id* is the atom ID.  *Mol* is the molecule ID, included in the data
+*Id* is the atom ID.  *Mol* is the molecule ID, included in the data file
-file for molecular systems.  *Proc* is the ID of the processor (0 to
+for molecular systems.  *Proc* is the ID of the processor (0 to
 :math:`N_\text{procs}-1`) that currently owns the atom.  *Procp1* is the
 proc ID+1, which can be convenient in place of a *type* attribute (1 to
 :math:`N_\text{types}`) for coloring atoms in a visualization program.
-*Type* is the atom type (1 to :math:`N_\text{types}`).  *Element* is
+*Type* is the atom type (1 to :math:`N_\text{types}`).  *Typelabel* is the
-typically the chemical name of an element, which you must assign to each
+atom :doc:`type label <Howto_type_labels>`.  *Element* is typically the
-type via the :doc:`dump_modify element <dump_modify>` command.  More
+chemical name of an element, which you must assign to each type via the
-generally, it can be any string you wish to associated with an atom
+:doc:`dump_modify element <dump_modify>` command.  More generally, it can
-type.  *Mass* is the atom mass. The quantities *vx*, *vy*, *vz*, *fx*,
+be any string you wish to associated with an atom type.  *Mass* is the atom
-*fy*, *fz*, and *q* are components of atom velocity and force and atomic
+mass.  The quantities *vx*, *vy*, *vz*, *fx*, *fy*, *fz*, and *q* are
-charge.
+components of atom velocity and force and atomic charge.
 
 There are several options for outputting atom coordinates.  The *x*,
 *y*, and *z* attributes write atom coordinates "unscaled", in the

doc/src/dump_modify.rst

@@ -17,7 +17,7 @@ Syntax
 * one or more keyword/value pairs may be appended
 
 * these keywords apply to various dump styles
-* keyword = *append* or *at* or *balance* or *buffer* or *colname* or *delay* or *element* or *every* or *every/time* or *fileper* or *first* or *flush* or *format* or *header* or *image* or *label* or *maxfiles* or *nfile* or *pad* or *pbc* or *precision* or *region* or *refresh* or *scale* or *sfactor* or *skip* or *sort* or *tfactor* or *thermo* or *thresh* or *time* or *triclinic/general* or *units* or *unwrap*
+* keyword = *append* or *at* or *balance* or *buffer* or *colname* or *delay* or *element* or *every* or *every/time* or *fileper* or *first* or *flush* or *format* or *header* or *image* or *label* or *maxfiles* or *nfile* or *pad* or *pbc* or *precision* or *region* or *refresh* or *scale* or *sfactor* or *skip* or *sort* or *tfactor* or *thermo* or *thresh* or *time* or *triclinic/general* or *types* or *units* or *unwrap*
 
   .. parsed-literal::
 
@@ -81,6 +81,7 @@ Syntax
          these 3 args can be replaced by the word "none" to turn off thresholding
        *time* arg = *yes* or *no*
        *triclinic/general* arg = *yes* or *no*
+       *types* value = *numeric* or *labels*
        *units* arg = *yes* or *no*
        *unwrap* arg = *yes* or *no*
 
@@ -849,6 +850,13 @@ The default setting is *no*\ .
 
 ----------
 
+The *types* keyword applies only to the dump xyz style. If this keyword is
+used with a value of *numeric*, then numeric atom types are printed in the
+xyz file (default). If the value *labels* is specified, then
+:doc:`type labels <Howto_type_labels>` are printed for atom types.
+
+----------
+
 The *triclinic/general* keyword only applies to the dump *atom* and
 *custom* styles.  It can only be used with a value of *yes* if the
 simulation box was created as a general triclinic box.  See the
@@ -960,11 +968,11 @@ The option defaults are
 * sort = id for dump styles *dcd*, *xtc*, and *xyz*
 * thresh = none
 * time = no
-* triclinic/general no
+* triclinic/general = no
+* types = numeric
 * units = no
 * unwrap = no
 
 * compression_level = 9 (gz variants)
 * compression_level = 0 (zstd variants)
 * checksum = yes (zstd variants)
-

doc/src/fix_adapt.rst

@@ -21,17 +21,17 @@ Syntax
        *pair* args = pstyle pparam I J v_name
          pstyle = pair style name (e.g., lj/cut)
          pparam = parameter to adapt over time
-         I,J = type pair(s) to set parameter for
+         I,J = type pair(s) to set parameter for (integer or type label)
          v_name = variable with name that calculates value of pparam
        *bond* args = bstyle bparam I v_name
          bstyle = bond style name (e.g., harmonic)
          bparam = parameter to adapt over time
-         I = type bond to set parameter for
+         I = type bond to set parameter for (integer or type label)
          v_name = variable with name that calculates value of bparam
        *angle* args = astyle aparam I v_name
          astyle = angle style name (e.g., harmonic)
          aparam = parameter to adapt over time
-         I = type angle to set parameter for
+         I = type angle to set parameter for (integer or type label)
          v_name = variable with name that calculates value of aparam
        *kspace* arg = v_name
          v_name = variable with name that calculates scale factor on :math:`k`-space terms
@@ -67,6 +67,9 @@ Examples
    variable ramp_up equal "ramp(0.01,0.5)"
    fix stretch all adapt 1 bond harmonic r0 1 v_ramp_up
 
+   labelmap atom 1 c1
+   fix 1 all adapt 1 pair soft a c1 c1 v_prefactor
+
 Description
 """""""""""
 
@@ -254,10 +257,12 @@ should be specified to indicate which type pairs to apply it to.  If a global
 parameter is specified, the :math:`I` and :math:`J` settings still need to be
 specified, but are ignored.
 
-Similar to the :doc:`pair_coeff command <pair_coeff>`, :math:`I` and :math:`J`
+Similar to the :doc:`pair_coeff command <pair_coeff>`, :math:`I` and
-can be specified in one of two ways.  Explicit numeric values can be used for
+:math:`J` can be specified in one of several ways.  Explicit numeric values
-each, as in the first example above.  :math:`I \le J` is required.  LAMMPS sets
+can be used for each, as in the first example above.  Or, one or both of
-the coefficients for the symmetric :math:`J,I` interaction to the same values.
+the types in the I,J pair can be a :doc:`type label <Howto_type_labels>`.
+LAMMPS sets the coefficients for the symmetric :math:`J,I` interaction to
+the same values.
 
 A wild-card asterisk can be used in place of or in conjunction with
 the :math:`I,J` arguments to set the coefficients for multiple pairs of atom
@@ -266,8 +271,9 @@ 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 m to :math:`N`
 (inclusive).  A middle asterisk means all types from m to n
-(inclusive).  Note that only type pairs with :math:`I \le J` are considered; if
+(inclusive).  For the asterisk syntax, note that only type pairs with
-asterisks imply type pairs where :math:`J < I`, they are ignored.
+:math:`I \le J` are considered; if asterisks imply type pairs where
+:math:`J < I`, they are ignored.
 
 IMPORTANT NOTE: If :doc:`pair_style hybrid or hybrid/overlay
 <pair_hybrid>` is being used, then the *pstyle* will be a sub-style

doc/src/fix_adapt_fep.rst

@@ -21,13 +21,13 @@ Syntax
        *pair* args = pstyle pparam I J v_name
          pstyle = pair style name (e.g., lj/cut)
          pparam = parameter to adapt over time
-         I,J = type pair(s) to set parameter for
+         I,J = type pair(s) to set parameter for (integer or type label)
          v_name = variable with name that calculates value of pparam
        *kspace* arg = v_name
          v_name = variable with name that calculates scale factor on K-space terms
        *atom* args = aparam v_name
          aparam = parameter to adapt over time
-         I = type(s) to set parameter for
+         I = type(s) to set parameter for (integer or type label)
          v_name = variable with name that calculates value of aparam
 
 * zero or more keyword/value pairs may be appended
@@ -56,6 +56,9 @@ Examples
    fix 1 all adapt/fep 1 pair lj/cut epsilon * * v_scale1 coul/cut scale 3 3 v_scale2 scale yes reset yes
    fix 1 all adapt/fep 10 atom diameter 1 v_size
 
+   labelmap atom 1 c1
+   fix 1 all adapt/fep 1 pair soft a c1 c1 v_prefactor
+
 
 Example input scripts available: examples/PACKAGES/fep
 
@@ -218,10 +221,12 @@ be specified to indicate which type pairs to apply it to.  If a global
 parameter is specified, the *I* and *J* settings still need to be
 specified, but are ignored.
 
-Similar to the :doc:`pair_coeff command <pair_coeff>`, I and J can be
+Similar to the :doc:`pair_coeff command <pair_coeff>`, :math:`I` and
-specified in one of two ways.  Explicit numeric values can be used for
+:math:`J` can be specified in one of several ways.  Explicit numeric values
-each, as in the first example above.  :math:`I \le J` is required.  LAMMPS sets
+can be used for each, as in the first example above.  Or, one or both of
-the coefficients for the symmetric J,I interaction to the same values.
+the types in the I,J pair can be a :doc:`type label <Howto_type_labels>`.
+LAMMPS sets the coefficients for the symmetric :math:`J,I` interaction to
+the same values.
 
 A wild-card asterisk can be used in place of or in conjunction with
 the :math:`I,J` arguments to set the coefficients for multiple pairs of atom
@@ -230,8 +235,9 @@ 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 m to :math:`N`
 (inclusive).  A middle asterisk means all types from m to n
-(inclusive).  Note that only type pairs with :math:`I \le J` are considered; if
+(inclusive).  For the asterisk syntax, note that only type pairs with
-asterisks imply type pairs where :math:`J < I`, they are ignored.
+:math:`I \le J` are considered; if asterisks imply type pairs where
+:math:`J < I`, they are ignored.
 
 IMPROTANT NOTE: If :doc:`pair_style hybrid or hybrid/overlay <pair_hybrid>` is
 being used, then the *pstyle* will be a sub-style name.  You must specify

doc/src/fix_amoeba_bitorsion.rst

@@ -35,7 +35,11 @@ the implementation 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.
+important to properly represent their conformational behavior.  Each
+bitorsion interaction is thus defined for a 5-tuple of atoms
+:math:`IJKLM` with bonds between successive atoms in the list,
+i.e. two overlapping dihedral interactions for atoms :math:`IJKL` and
+:math:`JKLM`.
 
 The examples/amoeba directory has a sample input script and data file
 for ubiquitin, which illustrates use of the fix amoeba/bitorsion
@@ -68,14 +72,15 @@ lines:
           [...]
           N       3     314     315     317      318    330
 
-The first column is an index from 1 to :math:`N` to enumerate the bitorsion
+The first column is an index from 1 to :math:`N` to enumerate the
-5-atom tuples; it is ignored by LAMMPS.  The second column is the
+bitorsion 5-atom tuples; it is ignored by LAMMPS.  The second column
-*type* of the interaction; it is an index into the bitorsion force
+is the *type* of the interaction; it is an index into the bitorsion
-field file.  The remaining 5 columns are the atom IDs of the atoms in
+force field file.  The remaining 5 columns are the atom IDs of the
-the two 4-atom dihedrals that overlap to create the bitorsion 5-body
+atoms (in order) for the 5-tuple :math:`IJKLM`, as described above.
-interaction.  Note that the *bitorsions* and *BiTorsions* keywords for
+
-the header and body sections match those specified in the
+Note that the *bitorsions* and *BiTorsions* keywords for the header
-:doc:`read_data <read_data>` command following the data file name.
+and body sections match those specified in the :doc:`read_data
+<read_data>` command following the data file name.
 
 The data file should be generated by using the
 tools/tinker/tinker2lmp.py conversion script which creates a LAMMPS

doc/src/fix_amoeba_pitorsion.rst

@@ -57,7 +57,7 @@ should have two lines like these in its header section:
    M pitorsion types
    N pitorsions
 
-where :math:`N` is the number of pitorsion 5-body interactions and :math:`M` is
+where :math:`N` is the number of pitorsion 6-body interactions and :math:`M` is
 the number of pitorsion types.  It should also have two sections in the body
 of the data file like these with :math:`M` and :math:`N` lines each:
 
@@ -74,21 +74,20 @@ of the data file like these with :math:`M` and :math:`N` lines each:
 
    PiTorsions
 
-          1       1       8      10      12      18      20
+          1       1       2      4      3      20      21      24
-          2       5      18      20      22      25      27
+          2       5      21     23     22      37      38      41
           [...]
-          N       3     314     315     317      318    330
+          N       7      27     29     28      30      35      36
 
-For PiTorsion Coeffs, the first column is an index from 1 to :math:`M` to
+For PiTorsion Coeffs, the first column is an index from 1 to :math:`M`
-enumerate the pitorsion types.  The second column is the single
+to enumerate the pitorsion types.  The second column is the single
 prefactor coefficient needed for each type.
 
-For PiTorsions, the first column is an index from 1 to :math:`N` to enumerate
+For PiTorsions, the first column is an index from 1 to :math:`N` to
-the pitorsion 5-atom tuples; it is ignored by LAMMPS.  The second
+enumerate the pitorsion 6-atom tuples; it is ignored by LAMMPS.  The
-column is the "type" of the interaction; it is an index into the
+second column is the "type" of the interaction; it is an index into
-PiTorsion Coeffs.  The remaining 5 columns are the atom IDs of the
+the PiTorsion Coeffs.  The remaining 6 columns are the atom IDs of the
-atoms in the two 4-atom dihedrals that overlap to create the pitorsion
+atoms (in order) for the 6-tuple :math:`IJKLMN`, as described above.
-5-body interaction.
 
 Note that the *pitorsion types* and *pitorsions* and *PiTorsion
 Coeffs* and *PiTorsions* keywords for the header and body sections of

doc/src/fix_atom_swap.rst

@@ -21,7 +21,7 @@ Syntax
 
   .. parsed-literal::
 
-       *types* values = two or more atom types
+       *types* values = two or more atom types (1-Ntypes or type label)
        *mu* values = chemical potential of swap types (energy units)
        *ke* value = *no* or *yes*
          *no* = no conservation of kinetic energy after atom swaps

doc/src/fix_bond_break.rst

@@ -13,7 +13,7 @@ Syntax
 * ID, group-ID are documented in :doc:`fix <fix>` command
 * bond/break = style name of this fix command
 * Nevery = attempt bond breaking every this many steps
-* bondtype = type of bonds to break
+* bondtype = type of bonds to break (integer or type label)
 * Rmax = bond longer than Rmax can break (distance units)
 * zero or more keyword/value pairs may be appended
 * keyword = *prob*

doc/src/fix_bond_create.rst

@@ -17,9 +17,9 @@ Syntax
 * ID, group-ID are documented in :doc:`fix <fix>` command
 * bond/create = style name of this fix command
 * Nevery = attempt bond creation every this many steps
-* itype,jtype = atoms of itype can bond to atoms of jtype
+* itype,jtype = atoms of itype can bond to atoms of jtype (1-Ntypes or type label)
 * Rmin = 2 atoms separated by less than Rmin can bond (distance units)
-* bondtype = type of created bonds
+* bondtype = type of created bonds (integer or type label)
 * zero or more keyword/value pairs may be appended to args
 * keyword = *iparam* or *jparam* or *prob* or *atype* or *dtype* or *itype* or *aconstrain*
 
@@ -27,19 +27,19 @@ Syntax
 
        *iparam* values = maxbond, newtype
          maxbond = max # of bonds of bondtype the itype atom can have
-         newtype = change the itype atom to this type when maxbonds exist
+         newtype = change the itype atom to this type when maxbonds exist (1-Ntypes or type label)
        *jparam* values = maxbond, newtype
          maxbond = max # of bonds of bondtype the jtype atom can have
-         newtype = change the jtype atom to this type when maxbonds exist
+         newtype = change the jtype atom to this type when maxbonds exist (1-Ntypes or type label)
        *prob* values = fraction seed
          fraction = create a bond with this probability if otherwise eligible
          seed = random number seed (positive integer)
        *atype* value = angletype
-         angletype = type of created angles
+         angletype = type of created angles (integer or type label)
        *dtype* value = dihedraltype
-         dihedraltype = type of created dihedrals
+         dihedraltype = type of created dihedrals (integer or type label)
        *itype* value = impropertype
-         impropertype = type of created impropers
+         impropertype = type of created impropers (integer or type label)
        *aconstrain* value = amin amax
          amin = minimal angle at which new bonds can be created
          amax = maximal angle at which new bonds can be created
@@ -54,6 +54,10 @@ Examples
    fix 5 all bond/create 1 3 3 0.8 1 prob 0.5 85784 iparam 2 3 atype 1 dtype 2
    fix 5 all bond/create/angle 10 1 2 1.122 1 aconstrain 120 180 prob 1 4928459 iparam 2 1 jparam 2 2
 
+   labelmap atom 1 c1 2 n2
+   labelmap bond 1 c1-n2
+   fix 5 all bond/create 10 c1 n2 0.8 c1-n2
+
 Description
 """""""""""
 

doc/src/fix_charge_regulation.rst

@@ -13,8 +13,8 @@ Syntax
 
 * ID, group-ID are documented in fix command
 * charge/regulation = style name of this fix command
-* cation_type = atom type of free cations
+* cation_type = atom type of free cations (integer or type label)
-* anion_type = atom type of free anions
+* anion_type = atom type of free anions (integer or type label)
 
 * zero or more keyword/value pairs may be appended
 
@@ -27,8 +27,8 @@ Syntax
      *pIp* value = activity (effective concentration) of free cations (in the -log10 representation)
      *pIm* value = activity (effective concentration) of free anions (in the -log10 representation)
      *pKs* value = solvent self-dissociation constant (in the -log10 representation)
-     *acid_type* = atom type of acid groups
+     *acid_type* = atom type of acid groups (integer or type label)
-     *base_type*  = atom type of base groups
+     *base_type*  = atom type of base groups (integer or type label)
      *lunit_nm* value = unit length used by LAMMPS (# in the units of nanometers)
      *temp* value = temperature
      *tempfixid* value = fix ID of temperature thermostat
@@ -51,6 +51,9 @@ Examples
     fix chareg all charge/regulation 1 2 acid_type 3 base_type 4 pKa 5.0 pKb 6.0 pH 7.0 pIp 3.0 pIm 3.0 nevery 200 nmc 200 seed 123 tempfixid fT
     fix chareg all charge/regulation 1 2 pIp 3 pIm 3 onlysalt yes 2 -1 seed 123 tag yes temp 1.0
 
+    labelmap atom 1 H+ 2 OH-
+    fix chareg all charge/regulation H+ OH- pIp 3 pIm 3 onlysalt yes 2 -1 seed 123 tag yes temp 1.0
+
 Description
 """""""""""
 

doc/src/fix_deform.rst

@@ -64,6 +64,8 @@ Syntax
                  effectively an engineering shear strain rate
            *erate* value = R
              R = engineering shear strain rate (1/time units)
+           *erate/rescale* value = R (ONLY available in :doc:`fix deform/pressure <fix_deform_pressure>` command)
+             R = engineering shear strain rate (1/time units)
            *trate* value = R
              R = true shear strain rate (1/time units)
            *wiggle* values = A Tp

doc/src/fix_deform_pressure.rst

@@ -311,6 +311,10 @@ This fix is not invoked during :doc:`energy minimization <minimize>`.
 Restrictions
 """"""""""""
 
+This fix is part of the EXTRA-FIX package.  It is only enabled if LAMMPS
+was built with that package.  See the :doc:`Build package <Build_package>`
+page for more info.
+
 You cannot apply x, y, or z deformations to a dimension that is
 shrink-wrapped via the :doc:`boundary <boundary>` command.
 

doc/src/fix_deposit.rst

@@ -13,7 +13,7 @@ Syntax
 * ID, group-ID are documented in :doc:`fix <fix>` command
 * deposit = style name of this fix command
 * N = # of atoms or molecules to insert
-* type = atom type to assign to inserted atoms (offset for molecule insertion)
+* type = atom type (1-Ntypes or type label) to assign to inserted atoms (offset for molecule insertion)
 * M = insert a single atom or molecule every M steps
 * seed = random # seed (positive integer)
 * one or more keyword/value pairs may be appended to args
@@ -76,6 +76,9 @@ Examples
    fix 4 sputter deposit 1000 2 500 12235 region sphere vz -1.0 -1.0 target 5.0 5.0 0.0 units lattice
    fix 5 insert deposit 200 2 100 777 region disk gaussian 5.0 5.0 9.0 1.0 units box
 
+   labelmap atom 1 Au
+   fix 4 sputter deposit 1000 Au 500 12235 region sphere vz -1.0 -1.0 target 5.0 5.0 0.0 units lattice
+
 Description
 """""""""""
 

doc/src/fix_gcmc.rst

@@ -15,7 +15,7 @@ Syntax
 * N = invoke this fix every N steps
 * X = average number of GCMC exchanges to attempt every N steps
 * M = average number of MC moves to attempt every N steps
-* type = atom type for inserted atoms (must be 0 if mol keyword used)
+* type = atom type (1-Ntypes or type label) for inserted atoms (must be 0 if mol keyword used)
 * seed = random # seed (positive integer)
 * T = temperature of the ideal gas reservoir (temperature units)
 * mu = chemical potential of the ideal gas reservoir (energy units)
@@ -45,7 +45,7 @@ Syntax
        *group* value = group-ID
          group-ID = group-ID for inserted atoms (string)
        *grouptype* values = type group-ID
-         type = atom type (int)
+         type = atom type (1-Ntypes or type label)
          group-ID = group-ID for inserted atoms (string)
        *intra_energy* value = intramolecular energy (energy units)
        *tfac_insert* value = scale up/down temperature of inserted atoms (unitless)
@@ -62,52 +62,47 @@ Examples
    fix 3 water gcmc 10 100 100 0 3456543 3.0 -2.5 0.1 mol my_one_water maxangle 180 full_energy
    fix 4 my_gas gcmc 1 10 10 1 123456543 300.0 -12.5 1.0 region disk
 
+   labelmap atom 1 Li
+   fix 2 ion gcmc 10 1000 1000 Li 29494 298.0 -0.5 0.01
+
 Description
 """""""""""
 
-This fix performs grand canonical Monte Carlo (GCMC) exchanges of
+This fix performs grand canonical Monte Carlo (GCMC) exchanges of atoms or
-atoms or molecules with an imaginary ideal gas
+molecules with an imaginary ideal gas reservoir at the specified T and
-reservoir at the specified T and chemical potential (mu) as discussed
+chemical potential (mu) as discussed in :ref:`(Frenkel) <Frenkel2>`.  It
-in :ref:`(Frenkel) <Frenkel2>`. It also
+also attempts Monte Carlo (MC) moves (translations and molecule rotations)
-attempts  Monte Carlo (MC) moves (translations and molecule
+within the simulation cell or region.  If used with the :doc:`fix nvt <fix_nh>`
-rotations) within the simulation cell or
-region. If used with the :doc:`fix nvt <fix_nh>`
 command, simulations in the grand canonical ensemble (muVT, constant
 chemical potential, constant volume, and constant temperature) can be
 performed.  Specific uses include computing isotherms in microporous
 materials, or computing vapor-liquid coexistence curves.
 
-Every N timesteps the fix attempts both GCMC exchanges
+Every N timesteps the fix attempts both GCMC exchanges (insertions or
-(insertions or deletions) and MC moves of gas atoms or molecules.
+deletions) and MC moves of gas atoms or molecules.  On those timesteps, the
-On those timesteps, the average number of attempted GCMC exchanges is X,
+average number of attempted GCMC exchanges is X, while the average number
-while the average number of attempted MC moves is M.
+of attempted MC moves is M.  For GCMC exchanges of either molecular or
-For GCMC exchanges of either molecular or atomic gasses,
+atomic gasses, these exchanges can be either deletions or insertions, with
-these exchanges can be either deletions or insertions,
+equal probability.
-with equal probability.
 
-The possible choices for MC moves are translation of an atom,
+The possible choices for MC moves are translation of an atom, translation
-translation of a molecule, and rotation of a molecule.
+of a molecule, and rotation of a molecule.  The relative amounts of each are
-The relative amounts of each are determined by the optional
+determined by the optional *mcmoves* keyword (see below).  The default
-*mcmoves* keyword (see below).
+behavior is as follows. If the *mol* keyword is used, only molecule
-The default behavior is as follows.
+translations and molecule rotations are performed with equal probability.
-If the *mol* keyword is used, only molecule translations
+Conversely, if the *mol* keyword is not used, only atom translations are
-and molecule rotations are performed with equal probability.
+performed.  M should typically be chosen to be approximately equal to the
-Conversely, if the *mol* keyword is not used, only atom
+expected number of gas atoms or molecules of the given type within the
-translations are performed.
+simulation cell or region, which will result in roughly one MC move per
-M should typically be
+atom or molecule per MC cycle.
-chosen to be approximately equal to the expected number of gas atoms
-or molecules of the given type within the simulation cell or region,
-which will result in roughly one MC move per atom or molecule
-per MC cycle.
 
-All inserted particles are always added to two groups: the default
+All inserted particles are always added to two groups: the default group
-group "all" and the fix group specified in the fix command.
+"all" and the fix group specified in the fix command.  In addition,
-In addition, particles are also added to any groups
+particles are also added to any groups specified by the *group* and
-specified by the *group* and *grouptype* keywords.  If inserted
+*grouptype* keywords.  If inserted particles are individual atoms, they are
-particles are individual atoms, they are assigned the atom type given
+assigned the atom type given by the type argument.  If they are molecules,
-by the type argument.  If they are molecules, the type argument has no
+the type argument has no effect and must be set to zero. Instead, the type
-effect and must be set to zero. Instead, the type of each atom in the
+of each atom in the inserted molecule is specified in the file read by the
-inserted molecule is specified in the file read by the
 :doc:`molecule <molecule>` command.
 
 .. note::

doc/src/fix_ipi.rst

@@ -35,23 +35,24 @@ Description
 """""""""""
 
 This fix enables LAMMPS to be run as a client for the i-PI Python
-wrapper :ref:`(IPI) <ipihome>` for performing a path integral molecular dynamics
+wrapper :ref:`(IPI) <ipihome>`. i-PI is a universal force engine,
-(PIMD) simulation.  The philosophy behind i-PI is described in the
+designed to perform advanced molecular simulations, with a special
-following publication :ref:`(IPI-CPC) <IPICPC>`.
+focus on path integral molecular dynamics (PIMD) simulation.
+The philosophy behind i-PI is to separate the evaluation of the
+energy and forces, which is delegated to the client, and the evolution
+of the dynamics, that is the responsibility of i-PI. This approach also
+simplifies combining energies computed from different codes, which
+can for instance be useful to mix first-principles calculations,
+empirical force fields or machine-learning potentials.
+The following publication :ref:`(IPI-CPC-2014) <IPICPC>` discusses the
+overall implementation of i-PI, and focuses on path-integral techniques,
+while a later release :ref:`(IPI-CPC-2019) <IPICPC2>` introduces several
+additional features and simulation schemes.
 
-A version of the i-PI package, containing only files needed for use
+The communication between i-PI and LAMMPS takes place using sockets,
-with LAMMPS, is provided in the tools/i-pi directory.  See the
+and is reduced to the bare minimum. All the parameters of the dynamics
-tools/i-pi/manual.pdf for an introduction to i-PI.  The
+are specified in the input of i-PI, and all the parameters of the force
-examples/PACKAGES/i-pi directory contains example scripts for using i-PI
+field must be specified as LAMMPS inputs, preceding the *fix ipi* command.
-with LAMMPS.
-
-In brief, the path integral molecular dynamics is performed by the
-Python wrapper, while the client (LAMMPS in this case) simply computes
-forces and energy for each configuration. The communication between
-the two components takes place using sockets, and is reduced to the
-bare minimum. All the parameters of the dynamics are specified in the
-input of i-PI, and all the parameters of the force field must be
-specified as LAMMPS inputs, preceding the *fix ipi* command.
 
 The server address must be specified by the *address* argument, and
 can be either the IP address, the fully-qualified name of the server,
@@ -75,6 +76,20 @@ If the cell varies too wildly, it may be advisable to re-initialize
 these interactions at each call. This behavior can be requested by
 setting the *reset* switch.
 
+Obtaining i-PI
+""""""""""""""
+
+Here are the commands to set up a virtual environment and install
+i-PI into it with all its dependencies via the PyPi repository and
+the pip package manager.
+
+.. code-block:: sh
+
+   python -m venv ipienv
+   source ipienv/bin/activate
+   pip install --upgrade pip
+   pip install ipi
+
 Restart, fix_modify, output, run start/stop, minimize info
 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
 
@@ -111,9 +126,14 @@ Related commands
 
 .. _IPICPC:
 
-**(IPI-CPC)** Ceriotti, More and Manolopoulos, Comp Phys Comm, 185,
+**(IPI-CPC-2014)** Ceriotti, More and Manolopoulos, Comp Phys Comm 185,
 1019-1026 (2014).
 
+.. _IPICPC2:
+
+**(IPI-CPC-2019)** Kapil et al., Comp Phys Comm 236, 214-223 (2019).
+
+
 .. _ipihome:
 
 **(IPI)**

doc/src/fix_mol_swap.rst

@@ -14,7 +14,7 @@ Syntax
 * atom/swap = style name of this fix command
 * N = invoke this fix every N steps
 * X = number of swaps to attempt every N steps
-* itype,jtype = two atom types to swap with each other
+* itype,jtype = two atom types (1-Ntypes or type label) to swap with each other
 * seed = random # seed (positive integer)
 * T = scaling temperature of the MC swaps (temperature units)
 * zero or more keyword/value pairs may be appended to args
@@ -34,6 +34,9 @@ Examples
    fix 2 all mol/swap 100 1 2 3 29494 300.0 ke no
    fix mySwap fluid mol/swap 500 10 1 2 482798 1.0
 
+   labelmap atom 1 A 2 B
+   fix mySwap fluid mol/swap 500 10 A B 482798 1.0
+
 Description
 """""""""""
 

doc/src/fix_reaxff_species.rst

@@ -134,36 +134,34 @@ value.  For example, AuO.pos.\* becomes AuO.pos.0, AuO.pos.1000, etc.
 
 .. versionadded:: 3Aug2022
 
-The optional keyword *delete* enables the periodic removal of
+The optional keyword *delete* enables the periodic removal of molecules
-molecules from the system.  Criteria for deletion can be either a list
+from the system :ref:`(Gissinger) <Delete>`.  Criteria for deletion can
-of specific chemical formulae or a range of molecular weights.
+be either a list of specific chemical formulae or a range of molecular
-Molecules are deleted every *Nfreq* timesteps, and bond connectivity
+weights.  Molecules are deleted every *Nfreq* timesteps, and bond
-is determined using the *Nevery* and *Nrepeat* keywords.  The
+connectivity is determined using the *Nevery* and *Nrepeat* keywords.  The
-*filedel* argument is the name of the output file that records the
+*filedel* argument is the name of the output file that records the species
-species that are removed from the system.  The *specieslist* keyword
+that are removed from the system.  The *specieslist* keyword permits
-permits specific chemical species to be deleted.  The *Nspecies*
+specific chemical species to be deleted.  The *Nspecies* argument specifies
-argument specifies how many species are eligible for deletion and is
+how many species are eligible for deletion and is followed by a list of
-followed by a list of chemical formulae, whose strings are compared to
+chemical formulae, whose strings are compared to species identified by this
-species identified by this fix.  For example, "specieslist 2 CO CO2"
+fix.  For example, "specieslist 2 CO CO2" deletes molecules that are
-deletes molecules that are identified as "CO" and "CO2" in the species
+identified as "CO" and "CO2" in the species output file.  When using the
-output file.  When using the *specieslist* keyword, the *filedel* file
+*specieslist* keyword, the *filedel* file has the following format: the
-has the following format: the first line lists the chemical formulae
+first line lists the chemical formulae eligible for deletion, and each
-eligible for deletion, and each additional line contains the timestep
+additional line contains the timestep on which a molecule deletion occurs
-on which a molecule deletion occurs and the number of each species
+and the number of each species deleted on that timestep.  The *masslimit*
-deleted on that timestep.  The *masslimit* keyword permits deletion of
+keyword permits deletion of molecules with molecular weights between
-molecules with molecular weights between *massmin* and *massmax*.
+*massmin* and *massmax*.  When using the *masslimit* keyword, each line of
-When using the *masslimit* keyword, each line of the *filedel* file
+the *filedel* file contains the timestep on which deletions occurs,
-contains the timestep on which deletions occurs, followed by how many
+followed by how many of each species are deleted (with quantities preceding
-of each species are deleted (with quantities preceding chemical
+chemical formulae).  The *specieslist* and *masslimit* keywords cannot both
-formulae).  The *specieslist* and *masslimit* keywords cannot both be
+be used in the same *reaxff/species* fix.  The *delete_rate_limit* keyword
-used in the same *reaxff/species* fix.  The *delete_rate_limit*
+can enforce an upper limit on the overall rate of molecule deletion.  The
-keyword can enforce an upper limit on the overall rate of molecule
+number of deletion occurrences is limited to Nlimit within an interval of
-deletion.  The number of deletion occurrences is limited to Nlimit
+Nsteps timesteps.  Nlimit can be specified with an equal-style
-within an interval of Nsteps timesteps.   Nlimit can be specified with
+:doc:`variable <variable>`.  When using the *delete_rate_limit* keyword, no
-an equal-style :doc:`variable <variable>`.  When using the
+deletions are permitted to occur within the first Nsteps timesteps of the
-*delete_rate_limit* keyword, no deletions are permitted to occur
+first run (after reading a either a data or restart file).
-within the first Nsteps timesteps of the first run (after reading a
-either a data or restart file).
 
 ----------
 
@@ -235,3 +233,7 @@ Default
 The default values for bond-order cutoffs are 0.3 for all I-J pairs.
 The default element symbols are taken from the ReaxFF pair_coeff command.
 Position files are not written by default.
+
+.. _Delete:
+
+**(Gissinger)** Jacob R. Gissinger, Scott R. Zavada, Joseph G. Smith, Josh Kemppainen, Ivan Gallegos, Gregory M. Odegard, Emilie J. Siochi, and Kristopher E. Wise, Carbon, 202, 336-347 (2023).

doc/src/fix_sgcmc.rst

@@ -15,7 +15,7 @@ Syntax
 * every_nsteps = number of MD steps between MC cycles
 * swap_fraction = fraction of a full MC cycle carried out at each call (a value of 1.0 will perform as many trial moves as there are atoms)
 * temperature = temperature that enters Boltzmann factor in Metropolis criterion (usually the same as MD temperature)
-* deltamu = chemical potential difference(s) (`N-1` values must be provided, with `N` being the number of elements)
+* deltamu = `N-1` chemical potential differences :math:`\mu_1-\mu_2, \ldots, \mu_1-\mu_N` (`N` is the number of atom types)
 * Zero or more keyword/value pairs may be appended to fix definition line:
 
   .. parsed-literal::
@@ -23,7 +23,7 @@ Syntax
      keyword = *variance* or *randseed* or *window_moves* or *window_size*
        *variance* kappa conc1 [conc2] ... [concN]
          kappa = variance constraint parameter
-         conc1,conc2,...  = target concentration(s) in the range 0.0-1.0 (*N-1* values must be provided, with *N* being the number of elements)
+         `c_2`, `c_3`,..., `c_N` = `N-1` target concentration fractions
        *randseed* N
          N = seed for pseudo random number generator
        *window_moves* N
@@ -90,11 +90,10 @@ the simulation, e.g., to speed up equilibration at low temperatures.
 
 ------------
 
-The parameter *deltamu* is used to set the chemical potential difference
+The parameter *deltamu* is used to set the chemical potential differences
-in the SGC MC algorithm (see Eq. 16 in :ref:`Sadigh1 <Sadigh1>`). By
+in the SGC MC algorithm (see Eq. 16 in :ref:`Sadigh1 <Sadigh1>`).
-convention it is the difference of the chemical potentials of elements
+The `N-1` differences are defined as :math:`\mu_1-\mu_2, \ldots, \mu_1-\mu_N`,
-`B`, `C` ..., with respect to element A. When the simulation includes
+where `N` is the number of atom types.
-`N` elements, `N-1` values must be specified.
 
 ------------
 
@@ -104,12 +103,12 @@ the effective average constraint in the parallel VC-SGC MC algorithm
 (parameter :math:`\delta\mu_0` in Eq. (20) of :ref:`Sadigh1
 <Sadigh1>`). The parameter *kappa* specifies the variance constraint
 (see Eqs. (20-21) in :ref:`Sadigh1 <Sadigh1>`).
-
+The parameter *conc* sets the `N-1` target atomic concentration
-The parameter *conc* sets the target concentration (parameter
+fractions (parameter :math:`c_0` in Eqs.  (20-21) of :ref:`Sadigh1 <Sadigh1>`)
-:math:`c_0` in Eqs.  (20-21) of :ref:`Sadigh1 <Sadigh1>`). The atomic
+:math:`0 \le c_2, \ldots, c_N \le 1`, with
-concentrations refer to components `B`, `C` ..., with `A` being set
+:math:`c_1 = 1 - \Sigma_{i=2}^N c_i`.
-automatically. When the simulation includes `N` elements, `N-1`
+When the simulation includes `N` atom types (elements),
-concentration values must be specified.
+`N-1` concentration values must be specified.
 
 ------------
 
@@ -143,10 +142,10 @@ components of the vector represent the following quantities:
 
 * 1 = The absolute number of accepted trial swaps during the last MC step
 * 2 = The absolute number of rejected trial swaps during the last MC step
-* 3 = The current global concentration of species *A* (= number of atoms of type 1 / total number of atoms)
+* 3 = Current global concentration `c_1` (= number of atoms of type 1 / total number of atoms)
-* 4 = The current global concentration of species *B* (= number of atoms of type 2 / total number of atoms)
+* 4 = Current global concentration `c_2` (= number of atoms of type 2 / total number of atoms)
 * ...
-* N+2: The current global concentration of species *X* (= number of atoms of type *N* / total number of atoms)
+* N+2: Current global concentration `c_N` (= number of atoms of type *N* / total number of atoms)
 
 The vector values calculated by this fix are "intensive".
 

doc/src/fix_wall_gran.rst

@@ -115,6 +115,18 @@ friction and twisting friction supported by the :doc:`pair_style granular <pair_
 supported for walls. These are discussed in greater detail on the doc
 page for :doc:`pair_style granular <pair_granular>`.
 
+.. note::
+   When *fstyle* *granular* is specified, the associated *fstyle_params* are taken as
+   those for a wall/particle interaction. For example, for the *hertz/material* normal
+   contact model with :math:`E = 960` and :math:`\nu = 0.2`, the effective Young's
+   modulus for a wall/particle interaction is computed as
+   :math:`E_{eff} = \frac{960}{2(1-0.2^2)} = 500`. Any pair coefficients defined by
+   :doc:`pair_style granular <pair_granular>` are not taken into consideration. To
+   model different wall/particle interactions for particles of different material
+   types, the user may define multiple fix wall/gran commands operating on separate
+   groups (e.g. based on particle type) each with a different wall/particle effective
+   Young's modulus.
+
 Note that you can choose a different force styles and/or different
 values for the wall/particle coefficients than for particle/particle
 interactions.  E.g. if you wish to model the wall as a different

doc/src/fix_widom.rst

@@ -14,7 +14,7 @@ Syntax
 * widom = style name of this fix command
 * N = invoke this fix every N steps
 * M = number of Widom insertions to attempt every N steps
-* type = atom type for inserted atoms (must be 0 if mol keyword used)
+* type = atom type (1-Ntypes or type label) for inserted atoms (must be 0 if mol keyword used)
 * seed = random # seed (positive integer)
 * T = temperature of the system (temperature units)
 * zero or more keyword/value pairs may be appended to args
@@ -38,6 +38,9 @@ Examples
    fix 2 gas widom 1 50000 1 19494 2.0
    fix 3 water widom 1000 100 0 29494 300.0 mol h2omol full_energy
 
+   labelmap atom 1 Li
+   fix 2 ion widom 1 50000 Li 19494 2.0
+
 Description
 """""""""""
 

doc/src/group.rst

@@ -20,13 +20,13 @@ Syntax
        *empty* = no args
        *region* args = region-ID
        *type* or *id* or *molecule*
-         args = list of one or more atom types, atom IDs, or molecule IDs
+         args = list of one or more atom types (1-Ntypes or type label), atom IDs, or molecule IDs
-           any entry in list can be a sequence formatted as A:B or A:B:C where
+           any numeric entry in list can be a sequence formatted as A:B or A:B:C where
            A = starting index, B = ending index,
            C = increment between indices, 1 if not specified
          args = logical value
            logical = "<" or "<=" or ">" or ">=" or "==" or "!="
-           value = an atom type or atom ID or molecule ID (depending on *style*\ )
+           value = an atom type (1-Ntypes or type label) or atom ID or molecule ID (depending on *style*\ )
          args = logical value1 value2
            logical = "<>"
            value1,value2 = atom types or atom IDs or molecule IDs (depending on *style*\ )
@@ -52,6 +52,7 @@ Examples
 
    group edge region regstrip
    group water type 3 4
+   group water type OW HT
    group sub id 10 25 50
    group sub id 10 25 50 500:1000
    group sub id 100:10000:10
@@ -119,7 +120,7 @@ three styles can use arguments specified in one of two formats.
 
 The first format is a list of values (types or IDs).  For example, the
 second command in the examples above puts all atoms of type 3 or 4 into
-the group named *water*\ .  Each entry in the list can be a
+the group named *water*\ .  Each numeric entry in the list can be a
 colon-separated sequence ``A:B`` or ``A:B:C``, as in two of the examples
 above.  A "sequence" generates a sequence of values (types or IDs),
 with an optional increment.  The first example with ``500:1000`` has the
@@ -135,7 +136,8 @@ except ``<>`` take a single argument.  The third example above adds all
 atoms with IDs from 1 to 150 to the group named *sub*\ .  The logical ``<>``
 means "between" and takes 2 arguments.  The fourth example above adds all
 atoms belonging to molecules with IDs from 50 to 250 (inclusive) to
-the group named polyA.
+the group named polyA.  For the *type* style, type labels are converted into
+numeric types before being evaluated.
 
 The *variable* style evaluates a variable to determine which atoms to
 add to the group.  It must be an :doc:`atom-style variable <variable>`

doc/src/group2ndx.rst

@@ -34,32 +34,66 @@ Description
 Write or read a Gromacs style index file in text format that associates
 atom IDs with the corresponding group definitions. This index file can be
 used with in combination with Gromacs analysis tools or to import group
-definitions into the :doc:`fix colvars <fix_colvars>` input file. It can
+definitions into the :doc:`fix colvars <fix_colvars>` input file.
-also be used to save and restore group definitions for static groups.
+
+It can also be used to save and restore group definitions for static groups
+using the individual atom IDs. This may be important if the original
+group definition depends on a region or otherwise on the geometry and thus
+cannot be easily recreated.
+
+Another application would be to import atom groups defined for Gromacs
+simulation into LAMMPS.  When translating Gromacs topology and geometry
+data to LAMMPS.
 
 The *group2ndx* command will write group definitions to an index file.
-Without specifying any group IDs, all groups will be written to the index
+Without specifying any group IDs, all groups will be written to the
-file. When specifying group IDs, only those groups will be written to the
+index file.  When specifying group IDs, only those groups will be
-index file. In order to follow the Gromacs conventions, the group *all*
+written to the index file.  In order to follow the Gromacs conventions,
-will be renamed to *System* in the index file.
+the group *all* will be renamed to *System* in the index file.
 
-The *ndx2group* command will create of update group definitions from those
+The *ndx2group* command will create of update group definitions from
-stored in an index file. Without specifying any group IDs, all groups except
+those stored in an index file.  Without specifying any group IDs, all
-*System* will be read from the index file and the corresponding groups
+groups except *System* will be read from the index file and the
-recreated. If a group of the same name already exists, it will be completely
+corresponding groups recreated.  If a group of the same name already
-reset. When specifying group IDs, those groups, if present, will be read
+exists, it will be completely reset.  When specifying group IDs, those
-from the index file and restored.
+groups, if present, will be read from the index file and restored.
+
+File Format
+"""""""""""
+
+The file format is equivalent and compatible with what is produced by
+the `Gromacs make_ndx command <https://manual.gromacs.org/current/onlinehelp/gmx-make_ndx.html>`_.
+and follows the `Gromacs definition of an ndx file <https://manual.gromacs.org/current/reference-manual/file-formats.html#ndx>`_
+
+Each group definition begins with the group name in square brackets with
+blanks, e.g. \[ water \] and is then followed by the list of atom
+indices, which may be spread over multiple lines.  Here is a small
+example file:
+
+.. code-block:: ini
+
+   [ Oxygen ]
+   1  4  7
+   [ Hydrogen ]
+   2  3  5  6
+   8  9
+   [ Water ]
+   1 2 3 4 5 6 7 8 9
+
+The index file defines 3 groups: Oxygen, Hydrogen, and Water and the
+latter happens to be the union of the first two.
 
 ----------
 
 Restrictions
 """"""""""""
 
-This command requires that atoms have atom IDs, since this is the
+These commands require that atoms have atom IDs, since this is the
 information that is written to the index file.
 
-These commands are part of the COLVARS package.  They are only
+These commands are part of the EXTRA-COMMAND package.  They are only
-enabled if LAMMPS was built with that package.  See the :doc:`Build package <Build_package>` page for more info.
+enabled if LAMMPS was built with that package.  See the
+:doc:`Build package <Build_package>` page for more info.
 
 Related commands
 """"""""""""""""

doc/src/improper_cossq.rst

@@ -64,8 +64,8 @@ Restrictions
 """"""""""""
 
 This improper style can only be used if LAMMPS was built with the
-MOLECULE package.  See the :doc:`Build package <Build_package>` doc
+EXTRA-MOLECULE package.  See the :doc:`Build package <Build_package>`
-page for more info.
+doc page for more info.
 
 Related commands
 """"""""""""""""

doc/src/improper_distance.rst

@@ -54,8 +54,8 @@ Restrictions
 """"""""""""
 
 This improper style can only be used if LAMMPS was built with the
-MOLECULE package.  See the :doc:`Build package <Build_package>` doc
+EXTRA-MOLECULE package.  See the :doc:`Build package <Build_package>`
-page for more info.
+doc page for more info.
 
 Related commands
 """"""""""""""""

doc/src/improper_fourier.rst

@@ -60,8 +60,8 @@ Restrictions
 """"""""""""
 
 This angle style can only be used if LAMMPS was built with the
-MOLECULE package.  See the :doc:`Build package <Build_package>` doc
+EXTRA-MOLECULE package.  See the :doc:`Build package <Build_package>`
-page for more info.
+doc page for more info.
 
 Related commands
 """"""""""""""""

doc/src/improper_ring.rst

@@ -72,8 +72,8 @@ Restrictions
 """"""""""""
 
 This improper style can only be used if LAMMPS was built with the
-MOLECULE package.  See the :doc:`Build package <Build_package>` doc
+EXTRA-MOLECULE package.  See the :doc:`Build package <Build_package>`
-page for more info.
+doc page for more info.
 
 Related commands
 """"""""""""""""

doc/src/labelmap.rst

@@ -24,6 +24,7 @@ Examples
 
 .. code-block:: LAMMPS
 
+   labelmap atom 1 c1 2 hc 3 cp 4 nt
    labelmap atom 3 carbon 4 'c3"' 5 "c1'" 6 "c#"
    labelmap atom $(label2type(atom,carbon)) C  # change type label from 'carbon' to 'C'
    labelmap clear

doc/src/pair_coul.rst

@@ -379,10 +379,11 @@ These pair styles can only be used via the *pair* keyword of the
 Restrictions
 """"""""""""
 
-The *coul/cut/global*, *coul/long*, *coul/msm*, *coul/streitz*, and *tip4p/long* styles
+The *coul/long*, *coul/msm*, *coul/streitz*, and *tip4p/long* styles are
-are part of the KSPACE package.  They are only enabled if LAMMPS was built
+part of the KSPACE package.  The *coul/cut/global* and *coul/exclude* are
-with that package.  See the :doc:`Build package <Build_package>` doc page
+part of the EXTRA-PAIR package.  A pair style is only enabled if LAMMPS was
-for more info.
+built with its corresponding package.  See the :doc:`Build package <Build_package>`
+doc page for more info.
 
 Related commands
 """"""""""""""""

doc/src/pair_dpd_coul_slater_long.rst

@@ -0,0 +1,179 @@
+.. index:: pair_style dpd/coul/slater/long
+.. index:: pair_style dpd/coul/slater/long/gpu
+
+pair_style dpd/coul/slater/long command
+=======================================
+
+Accelerator Variants: *dpd/coul/slater/long/gpu*
+
+Syntax
+""""""
+
+.. code-block:: LAMMPS
+
+   pair_style dpd/coul/slater/long T cutoff_DPD seed lambda cutoff_coul
+
+   pair_coeff I J a_IJ Gamma is_charged
+
+* T = temperature (temperature units) (dpd only)
+* cutoff_DPD = global cutoff for DPD interactions (distance units)
+* seed = random # seed (positive integer)
+* lambda = decay length of the charge (distance units)
+* cutoff_coul = real part cutoff for Coulombic interactions (distance units)
+* I,J = numeric atom types, or type labels
+* Gamma = DPD Gamma coefficient
+* is_charged (boolean) set to yes if I and J are charged beads
+
+Examples
+""""""""
+
+.. code-block:: LAMMPS
+
+   pair_style dpd/coul/slater/long 1.0 2.5 34387 0.25 3.0
+   pair_coeff 1 1 78.0 4.5 # not charged by default
+   pair_coeff 2 2 78.0 4.5 yes
+
+
+Description
+"""""""""""
+
+.. versionadded:: TBD
+
+Style *dpd/coul/slater/long* computes a force field for dissipative particle dynamics
+(DPD) following the exposition in :ref:`(Groot) <Groot5>` with the addition of
+electrostatic interactions. The coulombic forces in mesoscopic models
+employ potentials without explicit excluded-volume interactions.
+The goal is to prevent artificial ionic pair formation by including a charge
+distribution in the Coulomb potential, following the formulation of
+:ref:`(Melchor) <Melchor1>`:
+
+The force on bead I due to bead J is given as a sum
+of 4 terms
+
+.. math::
+
+   \vec{f}  = & (F^C + F^D + F^R + F^E) \hat{r_{ij}} \\
+   F^C      = & A w(r) \qquad \qquad \qquad \qquad \qquad r < r_c \\
+   F^D      = & - \gamma w^2(r) (\hat{r_{ij}} \bullet \vec{v}_{ij}) \qquad \qquad r < r_c \\
+   F^R      = & \sigma w(r) \alpha (\Delta t)^{-1/2} \qquad \qquad \qquad r < r_c \\
+   w(r)     = & 1 - \frac{r}{r_c} \\
+   F^E      = & \frac{Cq_iq_j}{\epsilon r^2} \left( 1- exp\left( \frac{2r_{ij}}{\lambda} \right) \left( 1 + \frac{2r_{ij}}{\lambda} \left( 1 + \frac{r_{ij}}{\lambda} \right)\right) \right)
+
+where :math:`F^C` is a conservative force, :math:`F^D` is a dissipative
+force, :math:`F^R` is a random force, and :math:`F^E` is an electrostatic force.
+:math:`\hat{r_{ij}}` is a unit vector in the direction
+:math:`r_i - r_j`, :math:`\vec{v}_{ij}` is
+the vector difference in velocities of the two atoms :math:`\vec{v}_i -
+\vec{v}_j`, :math:`\alpha` is a Gaussian random number with zero mean
+and unit variance, *dt* is the timestep size, and :math:`w(r)` is a
+weighting factor that varies between 0 and 1.  :math:`r_c` is the
+pairwise cutoff.  :math:`\sigma` is set equal to :math:`\sqrt{2 k_B T
+\gamma}`, where :math:`k_B` is the Boltzmann constant and *T* is the
+temperature parameter in the pair_style command.
+C is the same Coulomb conversion factor as in the pair_styles
+coul/cut and coul/long. In this way the Coulomb
+interaction between ions is corrected at small distances r, and
+the long-range interactions are computed either by the Ewald or the PPPM technique.
+
+
+The following parameters must be defined for each
+pair of atoms types via the :doc:`pair_coeff <pair_coeff>` command as in
+the examples above, or in the data file or restart files read by the
+:doc:`read_data <read_data>` or :doc:`read_restart <read_restart>`
+commands:
+
+* A (force units)
+* :math:`\gamma` (force/velocity units)
+* is_charged (boolean)
+
+
+.. note::
+
+   This style is the combination of :doc:`pair_style dpd <pair_dpd>` and :doc:`pair_style coul/slater/long <pair_coul_slater>`.
+
+----------
+
+.. include:: accel_styles.rst
+
+----------
+
+Mixing, shift, table, tail correction, restart, rRESPA info
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+This pair style does not support mixing.  Thus, coefficients for all
+I,J pairs must be specified explicitly.
+
+This pair style does not support the :doc:`pair_modify <pair_modify>`
+shift option for the energy of the pair interaction.
+
+The :doc:`pair_modify <pair_modify>` table option is not relevant
+for this pair style.
+
+This pair style does not support the :doc:`pair_modify <pair_modify>`
+tail option for adding long-range tail corrections to energy and
+pressure.
+
+This pair style writes its information to :doc:`binary restart files
+<restart>`, so pair_style and pair_coeff commands do not need to be
+specified in an input script that reads a restart file.  Note that the
+user-specified random number seed is stored in the restart file, so when
+a simulation is restarted, each processor will re-initialize its random
+number generator the same way it did initially.  This means the random
+forces will be random, but will not be the same as they would have been
+if the original simulation had continued past the restart time.
+
+This pair style can only be used via the *pair* keyword of the
+:doc:`run_style respa <run_style>` command.  They do not support the
+*inner*, *middle*, *outer* keywords.
+
+
+----------
+
+Restrictions
+""""""""""""
+
+This style is part of the DPD-BASIC package.  It is only enabled if
+LAMMPS was built with that package.  See the :doc:`Build package
+<Build_package>` page for more info.
+
+The default frequency for rebuilding neighbor lists is every 10 steps
+(see the :doc:`neigh_modify <neigh_modify>` command). This may be too
+infrequent since particles move rapidly and
+can overlap by large amounts.  If this setting yields a non-zero number
+of "dangerous" reneighborings (printed at the end of a simulation), you
+should experiment with forcing reneighboring more often and see if
+system energies/trajectories change.
+
+This pair style requires you to use the :doc:`comm_modify vel yes
+<comm_modify>` command so that velocities are stored by ghost atoms.
+
+This pair style also requires the long-range solvers included in the KSPACE package.
+
+
+This pair style will not restart exactly when using the
+:doc:`read_restart <read_restart>` command, though they should provide
+statistically similar results.  This is because the forces they compute
+depend on atom velocities.  See the :doc:`read_restart <read_restart>`
+command for more details.
+
+Related commands
+""""""""""""""""
+
+:doc:`pair_style dpd <pair_dpd>`, :doc:`pair_style coul/slater/long <pair_coul_slater>`,
+:doc:`pair_coeff <pair_coeff>`, :doc:`fix nvt <fix_nh>`, :doc:`fix langevin <fix_langevin>`,
+:doc:`pair_style srp <pair_srp>`, :doc:`fix mvv/dpd <fix_mvv_dpd>`.
+
+Default
+"""""""
+
+is_charged = no
+
+----------
+
+.. _Groot5:
+
+**(Groot)** Groot and Warren, J Chem Phys, 107, 4423-35 (1997).
+
+.. _Melchor1:
+
+**(Melchor)** Gonzalez-Melchor, Mayoral, Velazquez, and Alejandre, J Chem Phys, 125, 224107 (2006).

doc/src/pair_granular.rst

@@ -111,7 +111,7 @@ For the *hertz* model, the normal component of force is given by:
 
    \mathbf{F}_{ne, Hertz} = k_n R_{eff}^{1/2}\delta_{ij}^{3/2} \mathbf{n}
 
-Here, :math:`R_{eff} = \frac{R_i R_j}{R_i + R_j}` is the effective
+Here, :math:`R_{eff} = R = \frac{R_i R_j}{R_i + R_j}` is the effective
 radius, denoted for simplicity as *R* from here on.  For *hertz*, the
 units of the spring constant :math:`k_n` are *force*\ /\ *length*\ \^2, or
 equivalently *pressure*\ .
@@ -120,13 +120,14 @@ For the *hertz/material* model, the force is given by:
 
 .. math::
 
-   \mathbf{F}_{ne, Hertz/material} = \frac{4}{3} E_{eff} R_{eff}^{1/2}\delta_{ij}^{3/2} \mathbf{n}
+   \mathbf{F}_{ne, Hertz/material} = \frac{4}{3} E_{eff} R^{1/2}\delta_{ij}^{3/2} \mathbf{n}
 
-Here, :math:`E_{eff} = E = \left(\frac{1-\nu_i^2}{E_i} + \frac{1-\nu_j^2}{E_j}\right)^{-1}` is the effective Young's
+Here, :math:`E_{eff} = E = \left(\frac{1-\nu_i^2}{E_i} + \frac{1-\nu_j^2}{E_j}\right)^{-1}`
-modulus, with :math:`\nu_i, \nu_j` the Poisson ratios of the particles of
+is the effective Young's modulus, with :math:`\nu_i, \nu_j` the Poisson ratios
-types *i* and *j*\ . Note that if the elastic modulus and the shear
+of the particles of types *i* and *j*. :math:`E_{eff}` is denoted as *E* from here on.
-modulus of the two particles are the same, the *hertz/material* model
+Note that if the elastic modulus and the shear modulus of the two particles are the
-is equivalent to the *hertz* model with :math:`k_n = 4/3 E_{eff}`
+same, the *hertz/material* model is equivalent to the *hertz* model with
+:math:`k_n = 4/3 E`
 
 The *dmt* model corresponds to the
 :ref:`(Derjaguin-Muller-Toporov) <DMT1975>` cohesive model, where the force
@@ -270,7 +271,8 @@ where :math:`k_n = \frac{4}{3} E_{eff}` for the *hertz/material* model. Since
 *coeff_restitution* accounts for the effective mass, effective radius, and
 pairwise overlaps (except when used with the *hooke* normal model) when calculating
 the damping coefficient, it accurately reproduces the specified coefficient of
-restitution for both monodisperse and polydisperse particle pairs.
+restitution for both monodisperse and polydisperse particle pairs.  This damping
+model is not compatible with cohesive normal models such as *JKR* or *DMT*.
 
 The total normal force is computed as the sum of the elastic and
 damping components:
@@ -441,11 +443,11 @@ discussion above. To match the Mindlin solution, one should set
 
    G_{eff} = \left(\frac{2-\nu_i}{G_i} + \frac{2-\nu_j}{G_j}\right)^{-1}
 
-where :math:`G` is the shear modulus, related to Young's modulus :math:`E`
+where :math:`G_i` is the shear modulus of a particle of type :math:`i`, related to Young's
-and Poisson's ratio :math:`\nu` by :math:`G = E/(2(1+\nu))`. This can also be
+modulus :math:`E_i` and Poisson's ratio :math:`\nu_i` by :math:`G_i = E_i/(2(1+\nu_i))`.
-achieved by specifying *NULL* for :math:`k_t`, in which case a
+This can also be achieved by specifying *NULL* for :math:`k_t`, in which case a
-normal contact model that specifies material parameters :math:`E` and
+normal contact model that specifies material parameters :math:`E_i` and
-:math:`\nu` is required (e.g. *hertz/material*, *dmt* or *jkr*\ ). In this
+:math:`\nu_i` is required (e.g. *hertz/material*, *dmt* or *jkr*\ ). In this
 case, mixing of the shear modulus for different particle types *i* and
 *j* is done according to the formula above.
 
@@ -575,7 +577,7 @@ opposite torque on each particle, according to:
 
 .. math::
 
-   \tau_{roll,i} =  R_{eff} \mathbf{n} \times \mathbf{F}_{roll}
+   \tau_{roll,i} =  R \mathbf{n} \times \mathbf{F}_{roll}
 
 .. math::
 

doc/src/pair_hybrid.rst

@@ -1,28 +1,41 @@
 .. index:: pair_style hybrid
 .. index:: pair_style hybrid/kk
+.. index:: pair_style hybrid/omp
+.. index:: pair_style hybrid/molecular
+.. index:: pair_style hybrid/molecular/omp
 .. index:: pair_style hybrid/overlay
+.. index:: pair_style hybrid/overlay/omp
 .. index:: pair_style hybrid/overlay/kk
 .. index:: pair_style hybrid/scaled
+.. index:: pair_style hybrid/scaled/omp
 
 pair_style hybrid command
 =========================
 
-Accelerator Variants: *hybrid/kk*
+Accelerator Variants: *hybrid/kk*, *hybrid/omp*
+
+pair_style hybrid/molecular command
+===================================
+
+Accelerator Variant: *hybrid/molecular/omp*
 
 pair_style hybrid/overlay command
 =================================
 
-Accelerator Variants: *hybrid/overlay/kk*
+Accelerator Variants: *hybrid/overlay/kk*, *hybrid/overlay/omp*
 
 pair_style hybrid/scaled command
 ==================================
 
+Accelerator Variant: *hybrid/scaled/omp*
+
 Syntax
 """"""
 
 .. code-block:: LAMMPS
 
    pair_style hybrid style1 args style2 args ...
+   pair_style hybrid/molecular factor1 style1 args factor2 style 2 args
    pair_style hybrid/overlay style1 args style2 args ...
    pair_style hybrid/scaled factor1 style1 args factor2 style 2 args ...
 
@@ -47,6 +60,10 @@ Examples
    pair_coeff * * tersoff Si.tersoff Si
    pair_coeff * * sw Si.sw Si
 
+   pair_style hybrid/molecular lj/cut 2.5 lj/cut 2.5
+   pair_coeff * * lj/cut 1 1.0 1.0
+   pair_coeff * * lj/cut 2 1.5 1.0
+
    variable one equal ramp(1.0,0.0)
    variable two equal 1.0-v_one
    pair_style hybrid/scaled v_one lj/cut 2.5 v_two morse 2.5
@@ -56,17 +73,26 @@ Examples
 Description
 """""""""""
 
-The *hybrid*, *hybrid/overlay*, and *hybrid/scaled* styles enable the
+The *hybrid*, *hybrid/overlay*, *hybrid/molecular*, and *hybrid/scaled*
-use of multiple pair styles in one simulation.  With the *hybrid* style,
+styles enable the use of multiple pair styles in one simulation.  With
-exactly one pair style is assigned to each pair of atom types.  With the
+the *hybrid* style, exactly one pair style is assigned to each pair of
-*hybrid/overlay* and *hybrid/scaled* styles, one or more pair styles can
+atom types.  With the *hybrid/overlay* and *hybrid/scaled* styles, one
-be assigned to each pair of atom types.  The assignment of pair styles
+or more pair styles can be assigned to each pair of atom types.  With
-to type pairs is made via the :doc:`pair_coeff <pair_coeff>` command.
+the hybrid/molecular style, pair styles are assigned to either intra-
-The major difference between the *hybrid/overlay* and *hybrid/scaled*
+or inter-molecular interactions.
-styles is that the *hybrid/scaled* adds a scale factor for each
+
-sub-style contribution to forces, energies and stresses.  Because of the
+The assignment of pair styles to type pairs is made via the
-added complexity, the *hybrid/scaled* style has more overhead and thus
+:doc:`pair_coeff <pair_coeff>` command.  The major difference between
-may be slower than *hybrid/overlay*.
+the *hybrid/overlay* and *hybrid/scaled* styles is that the
+*hybrid/scaled* adds a scale factor for each sub-style contribution to
+forces, energies and stresses.  Because of the added complexity, the
+*hybrid/scaled* style has more overhead and thus may be slower than
+*hybrid/overlay*.
+
+The *hybrid/molecular* pair style accepts *only* two sub-styles: the
+first is assigned to intra-molecular interactions (i.e. both atoms
+have the same molecule ID), the second to inter-molecular interactions
+(i.e. interacting atoms have different molecule IDs).
 
 Here are two examples of hybrid simulations.  The *hybrid* style could
 be used for a simulation of a metal droplet on a LJ surface.  The metal
@@ -476,6 +502,8 @@ the same or else LAMMPS will generate an error.
 Pair style *hybrid/scaled* currently only works for non-accelerated
 pair styles and pair styles from the OPT package.
 
+Pair style *hybrid/molecular* is not compatible with manybody potentials.
+
 When using pair styles from the GPU package they must not be listed
 multiple times.  LAMMPS will detect this and abort.
 

doc/src/pair_oxdna.rst

@@ -37,18 +37,19 @@ Syntax
 
      *oxdna/stk* args = seq T xi kappa 6.0 0.4 0.9 0.32 0.75 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 2.0 0.65 2.0 0.65
        seq = seqav (for average sequence stacking strength) or seqdep (for sequence-dependent stacking strength)
-       T = temperature (oxDNA units, 0.1 = 300 K)
+       T = temperature (LJ units: 0.1 = 300 K, real units: 300 = 300 K)
-       xi = 1.3448 (temperature-independent coefficient in stacking strength)
+       xi = 1.3448 (LJ units) or 8.01727944817084 (real units), temperature-independent coefficient in stacking strength
-       kappa = 2.6568 (coefficient of linear temperature dependence in stacking strength)
+       kappa = 2.6568 (LJ units) or 0.005279604 (real units), coefficient of linear temperature dependence in stacking strength
      *oxdna/hbond* args = seq eps 8.0 0.4 0.75 0.34 0.7 1.5 0 0.7 1.5 0 0.7 1.5 0 0.7 0.46 3.141592653589793 0.7 4.0 1.5707963267948966 0.45 4.0 1.5707963267948966 0.45
        seq = seqav (for average sequence base-pairing strength) or seqdep (for sequence-dependent base-pairing strength)
-       eps = 1.077 (between base pairs A-T and C-G) or 0 (all other pairs)
+       eps = 1.077 (LJ units) or 6.42073911784652 (real units), average hydrogen bonding strength between A-T and C-G Watson-Crick base pairs, 0 between all other pairs
 
 Examples
 """"""""
 
 .. code-block:: LAMMPS
 
+   # LJ units
    pair_style hybrid/overlay oxdna/excv oxdna/stk oxdna/hbond oxdna/xstk oxdna/coaxstk
    pair_coeff * * oxdna/excv    2.0 0.7 0.675 2.0 0.515 0.5 2.0 0.33 0.32
    pair_coeff * * oxdna/stk     seqdep 0.1 1.3448 2.6568 6.0 0.4 0.9 0.32 0.75 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 2.0 0.65 2.0 0.65
@@ -58,55 +59,105 @@ Examples
    pair_coeff * * oxdna/xstk    47.5 0.575 0.675 0.495 0.655 2.25 0.791592653589793 0.58 1.7 1.0 0.68 1.7 1.0 0.68 1.5 0 0.65 1.7 0.875 0.68 1.7 0.875 0.68
    pair_coeff * * oxdna/coaxstk 46.0 0.4 0.6 0.22 0.58 2.0 2.541592653589793 0.65 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 2.0 -0.65 2.0 -0.65
 
+   pair_style hybrid/overlay oxdna/excv oxdna/stk oxdna/hbond oxdna/xstk oxdna/coaxstk
+   pair_coeff * * oxdna/excv    oxdna_lj.cgdna
+   pair_coeff * * oxdna/stk     seqav 0.1 1.3448 2.6568 oxdna_lj.cgdna
+   pair_coeff * * oxdna/hbond   seqav oxdna_lj.cgdna
+   pair_coeff 1 4 oxdna/hbond   seqav oxdna_lj.cgdna
+   pair_coeff 2 3 oxdna/hbond   seqav oxdna_lj.cgdna
+   pair_coeff * * oxdna/xstk    oxdna_lj.cgdna
+   pair_coeff * * oxdna/coaxstk oxdna_lj.cgdna
+
+   # Real units
+   pair_style hybrid/overlay oxdna/excv oxdna/stk oxdna/hbond oxdna/xstk oxdna/coaxstk
+   pair_coeff * * oxdna/excv    11.92337812042065 5.9626 5.74965 11.92337812042065 4.38677 4.259 11.92337812042065 2.81094 2.72576
+   pair_coeff * * oxdna/stk     seqdep 300.0 8.01727944817084 0.005279604 0.70439070204273 3.4072 7.6662 2.72576 6.3885 1.3 0.0 0.8 0.9 0.0 0.95 0.9 0.0 0.95 2.0 0.65 2.0 0.65
+   pair_coeff * * oxdna/hbond   seqdep 0.0 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592654 0.7 4.0 1.570796327 0.45 4.0 1.570796327 0.45
+   pair_coeff 1 4 oxdna/hbond   seqdep 6.42073911784652 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592654 0.7 4.0 1.570796327 0.45 4.0 1.570796327 0.45
+   pair_coeff 2 3 oxdna/hbond   seqdep 6.42073911784652 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592654 0.7 4.0 1.570796327 0.45 4.0 1.570796327 0.45
+   pair_coeff * * oxdna/xstk    3.9029021145006 4.89785 5.74965 4.21641 5.57929 2.25 0.791592654 0.58 1.7 1.0 0.68 1.7 1.0 0.68 1.5 0.0 0.65 1.7 0.875 0.68 1.7 0.875 0.68
+   pair_coeff * * oxdna/coaxstk 3.77965257404268 3.4072 5.1108 1.87396 4.94044 2.0 2.541592654 0.65 1.3 0.0 0.8 0.9 0.0 0.95 0.9 0.0 0.95 2.0 -0.65 2.0 -0.65
+
+   pair_style hybrid/overlay oxdna/excv oxdna/stk oxdna/hbond oxdna/xstk oxdna/coaxstk
+   pair_coeff * * oxdna/excv    oxdna_real.cgdna
+   pair_coeff * * oxdna/stk     seqav 300.0 8.01727944817084 0.005279604 oxdna_real.cgdna
+   pair_coeff * * oxdna/hbond   seqav oxdna_real.cgdna
+   pair_coeff 1 4 oxdna/hbond   seqav oxdna_real.cgdna
+   pair_coeff 2 3 oxdna/hbond   seqav oxdna_real.cgdna
+   pair_coeff * * oxdna/xstk    oxdna_real.cgdna
+   pair_coeff * * oxdna/coaxstk oxdna_real.cgdna
+
+.. note::
+
+   The coefficients in the above examples are provided in forms
+   compatible with both *units lj* and *units real* (see documentation
+   of :doc:`units <units>`).  These can also be read from a potential
+   file with correct unit style by specifying the name of the
+   file. Several potential files for each unit style are included in the
+   ``potentials`` directory of the LAMMPS distribution.
+
 Description
 """""""""""
 
-The *oxdna* pair styles compute the pairwise-additive parts of the oxDNA force field
+The *oxdna* pair styles compute the pairwise-additive parts of the oxDNA
-for coarse-grained modelling of DNA. The effective interaction between the nucleotides consists of potentials for the
+force field for coarse-grained modelling of DNA. The effective
-excluded volume interaction *oxdna/excv*, the stacking *oxdna/stk*, cross-stacking *oxdna/xstk*
+interaction between the nucleotides consists of potentials for the
-and coaxial stacking interaction *oxdna/coaxstk* as well
+excluded volume interaction *oxdna/excv*, the stacking *oxdna/stk*,
-as the hydrogen-bonding interaction *oxdna/hbond* between complementary pairs of nucleotides on
+cross-stacking *oxdna/xstk* and coaxial stacking interaction
-opposite strands. Average sequence or sequence-dependent stacking and base-pairing strengths
+*oxdna/coaxstk* as well as the hydrogen-bonding interaction
-are supported :ref:`(Sulc) <Sulc1>`. Quasi-unique base-pairing between nucleotides can be achieved by using
+*oxdna/hbond* between complementary pairs of nucleotides on opposite
-more complementary pairs of atom types like 5-8 and 6-7, 9-12 and 10-11, 13-16 and 14-15, etc.
+strands. Average sequence or sequence-dependent stacking and
-This prevents the hybridization of in principle complementary bases within Ntypes/4 bases
+base-pairing strengths are supported :ref:`(Sulc) <Sulc1>`. Quasi-unique
-up and down along the backbone.
+base-pairing between nucleotides can be achieved by using more
+complementary pairs of atom types like 5-8 and 6-7, 9-12 and 10-11,
+13-16 and 14-15, etc.  This prevents the hybridization of in principle
+complementary bases within Ntypes/4 bases up and down along the
+backbone.
 
-The exact functional form of the pair styles is rather complex.
+The exact functional form of the pair styles is rather complex.  The
-The individual potentials consist of products of modulation factors,
+individual potentials consist of products of modulation factors, which
-which themselves are constructed from a number of more basic potentials
+themselves are constructed from a number of more basic potentials
-(Morse, Lennard-Jones, harmonic angle and distance) as well as quadratic smoothing and modulation terms.
+(Morse, Lennard-Jones, harmonic angle and distance) as well as quadratic
-We refer to :ref:`(Ouldridge-DPhil) <Ouldridge-DPhil1>` and :ref:`(Ouldridge) <Ouldridge1>`
+smoothing and modulation terms.  We refer to :ref:`(Ouldridge-DPhil)
-for a detailed description of the oxDNA force field.
+<Ouldridge-DPhil1>` and :ref:`(Ouldridge) <Ouldridge1>` for a detailed
+description of the oxDNA force field.
 
 .. note::
 
-   These pair styles have to be used together with the related oxDNA bond style
+   These pair styles have to be used together with the related oxDNA
-   *oxdna/fene* for the connectivity of the phosphate backbone (see also documentation of
+   bond style *oxdna/fene* for the connectivity of the phosphate
-   :doc:`bond_style oxdna/fene <bond_oxdna>`). Most of the coefficients
+   backbone (see also documentation of :doc:`bond_style oxdna/fene
-   in the above example have to be kept fixed and cannot be changed without reparameterizing the entire model.
+   <bond_oxdna>`). Most of the coefficients in the above example have to
-   Exceptions are the first four coefficients after *oxdna/stk* (seq=seqdep, T=0.1, xi=1.3448 and kappa=2.6568 in the above example)
+   be kept fixed and cannot be changed without reparameterizing the
-   and the first coefficient after *oxdna/hbond* (seq=seqdep in the above example).
+   entire model.  Exceptions are the first four coefficients after
-   When using a Langevin thermostat, e.g. through :doc:`fix langevin <fix_langevin>`
+   *oxdna/stk* (seq=seqdep, T=0.1, xi=1.3448 and kappa=2.6568 and
-   or :doc:`fix nve/dotc/langevin <fix_nve_dotc_langevin>`
+   corresponding *real unit* equivalents in the above examples) and the
-   the temperature coefficients have to be matched to the one used in the fix.
+   first coefficient after *oxdna/hbond* (seq=seqdep in the above
+   example).  When using a Langevin thermostat, e.g. through :doc:`fix
+   langevin <fix_langevin>` or :doc:`fix nve/dotc/langevin
+   <fix_nve_dotc_langevin>` the temperature coefficients have to be
+   matched to the one used in the fix.
 
 .. note::
 
-   These pair styles have to be used with the *atom_style hybrid bond ellipsoid oxdna*
+   These pair styles have to be used with the *atom_style hybrid bond
-   (see documentation of :doc:`atom_style <atom_style>`). The *atom_style oxdna*
+   ellipsoid oxdna* (see documentation of :doc:`atom_style
-   stores the 3'-to-5' polarity of the nucleotide strand, which is set through
+   <atom_style>`). The *atom_style oxdna* stores the 3'-to-5' polarity
-   the bond topology in the data file. The first (second) atom in a bond definition
+   of the nucleotide strand, which is set through the bond topology in
-   is understood to point towards the 3'-end (5'-end) of the strand.
+   the data file. The first (second) atom in a bond definition is
+   understood to point towards the 3'-end (5'-end) of the strand.
 
-Example input and data files for DNA duplexes can be found in examples/PACKAGES/cgdna/examples/oxDNA/ and /oxDNA2/.
+Example input and data files for DNA duplexes can be found in
-A simple python setup tool which creates single straight or helical DNA strands,
+``examples/PACKAGES/cgdna/examples/oxDNA/`` and ``.../oxDNA2/``.  A
-DNA duplexes or arrays of DNA duplexes can be found in examples/PACKAGES/cgdna/util/.
+simple python setup tool which creates single straight or helical DNA
+strands, DNA duplexes or arrays of DNA duplexes can be found in
+``examples/PACKAGES/cgdna/util/``.
 
 Please cite :ref:`(Henrich) <Henrich1>` in any publication that uses
-this implementation. An updated documentation that contains general information
+this implementation. An updated documentation that contains general
-on the model, its implementation and performance as well as the structure of
+information on the model, its implementation and performance as well as
-the data and input file can be found `here <PDF/CG-DNA.pdf>`_.
+the structure of the data and input file can be found `here
+<PDF/CG-DNA.pdf>`_.
 
 Please cite also the relevant oxDNA publications
 :ref:`(Ouldridge) <Ouldridge1>`,
@@ -115,6 +166,57 @@ and :ref:`(Sulc) <Sulc1>`.
 
 ----------
 
+Potential file reading
+""""""""""""""""""""""
+
+For each pair style above the first non-modifiable argument can be a
+filename, and if it is, no further arguments should be
+supplied. Therefore the following command:
+
+.. code-block:: LAMMPS
+
+   pair_coeff 1 4 oxdna/hbond   seqav oxdna_lj.cgdna
+
+will be interpreted as a request to read the corresponding hydrogen
+bonding potential parameters from the file with the given name. The file
+can define multiple potential parameters for both bonded and pair
+interactions, but for the example pair interaction above there must
+exist in the file a line of the form:
+
+.. code-block:: LAMMPS
+
+   1 4 hbond     <coefficients>
+
+If potential customization is required, the potential file reading can
+be mixed with the manual specification of the potential parameters. For
+example, the following command:
+
+.. code-block:: LAMMPS
+
+   pair_style hybrid/overlay oxdna/excv oxdna/stk oxdna/hbond oxdna/xstk oxdna/coaxstk
+   pair_coeff * * oxdna/excv    oxdna_lj.cgdna
+   pair_coeff * * oxdna/stk     seqav 0.1 1.3448 2.6568 6.0 0.4 0.9 0.32 0.75 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 2.0 0.65 2.0 0.65
+   pair_coeff * * oxdna/hbond   seqav oxdna_lj.cgdna
+   pair_coeff 1 4 oxdna/hbond   seqav oxdna_lj.cgdna
+   pair_coeff 2 3 oxdna/hbond   seqav oxdna_lj.cgdna
+   pair_coeff * * oxdna/xstk    oxdna_lj.cgdna
+   pair_coeff * * oxdna/coaxstk 46.0 0.4 0.6 0.22 0.58 2.0 2.541592653589793 0.65 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 2.0 -0.65 2.0 -0.65
+
+will read the stacking and coaxial stacking potential parameters from
+the manual specification and all others from the potential file
+*oxdna_lj.cgdna*.
+
+There are sample potential files for each unit style in the
+``potentials`` directory of the LAMMPS distribution. The potential file
+unit system must align with the units defined via the :doc:`units
+<units>` command. For conversion between different *LJ* and *real* unit
+systems for oxDNA, the python tool *lj2real.py* located in the
+``examples/PACKAGES/cgdna/util/`` directory can be used. This tool
+assumes similar file structure to the examples found in
+``examples/PACKAGES/cgdna/examples/``.
+
+----------
+
 Restrictions
 """"""""""""
 

doc/src/pair_oxdna2.rst

@@ -41,14 +41,14 @@ Syntax
 
      *oxdna2/stk* args = seq T xi kappa 6.0 0.4 0.9 0.32 0.75 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 2.0 0.65 2.0 0.65
        seq = seqav (for average sequence stacking strength) or seqdep (for sequence-dependent stacking strength)
-       T = temperature (oxDNA units, 0.1 = 300 K)
+       T = temperature (LJ units: 0.1 = 300 K, real units: 300 = 300 K)
-       xi = 1.3523 (temperature-independent coefficient in stacking strength)
+       xi = 1.3523 (LJ units) or 8.06199211612242 (real units), temperature-independent coefficient in stacking strength
-       kappa = 2.6717 (coefficient of linear temperature dependence in stacking strength)
+       kappa = 2.6717 (LJ units) or 0.005309213 (real units), coefficient of linear temperature dependence in stacking strength
      *oxdna2/hbond* args = seq eps 8.0 0.4 0.75 0.34 0.7 1.5 0 0.7 1.5 0 0.7 1.5 0 0.7 0.46 3.141592653589793 0.7 4.0 1.5707963267948966 0.45 4.0 1.5707963267948966 0.45
        seq = seqav (for average sequence base-pairing strength) or seqdep (for sequence-dependent base-pairing strength)
-       eps = 1.0678 (between base pairs A-T and C-G) or 0 (all other pairs)
+       eps = 1.0678 (LJ units) or 6.36589157849259 (real units), average hydrogen bonding strength between A-T and C-G Watson-Crick base pairs, 0 between all other pairs
      *oxdna2/dh* args = T rhos qeff
-       T = temperature (oxDNA units, 0.1 = 300 K)
+       T = temperature (LJ units: 0.1 = 300 K, real units: 300 = 300 K)
        rhos = salt concentration (mole per litre)
        qeff = 0.815 (effective charge in elementary charges)
 
@@ -57,6 +57,7 @@ Examples
 
 .. code-block:: LAMMPS
 
+   # LJ units
    pair_style hybrid/overlay oxdna2/excv oxdna2/stk oxdna2/hbond oxdna2/xstk oxdna2/coaxstk oxdna2/dh
    pair_coeff * * oxdna2/excv    2.0 0.7 0.675 2.0 0.515 0.5 2.0 0.33 0.32
    pair_coeff * * oxdna2/stk     seqdep 0.1 1.3523 2.6717 6.0 0.4 0.9 0.32 0.75 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 2.0 0.65 2.0 0.65
@@ -67,61 +68,169 @@ Examples
    pair_coeff * * oxdna2/coaxstk 58.5 0.4 0.6 0.22 0.58 2.0 2.891592653589793 0.65 1.3 0 0.8 0.9 0 0.95 0.9 0 0.95 40.0 3.116592653589793
    pair_coeff * * oxdna2/dh      0.1 0.5 0.815
 
+   pair_style hybrid/overlay oxdna2/excv oxdna2/stk oxdna2/hbond oxdna2/xstk oxdna2/coaxstk oxdna2/dh
+   pair_coeff * * oxdna2/excv    oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/stk     seqdep 0.1 1.3523 2.6717 oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/hbond   seqdep oxdna2_lj.cgdna
+   pair_coeff 1 4 oxdna2/hbond   seqdep oxdna2_lj.cgdna
+   pair_coeff 2 3 oxdna2/hbond   seqdep oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/xstk    oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/coaxstk oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/dh      0.1 0.5 oxdna2_lj.cgdna
+
+   # Real units
+   pair_style hybrid/overlay oxdna2/excv oxdna2/stk oxdna2/hbond oxdna2/xstk oxdna2/coaxstk oxdna2/dh
+   pair_coeff * * oxdna2/excv    11.92337812042065 5.9626 5.74965 11.92337812042065 4.38677 4.259 11.92337812042065 2.81094 2.72576
+   pair_coeff * * oxdna2/stk     seqdep 300.0 8.06199211612242 0.005309213 0.70439070204273 3.4072 7.6662 2.72576 6.3885 1.3 0.0 0.8 0.9 0.0 0.95 0.9 0.0 0.95 2.0 0.65 2.0 0.65
+   pair_coeff * * oxdna2/hbond   seqdep 0.0 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592654 0.7 4.0 1.570796327 0.45 4.0 1.570796327 0.45
+   pair_coeff 1 4 oxdna2/hbond   seqdep 6.36589157849259 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592654 0.7 4.0 1.570796327 0.45 4.0 1.570796327 0.45
+   pair_coeff 2 3 oxdna2/hbond   seqdep 6.36589157849259 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592654 0.7 4.0 1.570796327 0.45 4.0 1.570796327 0.45
+   pair_coeff * * oxdna2/xstk    3.9029021145006 4.89785 5.74965 4.21641 5.57929 2.25 0.791592654 0.58 1.7 1.0 0.68 1.7 1.0 0.68 1.5 0.0 0.65 1.7 0.875 0.68 1.7 0.875 0.68
+   pair_coeff * * oxdna2/coaxstk 4.80673207785863 3.4072 5.1108 1.87396 4.94044 2.0 2.891592653589793 0.65 1.3 0.0 0.8 0.9 0.0 0.95 0.9 0.0 0.95 40.0 3.116592653589793
+   pair_coeff * * oxdna2/dh      300.0 0.5 0.815
+
+   pair_style hybrid/overlay oxdna2/excv oxdna2/stk oxdna2/hbond oxdna2/xstk oxdna2/coaxstk oxdna2/dh
+   pair_coeff * * oxdna2/excv    oxdna2_real.cgdna
+   pair_coeff * * oxdna2/stk     seqdep 300.0 8.06199211612242 0.005309213 oxdna2_real.cgdna
+   pair_coeff * * oxdna2/hbond   seqdep oxdna2_real.cgdna
+   pair_coeff 1 4 oxdna2/hbond   seqdep oxdna2_real.cgdna
+   pair_coeff 2 3 oxdna2/hbond   seqdep oxdna2_real.cgdna
+   pair_coeff * * oxdna2/xstk    oxdna2_real.cgdna
+   pair_coeff * * oxdna2/coaxstk oxdna2_real.cgdna
+   pair_coeff * * oxdna2/dh      300.0 0.5 oxdna2_real.cgdna
+
+.. note::
+
+   The coefficients in the above examples are provided in forms
+   compatible with both *units lj* and *units real* (see documentation
+   of :doc:`units <units>`).  These can also be read from a potential
+   file with correct unit style by specifying the name of the
+   file. Several potential files for each unit style are included in the
+   ``potentials`` directory of the LAMMPS distribution.
+
 Description
 """""""""""
 
-The *oxdna2* pair styles compute the pairwise-additive parts of the oxDNA force field
+The *oxdna2* pair styles compute the pairwise-additive parts of the
-for coarse-grained modelling of DNA. The effective interaction between the nucleotides consists of potentials for the
+oxDNA force field for coarse-grained modelling of DNA. The effective
-excluded volume interaction *oxdna2/excv*, the stacking *oxdna2/stk*, cross-stacking *oxdna2/xstk*
+interaction between the nucleotides consists of potentials for the
-and coaxial stacking interaction *oxdna2/coaxstk*, electrostatic Debye-Hueckel interaction *oxdna2/dh*
+excluded volume interaction *oxdna2/excv*, the stacking *oxdna2/stk*,
-as well as the hydrogen-bonding interaction *oxdna2/hbond* between complementary pairs of nucleotides on
+cross-stacking *oxdna2/xstk* and coaxial stacking interaction
-opposite strands. Average sequence or sequence-dependent stacking and base-pairing strengths
+*oxdna2/coaxstk*, electrostatic Debye-Hueckel interaction *oxdna2/dh* as
-are supported :ref:`(Sulc) <Sulc2>`. Quasi-unique base-pairing between nucleotides can be achieved by using
+well as the hydrogen-bonding interaction *oxdna2/hbond* between
-more complementary pairs of atom types like 5-8 and 6-7, 9-12 and 10-11, 13-16 and 14-15, etc.
+complementary pairs of nucleotides on opposite strands. Average sequence
-This prevents the hybridization of in principle complementary bases within Ntypes/4 bases
+or sequence-dependent stacking and base-pairing strengths are supported
+:ref:`(Sulc) <Sulc2>`. Quasi-unique base-pairing between nucleotides can
+be achieved by using more complementary pairs of atom types like 5-8 and
+6-7, 9-12 and 10-11, 13-16 and 14-15, etc.  This prevents the
+hybridization of in principle complementary bases within Ntypes/4 bases
 up and down along the backbone.
 
-The exact functional form of the pair styles is rather complex.
+The exact functional form of the pair styles is rather complex.  The
-The individual potentials consist of products of modulation factors,
+individual potentials consist of products of modulation factors, which
-which themselves are constructed from a number of more basic potentials
+themselves are constructed from a number of more basic potentials
-(Morse, Lennard-Jones, harmonic angle and distance) as well as quadratic smoothing and modulation terms.
+(Morse, Lennard-Jones, harmonic angle and distance) as well as quadratic
-We refer to :ref:`(Snodin) <Snodin2>` and the original oxDNA publications :ref:`(Ouldridge-DPhil) <Ouldridge-DPhil2>`
+smoothing and modulation terms.  We refer to :ref:`(Snodin) <Snodin2>`
-and  :ref:`(Ouldridge) <Ouldridge2>` for a detailed description of the oxDNA2 force field.
+and the original oxDNA publications :ref:`(Ouldridge-DPhil)
+<Ouldridge-DPhil2>` and :ref:`(Ouldridge) <Ouldridge2>` for a detailed
+description of the oxDNA2 force field.
 
 .. note::
 
-   These pair styles have to be used together with the related oxDNA2 bond style
+   These pair styles have to be used together with the related oxDNA2
-   *oxdna2/fene* for the connectivity of the phosphate backbone (see also documentation of
+   bond style *oxdna2/fene* for the connectivity of the phosphate
-   :doc:`bond_style oxdna2/fene <bond_oxdna>`). Most of the coefficients
+   backbone (see also documentation of :doc:`bond_style oxdna2/fene
-   in the above example have to be kept fixed and cannot be changed without reparameterizing the entire model.
+   <bond_oxdna>`). Most of the coefficients in the above example have to
-   Exceptions are the first four coefficients after *oxdna2/stk* (seq=seqdep, T=0.1, xi=1.3523 and kappa=2.6717 in the above example),
+   be kept fixed and cannot be changed without reparameterizing the
-   the first coefficient after *oxdna2/hbond* (seq=seqdep in the above example) and the three coefficients
+   entire model.  Exceptions are the first four coefficients after
-   after *oxdna2/dh* (T=0.1, rhos=0.5, qeff=0.815 in the above example). When using a Langevin thermostat
+   *oxdna2/stk* (seq=seqdep, T=0.1, xi=1.3523 and kappa=2.6717 and
-   e.g. through :doc:`fix langevin <fix_langevin>` or :doc:`fix nve/dotc/langevin <fix_nve_dotc_langevin>`
+   corresponding *real unit* equivalents in the above examples).  the
-   the temperature coefficients have to be matched to the one used in the fix.
+   first coefficient after *oxdna2/hbond* (seq=seqdep in the above
+   example) and the three coefficients after *oxdna2/dh* (T=0.1,
+   rhos=0.5, qeff=0.815 in the above example). When using a Langevin
+   thermostat e.g. through :doc:`fix langevin <fix_langevin>` or
+   :doc:`fix nve/dotc/langevin <fix_nve_dotc_langevin>` the temperature
+   coefficients have to be matched to the one used in the fix.
 
 .. note::
 
-   These pair styles have to be used with the *atom_style hybrid bond ellipsoid oxdna*
+   These pair styles have to be used with the *atom_style hybrid bond
-   (see documentation of :doc:`atom_style <atom_style>`). The *atom_style oxdna*
+   ellipsoid oxdna* (see documentation of :doc:`atom_style
-   stores the 3'-to-5' polarity of the nucleotide strand, which is set through
+   <atom_style>`). The *atom_style oxdna* stores the 3'-to-5' polarity
-   the bond topology in the data file. The first (second) atom in a bond definition
+   of the nucleotide strand, which is set through the bond topology in
-   is understood to point towards the 3'-end (5'-end) of the strand.
+   the data file. The first (second) atom in a bond definition is
+   understood to point towards the 3'-end (5'-end) of the strand.
 
-Example input and data files for DNA duplexes can be found in examples/PACKAGES/cgdna/examples/oxDNA/ and /oxDNA2/.
+Example input and data files for DNA duplexes can be found in
-A simple python setup tool which creates single straight or helical DNA strands,
+``examples/PACKAGES/cgdna/examples/oxDNA/`` and ``.../oxDNA2/``.  A
-DNA duplexes or arrays of DNA duplexes can be found in examples/PACKAGES/cgdna/util/.
+simple python setup tool which creates single straight or helical DNA
+strands, DNA duplexes or arrays of DNA duplexes can be found in
+``examples/PACKAGES/cgdna/util/``.
 
 Please cite :ref:`(Henrich) <Henrich2>` in any publication that uses
-this implementation. An updated documentation that contains general information
+this implementation. An updated documentation that contains general
-on the model, its implementation and performance as well as the structure of
+information on the model, its implementation and performance as well as
-the data and input file can be found `here <PDF/CG-DNA.pdf>`_.
+the structure of the data and input file can be found `here
+<PDF/CG-DNA.pdf>`_.
 
 Please cite also the relevant oxDNA2 publications
 :ref:`(Snodin) <Snodin2>` and :ref:`(Sulc) <Sulc2>`.
 
 ----------
 
+Potential file reading
+""""""""""""""""""""""
+
+For each pair style above the first non-modifiable argument can be a
+filename (with exception of Debye-Hueckel, for which the effective
+charge argument can be a filename), and if it is, no further arguments
+should be supplied.  Therefore the following command:
+
+.. code-block:: LAMMPS
+
+   pair_coeff 1 4 oxdna2/hbond   seqdep oxdna_real.cgdna
+
+will be interpreted as a request to read the corresponding hydrogen
+bonding potential parameters from the file with the given name.  The
+file can define multiple potential parameters for both bonded and pair
+interactions, but for the example pair interaction above there must
+exist in the file a line of the form:
+
+.. code-block:: LAMMPS
+
+  1 4 hbond     <coefficients>
+
+If potential customization is required, the potential file reading can
+be mixed with the manual specification of the potential parameters. For
+example, the following command:
+
+.. code-block:: LAMMPS
+
+   pair_style hybrid/overlay oxdna2/excv oxdna2/stk oxdna2/hbond oxdna2/xstk oxdna2/coaxstk oxdna2/dh
+   pair_coeff * * oxdna2/excv    2.0 0.7 0.675 2.0 0.515 0.5 2.0 0.33 0.32
+   pair_coeff * * oxdna2/stk     seqdep 0.1 1.3523 2.6717 oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/hbond   seqdep oxdna2_lj.cgdna
+   pair_coeff 1 4 oxdna2/hbond   seqdep oxdna2_lj.cgdna
+   pair_coeff 2 3 oxdna2/hbond   seqdep oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/xstk    oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/coaxstk oxdna2_lj.cgdna
+   pair_coeff * * oxdna2/dh      0.1 0.5 0.815
+
+will read the excluded volume and Debye-Hueckel effective charge *qeff*
+parameters from the manual specification and all others from the
+potential file *oxdna2_lj.cgdna*.
+
+There are sample potential files for each unit style in the ``potentials``
+directory of the LAMMPS distribution. The potential file unit system
+must align with the units defined via the :doc:`units <units>`
+command. For conversion between different *LJ* and *real* unit systems
+for oxDNA, the python tool *lj2real.py* located in the
+``examples/PACKAGES/cgdna/util/`` directory can be used. This tool assumes
+similar file structure to the examples found in
+``examples/PACKAGES/cgdna/examples/``.
+
+----------
+
 Restrictions
 """"""""""""
 

doc/src/pair_oxrna2.rst

@@ -41,14 +41,14 @@ Syntax
 
      *oxrna2/stk* args = seq T xi kappa 6.0 0.43 0.93 0.35 0.78 0.9 0 0.95 0.9 0 0.95 1.3 0 0.8 1.3 0 0.8 2.0 0.65 2.0 0.65
        seq = seqav (for average sequence stacking strength) or seqdep (for sequence-dependent stacking strength)
-       T = temperature (oxDNA units, 0.1 = 300 K)
+       T = temperature (LJ units: 0.1 = 300 K, real units: 300 = 300 K)
-       xi = 1.40206 (temperature-independent coefficient in stacking strength)
+       xi = 1.40206 (LJ units) or 8.35864576375849 (real units), temperature-independent coefficient in stacking strength
-       kappa = 2.77 (coefficient of linear temperature dependence in stacking strength)
+       kappa = 2.77 (LJ units) or 0.005504556 (real units), coefficient of linear temperature dependence in stacking strength
      *oxrna2/hbond* args = seq eps 8.0 0.4 0.75 0.34 0.7 1.5 0 0.7 1.5 0 0.7 1.5 0 0.7 0.46 3.141592653589793 0.7 4.0 1.5707963267948966 0.45 4.0 1.5707963267948966 0.45
        seq = seqav (for average sequence base-pairing strength) or seqdep (for sequence-dependent base-pairing strength)
-       eps = 0.870439 (between base pairs A-T, C-G and G-T) or 0 (all other pairs)
+       eps = 0.870439 (LJ units) or 5.18928666388042 (real units), average hydrogen bonding strength between A-U and C-G Watson-Crick and G-U wobble base pairs, 0 between all other pairs
      *oxrna2/dh* args = T rhos qeff
-       T = temperature (oxDNA units, 0.1 = 300 K)
+       T = temperature (LJ units: 0.1 = 300 K, real units: 300 = 300 K)
        rhos = salt concentration (mole per litre)
        qeff = 1.02455 (effective charge in elementary charges)
 
@@ -57,6 +57,7 @@ Examples
 
 .. code-block:: LAMMPS
 
+   # LJ units
    pair_style hybrid/overlay oxrna2/excv oxrna2/stk oxrna2/hbond oxrna2/xstk oxrna2/coaxstk oxrna2/dh
    pair_coeff * * oxrna2/excv    2.0 0.7 0.675 2.0 0.515 0.5 2.0 0.33 0.32
    pair_coeff * * oxrna2/stk     seqdep 0.1 1.40206 2.77 6.0 0.43 0.93 0.35 0.78 0.9 0 0.95 0.9 0 0.95 1.3 0 0.8 1.3 0 0.8 2.0 0.65 2.0 0.65
@@ -68,58 +69,168 @@ Examples
    pair_coeff * * oxrna2/coaxstk 80 0.5 0.6 0.42 0.58 2.0 2.592 0.65 1.3 0.151 0.8 0.9 0.685 0.95 0.9 0.685 0.95 2.0 -0.65 2.0 -0.65
    pair_coeff * * oxrna2/dh      0.1 0.5 1.02455
 
+   pair_style hybrid/overlay oxrna2/excv oxrna2/stk oxrna2/hbond oxrna2/xstk oxrna2/coaxstk oxrna2/dh
+   pair_coeff * * oxrna2/excv    oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/stk     seqdep 0.1 1.40206 2.77 oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff 1 4 oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff 2 3 oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff 3 4 oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/xstk    oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/coaxstk oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/dh      0.1 0.5 oxrna2_lj.cgdna
+
+   # Real units
+   pair_style hybrid/overlay oxrna2/excv oxrna2/stk oxrna2/hbond oxrna2/xstk oxrna2/coaxstk oxrna2/dh
+   pair_coeff * * oxrna2/excv    11.92337812042065 5.9626 5.74965 11.92337812042065 4.38677 4.259 11.92337812042065 2.81094 2.72576
+   pair_coeff * * oxrna2/stk     seqdep 300.0 8.35864576375849 0.005504556 0.70439070204273 3.66274 7.92174 2.9813 6.64404 0.9 0.0 0.95 0.9 0.0 0.95 1.3 0.0 0.8 1.3 0.0 0.8 2.0 0.65 2.0 0.65
+   pair_coeff * * oxrna2/hbond   seqdep 0.0 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592653589793 0.7 4.0 1.5707963267948966 0.45 4.0 1.5707963267948966 0.45
+   pair_coeff 1 4 oxrna2/hbond   seqdep 5.18928666388042 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592653589793 0.7 4.0 1.5707963267948966 0.45 4.0 1.5707963267948966 0.45
+   pair_coeff 2 3 oxrna2/hbond   seqdep 5.18928666388042 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592653589793 0.7 4.0 1.5707963267948966 0.45 4.0 1.5707963267948966 0.45
+   pair_coeff 3 4 oxrna2/hbond   seqdep 5.18928666388042 0.93918760272364 3.4072 6.3885 2.89612 5.9626 1.5 0.0 0.7 1.5 0.0 0.7 1.5 0.0 0.7 0.46 3.141592653589793 0.7 4.0 1.5707963267948966 0.45 4.0 1.5707963267948966 0.45
+   pair_coeff * * oxrna2/xstk    4.92690859644113 4.259 5.1108 3.57756 4.94044 2.25 0.505 0.58 1.7 1.266 0.68 1.7 1.266 0.68 1.7 0.309 0.68 1.7 0.309 0.68
+   pair_coeff * * oxrna2/coaxstk 6.57330882442206 4.259 5.1108 3.57756 4.94044 2.0 2.592 0.65 1.3 0.151 0.8 0.9 0.685 0.95 0.9 0.685 0.95 2.0 -0.65 2.0 -0.65
+   pair_coeff * * oxrna2/dh      300.0 0.5 1.02455
+
+   pair_style hybrid/overlay oxrna2/excv oxrna2/stk oxrna2/hbond oxrna2/xstk oxrna2/coaxstk oxrna2/dh
+   pair_coeff * * oxrna2/excv    oxrna2_real.cgdna
+   pair_coeff * * oxrna2/stk     seqdep 300.0 8.35864576375849 0.005504556 oxrna2_real.cgdna
+   pair_coeff * * oxrna2/hbond   seqdep oxrna2_real.cgdna
+   pair_coeff 1 4 oxrna2/hbond   seqdep oxrna2_real.cgdna
+   pair_coeff 2 3 oxrna2/hbond   seqdep oxrna2_real.cgdna
+   pair_coeff 3 4 oxrna2/hbond   seqdep oxrna2_real.cgdna
+   pair_coeff * * oxrna2/xstk    oxrna2_real.cgdna
+   pair_coeff * * oxrna2/coaxstk oxrna2_real.cgdna
+   pair_coeff * * oxrna2/dh      300.0 0.5 oxrna2_real.cgdna
+
+.. note::
+
+   The coefficients in the above examples are provided in forms
+   compatible with both *units lj* and *units real* (see documentation
+   of :doc:`units <units>`).  These can also be read from a potential
+   file with correct unit style by specifying the name of the
+   file. Several potential files for each unit style are included in the
+   ``potentials`` directory of the LAMMPS distribution.
+
 Description
 """""""""""
 
-The *oxrna2* pair styles compute the pairwise-additive parts of the oxDNA force field
+The *oxrna2* pair styles compute the pairwise-additive parts of the
-for coarse-grained modelling of DNA. The effective interaction between the nucleotides consists of potentials for the
+oxDNA force field for coarse-grained modelling of RNA. The effective
-excluded volume interaction *oxrna2/excv*, the stacking *oxrna2/stk*, cross-stacking *oxrna2/xstk*
+interaction between the nucleotides consists of potentials for the
-and coaxial stacking interaction *oxrna2/coaxstk*, electrostatic Debye-Hueckel interaction *oxrna2/dh*
+excluded volume interaction *oxrna2/excv*, the stacking *oxrna2/stk*,
-as well as the hydrogen-bonding interaction *oxrna2/hbond* between complementary pairs of nucleotides on
+cross-stacking *oxrna2/xstk* and coaxial stacking interaction
-opposite strands. Average sequence or sequence-dependent stacking and base-pairing strengths
+*oxrna2/coaxstk*, electrostatic Debye-Hueckel interaction *oxrna2/dh* as
-are supported :ref:`(Sulc2) <Sulc32>`. Quasi-unique base-pairing between nucleotides can be achieved by using
+well as the hydrogen-bonding interaction *oxrna2/hbond* between
-more complementary pairs of atom types like 5-8 and 6-7, 9-12 and 10-11, 13-16 and 14-15, etc.
+complementary pairs of nucleotides on opposite strands. Average sequence
-This prevents the hybridization of in principle complementary bases within Ntypes/4 bases
+or sequence-dependent stacking and base-pairing strengths are supported
+:ref:`(Sulc2) <Sulc32>`. Quasi-unique base-pairing between nucleotides
+can be achieved by using more complementary pairs of atom types like 5-8
+and 6-7, 9-12 and 10-11, 13-16 and 14-15, etc.  This prevents the
+hybridization of in principle complementary bases within Ntypes/4 bases
 up and down along the backbone.
 
-The exact functional form of the pair styles is rather complex.
+The exact functional form of the pair styles is rather complex.  The
-The individual potentials consist of products of modulation factors,
+individual potentials consist of products of modulation factors, which
-which themselves are constructed from a number of more basic potentials
+themselves are constructed from a number of more basic potentials
-(Morse, Lennard-Jones, harmonic angle and distance) as well as quadratic smoothing and modulation terms.
+(Morse, Lennard-Jones, harmonic angle and distance) as well as quadratic
-We refer to :ref:`(Sulc1) <Sulc31>` and the original oxDNA publications :ref:`(Ouldridge-DPhil) <Ouldridge-DPhil3>`
+smoothing and modulation terms.  We refer to :ref:`(Sulc1) <Sulc31>` and
-and  :ref:`(Ouldridge) <Ouldridge3>` for a detailed description of the oxRNA2 force field.
+the original oxDNA publications :ref:`(Ouldridge-DPhil)
+<Ouldridge-DPhil3>` and :ref:`(Ouldridge) <Ouldridge3>` for a detailed
+description of the oxRNA2 force field.
 
 .. note::
 
-   These pair styles have to be used together with the related oxDNA2 bond style
+   These pair styles have to be used together with the related oxDNA2
-   *oxrna2/fene* for the connectivity of the phosphate backbone (see also documentation of
+   bond style *oxrna2/fene* for the connectivity of the phosphate
-   :doc:`bond_style oxrna2/fene <bond_oxdna>`). Most of the coefficients
+   backbone (see also documentation of :doc:`bond_style oxrna2/fene
-   in the above example have to be kept fixed and cannot be changed without reparameterizing the entire model.
+   <bond_oxdna>`). Most of the coefficients in the above example have to
-   Exceptions are the first four coefficients after *oxrna2/stk* (seq=seqdep, T=0.1, xi=1.40206 and kappa=2.77 in the above example),
+   be kept fixed and cannot be changed without reparameterizing the
-   the first coefficient after *oxrna2/hbond* (seq=seqdep in the above example) and the three coefficients
+   entire model.  Exceptions are the first four coefficients after
-   after *oxrna2/dh* (T=0.1, rhos=0.5, qeff=1.02455 in the above example). When using a Langevin thermostat
+   *oxrna2/stk* (seq=seqdep, T=0.1, xi=1.40206 and kappa=2.77 and
-   e.g. through :doc:`fix langevin <fix_langevin>` or :doc:`fix nve/dotc/langevin <fix_nve_dotc_langevin>`
+   corresponding *real unit* equivalents in the above examples), the
-   the temperature coefficients have to be matched to the one used in the fix.
+   first coefficient after *oxrna2/hbond* (seq=seqdep in the above
+   example) and the three coefficients after *oxrna2/dh* (T=0.1,
+   rhos=0.5, qeff=1.02455 in the above example). When using a Langevin
+   thermostat e.g. through :doc:`fix langevin <fix_langevin>` or
+   :doc:`fix nve/dotc/langevin <fix_nve_dotc_langevin>` the temperature
+   coefficients have to be matched to the one used in the fix.
 
 .. note::
 
-   These pair styles have to be used with the *atom_style hybrid bond ellipsoid oxdna*
+   These pair styles have to be used with the *atom_style hybrid bond
-   (see documentation of :doc:`atom_style <atom_style>`). The *atom_style oxdna*
+   ellipsoid oxdna* (see documentation of :doc:`atom_style
-   stores the 3'-to-5' polarity of the nucleotide strand, which is set through
+   <atom_style>`). The *atom_style oxdna* stores the 3'-to-5' polarity
-   the bond topology in the data file. The first (second) atom in a bond definition
+   of the nucleotide strand, which is set through the bond topology in
-   is understood to point towards the 3'-end (5'-end) of the strand.
+   the data file. The first (second) atom in a bond definition is
+   understood to point towards the 3'-end (5'-end) of the strand.
 
-Example input and data files for DNA duplexes can be found in examples/PACKAGES/cgdna/examples/oxDNA/ and /oxDNA2/.
+Example input and data files for DNA duplexes can be found in
-A simple python setup tool which creates single straight or helical DNA strands,
+``examples/PACKAGES/cgdna/examples/oxDNA/`` and ``.../oxDNA2/``.  A simple python
-DNA duplexes or arrays of DNA duplexes can be found in examples/PACKAGES/cgdna/util/.
+setup tool which creates single straight or helical DNA strands, DNA
+duplexes or arrays of DNA duplexes can be found in
+``examples/PACKAGES/cgdna/util/``.
 
 Please cite :ref:`(Henrich) <Henrich3>` in any publication that uses
-this implementation.  The article contains general information
+this implementation.  The article contains general information on the
-on the model, its implementation and performance as well as the structure of
+model, its implementation and performance as well as the structure of
-the data and input file. The preprint version of the article can be found
+the data and input file. The preprint version of the article can be
-`here <PDF/CG-DNA.pdf>`_.
+found `here <PDF/CG-DNA.pdf>`_.  Please cite also the relevant oxRNA2
-Please cite also the relevant oxRNA2 publications
+publications :ref:`(Sulc1) <Sulc31>` and :ref:`(Sulc2) <Sulc32>`.
-:ref:`(Sulc1) <Sulc31>` and :ref:`(Sulc2) <Sulc32>`.
+
+----------
+
+Potential file reading
+""""""""""""""""""""""
+
+For each pair style above the first non-modifiable argument can be a
+filename (with exception of Debye-Hueckel, for which the effective
+charge argument can be a filename), and if it is, no further arguments
+should be supplied.  Therefore the following command:
+
+.. code-block:: LAMMPS
+
+   pair_coeff 3 4 oxrna2/hbond   seqdep oxrna2_lj.cgdna
+
+will be interpreted as a request to read the corresponding hydrogen
+bonding potential parameters from the file with the given name.  The
+file can define multiple potential parameters for both bonded and pair
+interactions, but for the example pair interaction above there must
+exist in the file a line of the form:
+
+.. code-block:: LAMMPS
+
+  3 4 hbond     <coefficients>
+
+If potential customization is required, the potential file reading can
+be mixed with the manual specification of the potential parameters. For
+example, the following command:
+
+.. code-block:: LAMMPS
+
+   pair_style hybrid/overlay oxrna2/excv oxrna2/stk oxrna2/hbond oxrna2/xstk oxrna2/coaxstk oxrna2/dh
+   pair_coeff * * oxrna2/excv    2.0 0.7 0.675 2.0 0.515 0.5 2.0 0.33 0.32
+   pair_coeff * * oxrna2/stk     seqdep 0.1 1.40206 2.77 oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff 1 4 oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff 2 3 oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff 3 4 oxrna2/hbond   seqdep oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/xstk    oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/coaxstk oxrna2_lj.cgdna
+   pair_coeff * * oxrna2/dh      0.1 0.5 1.02455
+
+will read the excluded volume and Debye-Hueckel effective charge *qeff*
+parameters from the manual specification and all others from the
+potential file *oxrna2_lj.cgdna*.
+
+There are sample potential files for each unit style in the
+``potentials`` directory of the LAMMPS distribution. The potential file
+unit system must align with the units defined via the :doc:`units
+<units>` command. For conversion between different *LJ* and *real* unit
+systems for oxDNA, the python tool *lj2real.py* located in the
+``examples/PACKAGES/cgdna/util/`` directory can be used. This tool
+assumes similar file structure to the examples found in
+``examples/PACKAGES/cgdna/examples/``.
 
 ----------
 

doc/src/pair_soft.rst

@@ -1,11 +1,12 @@
 .. index:: pair_style soft
 .. index:: pair_style soft/gpu
+.. index:: pair_style soft/kk
 .. index:: pair_style soft/omp
 
 pair_style soft command
 =======================
 
-Accelerator Variants: *soft/gpu*, *soft/omp*
+Accelerator Variants: *soft/gpu*, *soft/kk*, *soft/omp*
 
 Syntax
 """"""

doc/src/pair_style.rst

@@ -108,6 +108,7 @@ accelerated styles exist.
 
 * :doc:`none <pair_none>` - turn off pairwise interactions
 * :doc:`hybrid <pair_hybrid>` - multiple styles of pairwise interactions
+* :doc:`hybrid/molecular <pair_hybrid>` - different pair styles for intra- and inter-molecular interactions
 * :doc:`hybrid/overlay <pair_hybrid>` - multiple styles of superposed pairwise interactions
 * :doc:`hybrid/scaled <pair_hybrid>` - multiple styles of scaled superposed pairwise interactions
 * :doc:`zero <pair_zero>` - neighbor list but no interactions
@@ -171,6 +172,7 @@ accelerated styles exist.
 * :doc:`coul/wolf <pair_coul>` - Coulomb via Wolf potential
 * :doc:`coul/wolf/cs <pair_cs>` - Coulomb via Wolf potential with core/shell adjustments
 * :doc:`dpd <pair_dpd>` - dissipative particle dynamics (DPD)
+* :doc:`dpd/coul/slater/long <pair_dpd_coul_slater_long>` - dissipative particle dynamics (DPD) with electrostatic interactions
 * :doc:`dpd/ext <pair_dpd_ext>` - generalized force field for DPD
 * :doc:`dpd/ext/tstat <pair_dpd_ext>` - pairwise DPD thermostatting  with generalized force field
 * :doc:`dpd/fdt <pair_dpd_fdt>` - DPD for constant temperature and pressure
@@ -382,6 +384,7 @@ accelerated styles exist.
 * :doc:`tracker <pair_tracker>` - monitor information about pairwise interactions
 * :doc:`tri/lj <pair_tri_lj>` - LJ potential between triangles
 * :doc:`ufm <pair_ufm>` -
+* :doc:`uf3 <pair_uf3>` - UF3 machine-learning potential
 * :doc:`vashishta <pair_vashishta>` - Vashishta 2-body and 3-body potential
 * :doc:`vashishta/table <pair_vashishta>` -
 * :doc:`wf/cut <pair_wf_cut>` - Wang-Frenkel Potential for short-ranged interactions

doc/src/pair_uf3.rst

@@ -0,0 +1,213 @@
+.. index:: pair_style uf3
+.. index:: pair_style uf3/kk
+
+pair_style uf3 command
+======================
+
+Accelerator Variants: *uf3/kk*
+
+Syntax
+""""""
+
+.. code-block:: LAMMPS
+
+    pair_style style BodyFlag
+
+* style = *uf3* or *uf3/kk*
+
+  .. parsed-literal::
+
+       BodyFlag = Indicates whether to calculate only 2-body or 2 and 3-body interactions. Possible values: 2 or 3
+
+Examples
+""""""""
+
+.. code-block:: LAMMPS
+
+    pair_style uf3 3
+    pair_coeff * * Nb.uf3 Nb
+
+    pair_style uf3 2
+    pair_coeff * * NbSn.uf3 Nb Sn
+
+    pair_style uf3 3
+    pair_coeff * * NbSn.uf3 Nb Sn
+
+Description
+"""""""""""
+
+.. versionadded:: TBD
+
+The *uf3* style computes the :ref:`Ultra-Fast Force Fields (UF3)
+<Xie23>` potential, a machine-learning interatomic potential. In UF3,
+the total energy of the system is defined via two- and three-body
+interactions:
+
+.. math::
+
+   E & = \sum_{i,j} V_2(r_{ij}) + \sum_{i,j,k} V_3 (r_{ij},r_{ik},r_{jk}) \\
+   V_2(r_{ij}) & = \sum_{n=0}^N c_n B_n(r_{ij}) \\
+   V_3 (r_{ij},r_{ik},r_{jk}) & = \sum_{l=0}^{N_l} \sum_{m=0}^{N_m} \sum_{n=0}^{N_n} c_{l,m,n} B_l(r_{ij}) B_m(r_{ik}) B_n(r_{jk})
+
+where :math:`V_2(r_{ij})` and :math:`V_3 (r_{ij},r_{ik},r_{jk})` are the
+two- and three-body interactions, respectively. For the two-body the
+summation is over all neighbors J and for the three-body the summation
+is over all neighbors J and K of atom I within a cutoff distance
+determined from the potential files. :math:`B_n(r_{ij})` are the cubic
+b-spline basis, :math:`c_n` and :math:`c_{l,m,n}` are the machine-learned
+interaction parameters and :math:`N`, :math:`N_l`, :math:`N_m`, and
+:math:`N_n` denote the number of basis functions per spline or tensor
+spline dimension.
+
+With *uf3* style only a single pair_coeff command is used to indicate the
+UF3 LAMMPS potential file containing all the two- and three-body interactions
+followed by N additional arguments specifying the mapping of UF3 elements to
+LAMMPS atom types, where N is the number of LAMMPS atom types:
+
+* UF3 LAMMPS potential file
+* N elements names = mapping of UF3 elements to atom types
+
+As an example, if a LAMMPS simulation contains 2 atom types (elements
+'A' and 'B'), the pair_coeff command will be:
+
+.. code-block:: LAMMPS
+
+   pair_style uf3 3
+   pair_coeff * * AB.uf3 A B
+
+The AB.uf3 file should contain all two-body (A-A, A-B, B-B) and three-body
+(A-A-A, A-A-B, A-B-B, B-A-A, B-A-B, B-B-B).
+
+If a value of "2" is specified in the :code:`pair_style uf3` command,
+only the two-body potentials are needed. For 3-body interaction the
+first atom type is the central atom. We recommend using the
+:code:`generate_uf3_lammps_pots.py` script (found `here
+<https://github.com/uf3/uf3/tree/develop/lammps_plugin/scripts>`_) for
+generating the UF3 LAMMPS potential file from the UF3 JSON potentials.
+
+----------
+
+UF3 LAMMPS potential file in the *potentials* directory of the LAMMPS
+distribution have a ".uf3" suffix. The interaction block in UF3 LAMMPS potential
+file should start with :code:`#UF3 POT` and end with :code:`#` characters.
+Following shows the format of a generic 2-body and 3-body potential block in
+UF3 LAMMPS potential file-
+
+.. code-block:: LAMMPS
+
+   #UF3 POT UNITS: units DATE: POT_GEN_DATE AUTHOR: AUTHOR_NAME CITATION: CITE
+   2B ELEMENT1 ELEMENT2 LEADING_TRIM TRAILING_TRIM
+   Rij_CUTOFF NUM_OF_KNOTS
+   BSPLINE_KNOTS
+   NUM_OF_COEFF
+   COEFF
+   #
+   #UF3 POT UNITS: units DATE: POT_GEN_DATE AUTHOR: AUTHOR_NAME CITATION: CITE
+   3B ELEMENT1 ELEMENT2 ELEMENT3 LEADING_TRIM TRAILING_TRIM
+   Rjk_CUTOFF Rik_CUTOFF Rij_CUTOFF NUM_OF_KNOTS_JK NUM_OF_KNOTS_IK NUM_OF_KNOTS_IJ
+   BSPLINE_KNOTS_FOR_JK
+   BSPLINE_KNOTS_FOR_IK
+   BSPLINE_KNOTS_FOR_IJ
+   SHAPE_OF_COEFF_MATRIX[I][J][K]
+   COEFF_MATRIX[0][0][K]
+   COEFF_MATRIX[0][1][K]
+   COEFF_MATRIX[0][2][K]
+   .
+   .
+   .
+   COEFF_MATRIX[1][0][K]
+   COEFF_MATRIX[1][1][K]
+   COEFF_MATRIX[1][2][K]
+   .
+   .
+   .
+   #
+
+The second line indicates whether the block contains data for 2-body
+(:code:`2B`) or 3-body (:code:`3B`) interaction. This is followed by element
+combination interaction, :code:`LEADING_TRIM` and :code:`TRAILING_TRIM`
+number on the same line. The current implementation is only tested for
+:code:`LEADING_TRIM=0` and :code:`TRAILING_TRIM=3`.
+If other values are used LAMMPS is terminated after issuing an error message.
+The :code:`Rij_CUTOFF` sets the 2-body cutoff for the interaction described
+by the potential block. :code:`NUM_OF_KNOTS` is the number of knots
+(or the length of the knot vector) present on the very next line. The
+:code:`BSPLINE_KNOTS` line should contain all the knots in ascending order.
+:code:`NUM_OF_COEFF` is the number of coefficients in the :code:`COEFF` line.
+All the numbers in the BSPLINE_KNOTS and COEFF line should be space-separated.
+Similar to the 2-body potential block, the third line sets the cutoffs and
+length of the knots. The cutoff distance between atom-type I and J is
+:code:`Rij_CUTOFF`, atom-type I and K is :code:`Rik_CUTOFF` and between
+J and K is :code:`Rjk_CUTOFF`.
+
+.. note::
+
+   The current implementation only works for UF3 potentials with cutoff
+   distances for 3-body interactions that follows
+   :code:`2Rij_CUTOFF=2Rik_CUTOFF=Rjk_CUTOFF` relation.
+
+The :code:`BSPLINE_KNOTS_FOR_JK`, :code:`BSPLINE_KNOTS_FOR_IK`, and
+:code:`BSPLINE_KNOTS_FOR_IJ` lines (note the order) contain the knots in
+increasing order for atoms J and K, I and K, and atoms I and J
+respectively. The number of knots is defined by the
+:code:`NUM_OF_KNOTS_*` characters in the previous line.  The shape of
+the coefficient matrix is defined on the
+:code:`SHAPE_OF_COEFF_MATRIX[I][J][K]` line followed by the columns of
+the coefficient matrix, one per line, as shown above. For example, if
+the coefficient matrix has the shape of 8x8x13, then
+:code:`SHAPE_OF_COEFF_MATRIX[I][J][K]` will be :code:`8 8 13` followed
+by 64 (8x8) lines each containing 13 coefficients separated by space.
+
+----------
+
+.. include:: accel_styles.rst
+
+----------
+
+Mixing, shift, table, tail correction, restart, rRESPA info
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+For atom type pairs I,J and I != J, where types I and J correspond to
+two different element types, mixing is performed by LAMMPS as described
+above from values in the potential file.
+
+This pair style does not support the :doc:`pair_modify <pair_modify>`
+shift, table, and tail options.
+
+This pair style does not write its information to :doc:`binary restart
+files <restart>`, since it is stored in potential file.
+
+This pair style can only be used via the *pair* keyword of the
+:doc:`run_style respa <run_style>` command.  It does not support the
+*inner*, *middle*, *outer* keywords.
+
+Restrictions
+""""""""""""
+
+The 'uf3' pair style is part of the ML-UF3 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 :doc:`newton <newton>` setting to be "on".
+
+The UF3 LAMMPS potential file provided with LAMMPS (see the potentials
+directory) are parameterized for metal :doc:`units <units>`.
+
+The single() function of 'uf3' pair style only return the 2-body
+interaction energy.
+
+Related commands
+""""""""""""""""
+
+:doc:`pair_coeff <pair_coeff>`
+
+Default
+"""""""
+
+none
+
+----------
+
+.. _Xie23:
+
+**(Xie23)** Xie, S.R., Rupp, M. & Hennig, R.G. Ultra-fast interpretable machine-learning potentials. npj Comput Mater 9, 162 (2023). https://doi.org/10.1038/s41524-023-01092-7

doc/src/run_style.rst

@@ -329,8 +329,10 @@ Restrictions
 The *verlet/split* style can only be used if LAMMPS was built with the
 REPLICA package.  Correspondingly the *respa/omp* style is available
 only if the OPENMP package was included. See the :doc:`Build package
-<Build_package>` page for more info.  It is not compatible with
+<Build_package>` page for more info.
-kspace styles from the INTEL package.
+
+Run style *verlet/split* is not compatible with kspace styles from
+the INTEL package and it is not compatible with any TIP4P styles.
 
 Whenever using rRESPA, the user should experiment with trade-offs in
 speed and accuracy for their system, and verify that they are

doc/src/variable.rst

@@ -67,7 +67,7 @@ Syntax
                            bound(group,dir,region), gyration(group,region), ke(group,reigon),
                            angmom(group,dim,region), torque(group,dim,region),
                            inertia(group,dimdim,region), omega(group,dim,region)
-         special functions = sum(x), min(x), max(x), ave(x), trap(x), slope(x), gmask(x), rmask(x), grmask(x,y), next(x), is_file(name), is_os(name), extract_setting(name), label2type(kind,label), is_typelabel(kind,label)
+         special functions = sum(x), min(x), max(x), ave(x), trap(x), slope(x), sort(x), rsort(x), gmask(x), rmask(x), grmask(x,y), next(x), is_file(name), is_os(name), extract_setting(name), label2type(kind,label), is_typelabel(kind,label)
          feature functions = is_available(category,feature), is_active(category,feature), is_defined(category,id)
          atom value = id[i], mass[i], type[i], mol[i], x[i], y[i], z[i], vx[i], vy[i], vz[i], fx[i], fy[i], fz[i], q[i]
          atom vector = id, mass, type, mol, radius, q, x, y, z, vx, vy, vz, fx, fy, fz
@@ -547,7 +547,7 @@ variables.
 +------------------------+---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------+
 | Region functions       | count(ID,IDR), mass(ID,IDR), charge(ID,IDR), xcm(ID,dim,IDR), vcm(ID,dim,IDR), fcm(ID,dim,IDR), bound(ID,dir,IDR), gyration(ID,IDR), ke(ID,IDR), angmom(ID,dim,IDR), torque(ID,dim,IDR), inertia(ID,dimdim,IDR), omega(ID,dim,IDR)                                                                                                                |
 +------------------------+---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------+
-| Special functions      | sum(x), min(x), max(x), ave(x), trap(x), slope(x), gmask(x), rmask(x), grmask(x,y), next(x), is_file(name), is_os(name), extract_setting(name), label2type(kind,label), is_typelabel(kind,label)                                                                                                                                                  |
+| Special functions      | sum(x), min(x), max(x), ave(x), trap(x), slope(x), sort(x), rsort(x), gmask(x), rmask(x), grmask(x,y), next(x), is_file(name), is_os(name), extract_setting(name), label2type(kind,label), is_typelabel(kind,label)                                                                                                                               |
 +------------------------+---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------+
 | Feature functions      | is_available(category,feature), is_active(category,feature), is_defined(category,id)                                                                                                                                                                                                                                                              |
 +------------------------+---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------+
@@ -913,23 +913,27 @@ Special Functions
 Special functions take specific kinds of arguments, meaning their
 arguments cannot be formulas themselves.
 
-The sum(x), min(x), max(x), ave(x), trap(x), and slope(x) functions
+The sum(x), min(x), max(x), ave(x), trap(x), slope(x), sort(x), and
-each take 1 argument which is of the form "c_ID" or "c_ID[N]" or
+rsort(x) functions each take 1 argument which is of the form "c_ID" or
-"f_ID" or "f_ID[N]" or "v_name".  The first two are computes and the
+"c_ID[N]" or "f_ID" or "f_ID[N]" or "v_name".  The first two are
-second two are fixes; the ID in the reference should be replaced by
+computes and the second two are fixes; the ID in the reference should be
-the ID of a compute or fix defined elsewhere in the input script.  The
+replaced by the ID of a compute or fix defined elsewhere in the input
-compute or fix must produce either a global vector or array.  If it
+script.  The compute or fix must produce either a global vector or
-produces a global vector, then the notation without "[N]" should be
+array.  If it produces a global vector, then the notation without "[N]"
-used.  If it produces a global array, then the notation with "[N]"
+should be used.  If it produces a global array, then the notation with
-should be used, when N is an integer, to specify which column of the
+"[N]" should be used, where N is an integer, to specify which column of
-global array is being referenced.  The last form of argument "v_name"
+the global array is being referenced.  The last form of argument
-is for a vector-style variable where "name" is replaced by the name of
+"v_name" is for a vector-style variable where "name" is replaced by the
-the variable.
+name of the variable.
 
-These functions operate on a global vector of inputs and reduce it to
+The sum(x), min(x), max(x), ave(x), trap(x), and slope(x) functions
-a single scalar value.  This is analogous to the operation of the
+operate on a global vector of inputs and reduce it to a single scalar
-:doc:`compute reduce <compute_reduce>` command, which performs similar
+value.  This is analogous to the operation of the :doc:`compute reduce
-operations on per-atom and local vectors.
+<compute_reduce>` command, which performs similar operations on per-atom
+and local vectors.
+
+The sort(x) and rsort(x) functions operate on a global vector of inputs
+and return a global vector of the same length.
 
 The sum() function calculates the sum of all the vector elements.  The
 min() and max() functions find the minimum and maximum element
@@ -953,6 +957,12 @@ of points, equally spaced by 1 in their x coordinate: (1,V1), (2,V2),
 length N.  The returned value is the slope of the line.  If the line
 has a single point or is vertical, it returns 1.0e20.
 
+.. versionadded:: TBD
+
+The sort(x) and rsort(x) functions sort the data of the input vector by
+their numeric value: sort(x) sorts in ascending order, rsort(x) sorts
+in descending order.
+
 The gmask(x) function takes 1 argument which is a group ID.  It
 can only be used in atom-style variables.  It returns a 1 for
 atoms that are in the group, and a 0 for atoms that are not.

doc/src/write_data.rst

@@ -189,4 +189,4 @@ Related commands
 Default
 """""""
 
-The option defaults are pair = ii and types_style = numeric.
+The option defaults are pair = ii and types = numeric.

doc/utils/requirements.txt

@@ -1,4 +1,4 @@
-Sphinx >= 5.3.0, <7.3
+Sphinx >= 5.3.0, <8.0
 sphinxcontrib-spelling
 sphinxcontrib-jquery
 git+https://github.com/akohlmey/sphinx-fortran@parallel-read

doc/utils/sphinx-config/false_positives.txt

@@ -992,6 +992,7 @@ emax
 Emax
 Embt
 emi
+Emilie
 Emmrich
 emol
 eN
@@ -1433,6 +1434,7 @@ Hendrik
 Henin
 Henkelman
 Henkes
+Hennig
 henrich
 Henrich
 Hermitian
@@ -1594,6 +1596,7 @@ interlayer
 intermolecular
 interoperable
 Interparticle
+interpretable
 interstitials
 intertube
 Intr
@@ -1732,6 +1735,7 @@ Kalia
 Kamberaj
 Kantorovich
 Kapfer
+Kapil
 Karhunen
 Karls
 Karlsruhe
@@ -1763,8 +1767,10 @@ keflag
 Keir
 Kelchner
 Kelkar
+Kemppainen
 Kemper
 kepler
+Kemppainen
 keV
 Keyes
 keyfile
@@ -1813,6 +1819,7 @@ Koziol
 Kp
 kradius
 Kraker
+Krass
 Kraus
 Kremer
 Kress
@@ -2483,6 +2490,7 @@ Nevery
 newfile
 Newns
 newtype
+nextsort
 Neyts
 Nf
 nfft
@@ -2672,6 +2680,7 @@ nzlo
 ocl
 octahedral
 octants
+Odegard
 Ohara
 O'Hearn
 ohenrich
@@ -2953,6 +2962,7 @@ Priya
 proc
 Proc
 procs
+procgrid
 progguide
 Prony
 ps
@@ -3252,6 +3262,7 @@ rRESPA
 Rsi
 Rso
 Rspace
+rsort
 rsq
 rst
 rstyle
@@ -3264,6 +3275,7 @@ Rudranarayan
 Rudzinski
 Runge
 runtime
+Rupp
 Rutuparna
 rx
 rxd
@@ -3287,6 +3299,7 @@ Saidi
 saip
 Salanne
 Salles
+sametag
 sandia
 Sandia
 sandybrown
@@ -3404,6 +3417,7 @@ sinh
 sinusoid
 sinusoidally
 SiO
+Siochi
 Sirk
 Sival
 sizeI
@@ -3451,6 +3465,7 @@ solvated
 solvation
 someuser
 Sorensen
+sortfreq
 soundspeed
 sourceforge
 Souza
@@ -3810,6 +3825,8 @@ uChem
 uCond
 uef
 UEF
+uf
+uf3
 ufm
 Uhlenbeck
 Ui
@@ -4149,6 +4166,7 @@ yy
 yz
 Zagaceta
 Zannoni
+Zavada
 Zavattieri
 zbl
 ZBL

examples/COUPLE/plugin/liblammpsplugin.c

@@ -101,6 +101,7 @@ liblammpsplugin_t *liblammpsplugin_load(const char *lib)
   ADDSYM(extract_setting);
   ADDSYM(extract_global_datatype);
   ADDSYM(extract_global);
+  ADDSYM(map_atom);
 
   ADDSYM(extract_atom_datatype);
   ADDSYM(extract_atom);

examples/COUPLE/plugin/liblammpsplugin.h

@@ -146,6 +146,7 @@ struct _liblammpsplugin {
   int (*extract_setting)(void *, const char *);
   int *(*extract_global_datatype)(void *, const char *);
   void *(*extract_global)(void *, const char *);
+  void *(*map_atom)(void *, const void *);
 
   int *(*extract_atom_datatype)(void *, const char *);
   void *(*extract_atom)(void *, const char *);

examples/PACKAGES/cgdna/examples/lj_units/oxDNA/potential_file/in.duplex1

@@ -0,0 +1,70 @@
+variable number	equal 1
+variable ofreq	equal 1000
+variable efreq	equal 1000
+variable T      equal 0.1
+
+units lj
+
+dimension 3
+
+newton on
+
+boundary  p p p
+
+atom_style hybrid bond ellipsoid oxdna
+atom_modify sort 0 1.0
+
+# Pair interactions require lists of neighbours to be calculated
+neighbor 2.0 bin
+neigh_modify every 1 delay 0 check yes
+
+read_data data.duplex1
+
+set atom * mass 3.1575
+
+group all type 1 4
+
+# oxDNA bond interactions - FENE backbone
+bond_style oxdna/fene
+bond_coeff * oxdna_lj.cgdna
+special_bonds lj 0 1 1
+
+# oxDNA pair interactions
+pair_style hybrid/overlay oxdna/excv oxdna/stk oxdna/hbond oxdna/xstk oxdna/coaxstk
+pair_coeff * * oxdna/excv    oxdna_lj.cgdna
+pair_coeff * * oxdna/stk     seqav 0.1 1.3448 2.6568 oxdna_lj.cgdna
+pair_coeff * * oxdna/hbond   seqav oxdna_lj.cgdna
+pair_coeff 1 4 oxdna/hbond   seqav oxdna_lj.cgdna
+pair_coeff 2 3 oxdna/hbond   seqav oxdna_lj.cgdna
+pair_coeff * * oxdna/xstk    oxdna_lj.cgdna
+pair_coeff * * oxdna/coaxstk oxdna_lj.cgdna
+
+# NVE ensemble
+fix 1 all nve/asphere
+#fix 2 all langevin ${T} ${T} 2.5 457145 angmom 10
+
+timestep 1e-5
+
+#comm_style tiled
+fix 3 all balance 1000 1.03 shift xyz 10 1.03
+comm_modify cutoff 3.8
+
+compute quat all property/atom quatw quati quatj quatk
+
+compute erot all erotate/asphere
+compute ekin all ke
+compute epot all pe
+variable erot equal c_erot
+variable ekin equal c_ekin
+variable epot equal c_epot
+variable etot equal c_erot+c_ekin+c_epot
+fix 5 all print ${efreq} "$(step)  ekin = ${ekin} |  erot = ${erot} | epot = ${epot} | etot = ${etot}" screen yes
+
+dump out all custom ${ofreq} out.${number}.lammpstrj id mol type x y z ix iy iz vx vy vz c_quat[1] c_quat[2] c_quat[3] c_quat[4] angmomx angmomy angmomz
+dump_modify out sort id
+dump_modify out format line "%d %d %d  %22.15le %22.15le %22.15le  %d %d %d  %22.15le %22.15le %22.15le  %22.15le %22.15le %22.15le %22.15le  %22.15le %22.15le %22.15le"
+
+run 1000000
+
+write_data last_config.${number}.* nocoeff
+#write_restart last_config.${number}.*

examples/PACKAGES/cgdna/examples/lj_units/oxDNA/potential_file/oxdna_lj.cgdna

@@ -0,0 +1 @@
+../../../../../../../potentials/oxdna_lj.cgdna

examples/PACKAGES/cgdna/examples/lj_units/oxDNA2/duplex1/data.duplex1

@@ -0,0 +1,68 @@
+LAMMPS data file via write_data, version 27 May 2021
+
+10 atoms
+4 atom types
+8 bonds
+1 bond types
+10 ellipsoids
+
+-20 20 xlo xhi
+-20 20 ylo yhi
+-20 20 zlo zhi
+
+Masses
+
+1 3.1575
+2 3.1575
+3 3.1575
+4 3.1575
+
+Atoms # hybrid
+
+1 1 -0.33741452300167507 -0.43708835412476305 0.6450685042019271 1 1 3.7269849963023267 0 0 0
+2 2 -0.32142606102826937 -0.7137743037592722 1.1817366147004618 1 1 3.7269849963023267 0 0 0
+3 3 -0.130363628207774 -0.9147144801536078 1.62581312195109 1 1 3.7269849963023267 0 0 0
+4 4 0.16795127962282844 -0.9808507459807022 2.0894908590909003 1 1 3.7269849963023267 0 0 0
+5 1 0.46370423490634166 -0.7803347954883079 2.4251986815515827 1 1 3.7269849963023267 0 0 0
+6 4 -0.4462950185476711 0.09062163051035639 2.4668941268777607 2 1 3.7269849963023267 0 0 0
+7 1 -0.03377054097560965 0.20979847489755046 2.078208732038921 2 1 3.7269849963023267 0 0 0
+8 2 0.3297325391466579 0.17657587120899895 1.7206328374934152 2 1 3.7269849963023267 0 0 0
+9 3 0.6063699309305985 0.04682595158675571 1.2335049647817748 2 1 3.7269849963023267 0 0 0
+10 4 0.8003979559814726 -0.364393011459011 0.9884025318908612 2 1 3.7269849963023267 0 0 0
+
+Velocities
+
+1 0.320321385294804 -0.13632815939410442 -0.029398292452023418 0.3064009492028237 -0.15808560233691588 0.35878007201886397
+2 0.16868594667473025 -0.04950805613064363 0.15811007290373785 -0.07666583909321756 -0.0008074676325318194 -0.21475821019816385
+3 -0.22924557018300165 0.08381876748892438 -0.0919832851533896 0.4039387481683193 0.6123610642545824 -0.11063432848545783
+4 -0.22186204313310393 0.04952817499985707 -0.0693642101605732 -0.1358248430264938 0.4118493572385653 -0.056529305922687775
+5 0.08156431270087049 -0.2564594759800144 0.1724544416027875 0.05439894663158808 0.09338481510384318 0.3205408219238883
+6 0.03598698404367743 -0.04868642973674152 0.02860105256592922 0.04007709957283992 -0.317943400069374 0.36438025397586354
+7 -0.00822868972307372 0.047514932936351305 -0.027726409420297023 0.18356392696891796 -0.49877294396308003 0.06993313839189567
+8 -0.07177147672242379 0.1718272727853115 0.39056151182616994 -0.16728362538690794 -0.47839708820957955 -0.17897249005947627
+9 -0.1748638855727651 -0.0781638161351808 0.0560181565271157 -0.28062568580131014 0.2405396522734162 -0.4311598357169048
+10 0.18870318272756448 -0.1066780134639517 0.12610657946741227 -0.05740397100183697 0.36748833227892685 0.1498775724372025
+
+Bonds
+
+1 1 1 2
+2 1 2 3
+3 1 3 4
+4 1 4 5
+5 1 6 7
+6 1 7 8
+7 1 8 9
+8 1 9 10
+
+Ellipsoids
+
+1 1.173984503142341 1.173984503142341 1.173984503142341 0.9890278201757743 0.01779228232037064 -0.14337734159225404 0.030827642240801516
+2 1.173984503142341 1.173984503142341 1.173984503142341 0.939687458852748 0.04174166924055095 -0.023337773785056866 0.338674565089608
+3 1.173984503142341 1.173984503142341 1.173984503142341 0.8210113150655425 0.03012140921736572 0.017666019956944813 0.5698429897612057
+4 1.173984503142341 1.173984503142341 1.173984503142341 0.6623662858285051 -0.028186343967346823 0.022942552517501488 0.7482981175276918
+5 1.173984503142341 1.173984503142341 1.173984503142341 0.3601488726765216 0.0513614985821682 0.0724224158335286 0.9286602067807472
+6 1.173984503142341 1.173984503142341 1.173984503142341 0.11941234710084649 0.9244660117493703 -0.35317942248051865 -0.07979711784524246
+7 1.173984503142341 1.173984503142341 1.173984503142341 -0.17949125421205164 0.7412884899431119 -0.6379094464220707 0.1065166771202199
+8 1.173984503142341 1.173984503142341 1.173984503142341 -0.10483691088405202 0.5508895999584645 -0.8250090480220789 0.06992811634525403
+9 1.173984503142341 1.173984503142341 1.173984503142341 0.07777239911646 -0.3724087549185288 0.9103052384821374 -0.1631181963720798
+10 1.173984503142341 1.173984503142341 1.173984503142341 0.16279109707978262 0.027148630125149613 0.9849325709665359 -0.0516705065113425

examples/PACKAGES/cgdna/examples/lj_units/oxDNA2/potential_file/data.duplex1

@@ -0,0 +1,68 @@
+LAMMPS data file via write_data, version 27 May 2021
+
+10 atoms
+4 atom types
+8 bonds
+1 bond types
+10 ellipsoids
+
+-20 20 xlo xhi
+-20 20 ylo yhi
+-20 20 zlo zhi
+
+Masses
+
+1 3.1575
+2 3.1575
+3 3.1575
+4 3.1575
+
+Atoms # hybrid
+
+1 1 -0.33741452300167507 -0.43708835412476305 0.6450685042019271 1 1 3.7269849963023267 0 0 0
+2 2 -0.32142606102826937 -0.7137743037592722 1.1817366147004618 1 1 3.7269849963023267 0 0 0
+3 3 -0.130363628207774 -0.9147144801536078 1.62581312195109 1 1 3.7269849963023267 0 0 0
+4 4 0.16795127962282844 -0.9808507459807022 2.0894908590909003 1 1 3.7269849963023267 0 0 0
+5 1 0.46370423490634166 -0.7803347954883079 2.4251986815515827 1 1 3.7269849963023267 0 0 0
+6 4 -0.4462950185476711 0.09062163051035639 2.4668941268777607 2 1 3.7269849963023267 0 0 0
+7 1 -0.03377054097560965 0.20979847489755046 2.078208732038921 2 1 3.7269849963023267 0 0 0
+8 2 0.3297325391466579 0.17657587120899895 1.7206328374934152 2 1 3.7269849963023267 0 0 0
+9 3 0.6063699309305985 0.04682595158675571 1.2335049647817748 2 1 3.7269849963023267 0 0 0
+10 4 0.8003979559814726 -0.364393011459011 0.9884025318908612 2 1 3.7269849963023267 0 0 0
+
+Velocities
+
+1 0.320321385294804 -0.13632815939410442 -0.029398292452023418 0.3064009492028237 -0.15808560233691588 0.35878007201886397
+2 0.16868594667473025 -0.04950805613064363 0.15811007290373785 -0.07666583909321756 -0.0008074676325318194 -0.21475821019816385
+3 -0.22924557018300165 0.08381876748892438 -0.0919832851533896 0.4039387481683193 0.6123610642545824 -0.11063432848545783
+4 -0.22186204313310393 0.04952817499985707 -0.0693642101605732 -0.1358248430264938 0.4118493572385653 -0.056529305922687775
+5 0.08156431270087049 -0.2564594759800144 0.1724544416027875 0.05439894663158808 0.09338481510384318 0.3205408219238883
+6 0.03598698404367743 -0.04868642973674152 0.02860105256592922 0.04007709957283992 -0.317943400069374 0.36438025397586354
+7 -0.00822868972307372 0.047514932936351305 -0.027726409420297023 0.18356392696891796 -0.49877294396308003 0.06993313839189567
+8 -0.07177147672242379 0.1718272727853115 0.39056151182616994 -0.16728362538690794 -0.47839708820957955 -0.17897249005947627
+9 -0.1748638855727651 -0.0781638161351808 0.0560181565271157 -0.28062568580131014 0.2405396522734162 -0.4311598357169048
+10 0.18870318272756448 -0.1066780134639517 0.12610657946741227 -0.05740397100183697 0.36748833227892685 0.1498775724372025
+
+Bonds
+
+1 1 1 2
+2 1 2 3
+3 1 3 4
+4 1 4 5
+5 1 6 7
+6 1 7 8
+7 1 8 9
+8 1 9 10
+
+Ellipsoids
+
+1 1.173984503142341 1.173984503142341 1.173984503142341 0.9890278201757743 0.01779228232037064 -0.14337734159225404 0.030827642240801516
+2 1.173984503142341 1.173984503142341 1.173984503142341 0.939687458852748 0.04174166924055095 -0.023337773785056866 0.338674565089608
+3 1.173984503142341 1.173984503142341 1.173984503142341 0.8210113150655425 0.03012140921736572 0.017666019956944813 0.5698429897612057
+4 1.173984503142341 1.173984503142341 1.173984503142341 0.6623662858285051 -0.028186343967346823 0.022942552517501488 0.7482981175276918
+5 1.173984503142341 1.173984503142341 1.173984503142341 0.3601488726765216 0.0513614985821682 0.0724224158335286 0.9286602067807472
+6 1.173984503142341 1.173984503142341 1.173984503142341 0.11941234710084649 0.9244660117493703 -0.35317942248051865 -0.07979711784524246
+7 1.173984503142341 1.173984503142341 1.173984503142341 -0.17949125421205164 0.7412884899431119 -0.6379094464220707 0.1065166771202199
+8 1.173984503142341 1.173984503142341 1.173984503142341 -0.10483691088405202 0.5508895999584645 -0.8250090480220789 0.06992811634525403
+9 1.173984503142341 1.173984503142341 1.173984503142341 0.07777239911646 -0.3724087549185288 0.9103052384821374 -0.1631181963720798
+10 1.173984503142341 1.173984503142341 1.173984503142341 0.16279109707978262 0.027148630125149613 0.9849325709665359 -0.0516705065113425

examples/PACKAGES/cgdna/examples/lj_units/oxDNA2/potential_file/in.duplex1

@@ -0,0 +1,72 @@
+variable number equal 1
+variable ofreq  equal 1000
+variable efreq  equal 1000
+variable T      equal 0.1
+variable rhos   equal 0.2
+
+units lj
+
+dimension 3
+
+newton on
+
+boundary  p p p
+
+atom_style hybrid bond ellipsoid oxdna
+atom_modify sort 0 1.0
+
+# Pair interactions require lists of neighbours to be calculated
+neighbor 2.0 bin
+neigh_modify every 1 delay 0 check yes
+
+read_data data.duplex1
+
+set atom * mass 3.1575
+
+group all type 1 4
+
+# oxDNA2 bond interactions - FENE backbone
+bond_style oxdna2/fene
+bond_coeff * oxdna2_lj.cgdna
+special_bonds lj 0 1 1
+
+# oxDNA2 pair interactions
+pair_style hybrid/overlay oxdna2/excv oxdna2/stk oxdna2/hbond oxdna2/xstk oxdna2/coaxstk oxdna2/dh
+pair_coeff * * oxdna2/excv    oxdna2_lj.cgdna
+pair_coeff * * oxdna2/stk     seqdep 0.1 1.3523 2.6717 oxdna2_lj.cgdna
+pair_coeff * * oxdna2/hbond   seqdep oxdna2_lj.cgdna
+pair_coeff 1 4 oxdna2/hbond   seqdep oxdna2_lj.cgdna
+pair_coeff 2 3 oxdna2/hbond   seqdep oxdna2_lj.cgdna
+pair_coeff * * oxdna2/xstk    oxdna2_lj.cgdna
+pair_coeff * * oxdna2/coaxstk oxdna2_lj.cgdna
+pair_coeff * * oxdna2/dh      0.1 0.5 oxdna2_lj.cgdna
+
+# NVE ensemble
+fix 1 all nve/asphere
+#fix 2 all langevin ${T} ${T} 2.5 457145 angmom 10
+
+timestep 1e-5
+
+#comm_style tiled
+fix 3 all balance 1000 1.03 shift xyz 10 1.03
+comm_modify cutoff 3.8
+
+compute quat all property/atom quatw quati quatj quatk
+
+compute erot all erotate/asphere
+compute ekin all ke
+compute epot all pe
+variable erot equal c_erot
+variable ekin equal c_ekin
+variable epot equal c_epot
+variable etot equal c_erot+c_ekin+c_epot
+fix 5 all print ${efreq} "$(step)  ekin = ${ekin} |  erot = ${erot} | epot = ${epot} | etot = ${etot}" screen yes
+
+dump out all custom ${ofreq} out.${number}.lammpstrj id mol type x y z ix iy iz vx vy vz c_quat[1] c_quat[2] c_quat[3] c_quat[4] angmomx angmomy angmomz
+dump_modify out sort id
+dump_modify out format line "%d %d %d  %22.15le %22.15le %22.15le  %d %d %d  %22.15le %22.15le %22.15le  %22.15le %22.15le %22.15le %22.15le  %22.15le %22.15le %22.15le"
+
+run 1000000
+
+write_data last_config.${number}.* nocoeff
+#write_restart last_config.${number}.*