Merge pull request #3375 from akohlmey/next_patch_release

Step version strings for the next patch release

.gitattributes

@@ -3,6 +3,7 @@
 .github export-ignore
 .lgtm.yml export-ignore
 SECURITY.md export-ignore
+CITATION.cff export-ignore
 *       text=auto
 *.jpg   -text
 *.pdf   -text

.github/codecov.yml

@@ -7,7 +7,7 @@ coverage:
         threshold: 10%
         only_pulls: false
         branches:
-          - "master"
+          - "develop"
         flags:
           - "unit"
         paths:
@@ -16,14 +16,14 @@ coverage:
     project:
         default:
             branches:
-              - "master"
+              - "develop"
             paths:
               - "src"
             informational: true
     patch:
         default:
             branches:
-              - "master"
+              - "develop"
             paths:
               - "src"
             informational: true

CITATION.cff

@@ -0,0 +1,91 @@
+# YAML 1.2
+---
+cff-version: 1.2.0
+title: "LAMMPS: Large-scale Atomic/Molecular Massively Parallel Simulator"
+type: software
+authors:
+  - family-names: "Plimpton"
+    given-names: "Steven J."
+  - family-names: "Kohlmeyer"
+    given-names: "Axel"
+    orcid: "https://orcid.org/0000-0001-6204-6475"
+  - family-names: "Thompson"
+    given-names: "Aidan P."
+    orcid: "https://orcid.org/0000-0002-0324-9114"
+  - family-names: "Moore"
+    given-names: "Stan G."
+  - family-names: "Berger"
+    given-names: "Richard"
+    orcid: "https://orcid.org/0000-0002-3044-8266"
+doi: 10.5281/zenodo.3726416
+license: GPL-2.0-only
+url: https://www.lammps.org
+repository-code: https://github.com/lammps/lammps/
+keywords:
+  - "Molecular Dynamics"
+  - "Materials Modeling"
+message: "If you are referencing LAMMPS in a publication, please cite the paper below."
+preferred-citation:
+  type: article
+  doi: "10.1016/j.cpc.2021.108171"
+  url: "https://www.sciencedirect.com/science/article/pii/S0010465521002836"
+  authors:
+  - family-names: "Thompson"
+    given-names: "Aidan P."
+    orcid: "https://orcid.org/0000-0002-0324-9114"
+  - family-names: "Aktulga"
+    given-names: "H. Metin"
+  - family-names: "Berger"
+    given-names: "Richard"
+    orcid: "https://orcid.org/0000-0002-3044-8266"
+  - family-names: "Bolintineanu"
+    given-names: "Dan S."
+  - family-names: "Brown"
+    given-names: "W. Michael"
+  - family-names: "Crozier"
+    given-names: "Paul S."
+  - family-names: "in 't Veld"
+    given-names: "Pieter J."
+  - family-names: "Kohlmeyer"
+    given-names: "Axel"
+    orcid: "https://orcid.org/0000-0001-6204-6475"
+  - family-names: "Moore"
+    given-names: "Stan G."
+  - family-names: "Nguyen"
+    given-names: "Trung Dac"
+  - family-names: "Shan"
+    given-names: "Ray"
+  - family-names: "Stevens"
+    given-names: "Mark J."
+  - family-names: "Tranchida"
+    given-names: "Julien"
+  - family-names: "Trott"
+    given-names: "Christian"
+  - family-names: "Plimpton"
+    given-names: "Steven J."
+  title: "LAMMPS - a flexible simulation tool for particle-based materials modeling at the atomic, meso, and continuum scales"
+  journal: "Computer Physics Communications"
+  keywords:
+   - Molecular dynamics
+   - Materials modeling
+   - Parallel algorithms
+   - LAMMPS
+  month: 2
+  volume: 271
+  issn: 0010-4655
+  pages: 108171
+  year: 2022
+references:
+ - title: "Fast Parallel Algorithms for Short-Range Molecular Dynamics"
+   type: article
+   journal: Journal of Computational Physics
+   volume: 117
+   number: 1
+   pages: "1-19"
+   year: 1995
+   issn: 0021-9991
+   doi: 10.1006/jcph.1995.1039
+   url: https://www.sciencedirect.com/science/article/pii/S002199918571039X
+   authors:
+   - family-names: "Plimpton"
+     given-names: "Steve"

cmake/CMakeLists.txt

@@ -154,6 +154,19 @@ endif()
 ########################################################################
 # User input options                                                   #
 ########################################################################
+# set path to python interpreter and thus enforcing python version if
+# when in a virtual environment and PYTHON_EXECUTABLE is not set on command line
+if(DEFINED ENV{VIRTUAL_ENV} AND NOT PYTHON_EXECUTABLE)
+  if(CMAKE_HOST_SYSTEM_NAME STREQUAL "Windows")
+    set(PYTHON_EXECUTABLE "$ENV{VIRTUAL_ENV}/Scripts/python.exe")
+  else()
+    set(PYTHON_EXECUTABLE "$ENV{VIRTUAL_ENV}/bin/python")
+  endif()
+  set(Python_EXECUTABLE "${PYTHON_EXECUTABLE}")
+  message(STATUS "Running in virtual environment: $ENV{VIRTUAL_ENV}\n"
+    "   Setting Python interpreter to: ${PYTHON_EXECUTABLE}")
+endif()
+
 set(LAMMPS_MACHINE "" CACHE STRING "Suffix to append to lmp binary (WON'T enable any features automatically")
 mark_as_advanced(LAMMPS_MACHINE)
 if(LAMMPS_MACHINE)
@@ -194,6 +207,7 @@ option(CMAKE_VERBOSE_MAKEFILE "Generate verbose Makefiles" OFF)
 
 set(STANDARD_PACKAGES
   ADIOS
+  AMOEBA
   ASPHERE
   ATC
   AWPMD
@@ -202,7 +216,7 @@ set(STANDARD_PACKAGES
   BPM
   BROWNIAN
   CG-DNA
-  CG-SDK
+  CG-SPICA
   CLASS2
   COLLOID
   COLVARS
@@ -358,6 +372,7 @@ pkg_depends(MPIIO MPI)
 pkg_depends(ATC MANYBODY)
 pkg_depends(LATBOLTZ MPI)
 pkg_depends(SCAFACOS MPI)
+pkg_depends(AMOEBA KSPACE)
 pkg_depends(DIELECTRIC KSPACE)
 pkg_depends(DIELECTRIC EXTRA-PAIR)
 pkg_depends(CG-DNA MOLECULE)
@@ -402,9 +417,11 @@ endif()
 
 if(PKG_MSCG OR PKG_ATC OR PKG_AWPMD OR PKG_ML-QUIP OR PKG_LATTE OR PKG_ELECTRODE)
   enable_language(C)
-  find_package(LAPACK)
-  find_package(BLAS)
-  if(NOT LAPACK_FOUND OR NOT BLAS_FOUND)
+  if (NOT USE_INTERNAL_LINALG)
+    find_package(LAPACK)
+    find_package(BLAS)
+  endif()
+  if(NOT LAPACK_FOUND OR NOT BLAS_FOUND OR USE_INTERNAL_LINALG)
     include(CheckGeneratorSupport)
     if(NOT CMAKE_GENERATOR_SUPPORT_FORTRAN)
       status(FATAL_ERROR "Cannot build internal linear algebra library as CMake build tool lacks Fortran support")
@@ -633,7 +650,7 @@ endif()
 # packages which selectively include variants based on enabled styles
 # e.g. accelerator packages
 ######################################################################
-foreach(PKG_WITH_INCL CORESHELL DPD-SMOOTH PHONON QEQ OPENMP KOKKOS OPT INTEL GPU)
+foreach(PKG_WITH_INCL CORESHELL DPD-SMOOTH MISC PHONON QEQ OPENMP KOKKOS OPT INTEL GPU)
   if(PKG_${PKG_WITH_INCL})
     include(Packages/${PKG_WITH_INCL})
   endif()
@@ -780,9 +797,13 @@ if(BUILD_SHARED_LIBS)
     set(Python_ADDITIONAL_VERSIONS 3.12 3.11 3.10 3.9 3.8 3.7 3.6)
     find_package(PythonInterp) # Deprecated since version 3.12
     if(PYTHONINTERP_FOUND)
-        set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
+      set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
     endif()
   else()
+    # backward compatibility
+    if(PYTHON_EXECUTABLE)
+      set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
+    endif()
     find_package(Python COMPONENTS Interpreter)
   endif()
   if(BUILD_IS_MULTI_CONFIG)
@@ -815,11 +836,17 @@ endif()
 ###############################################################################
 if(BUILD_SHARED_LIBS OR PKG_PYTHON)
   if(CMAKE_VERSION VERSION_LESS 3.12)
+    # adjust so we find Python 3 versions before Python 2 on old systems with old CMake
+    set(Python_ADDITIONAL_VERSIONS 3.12 3.11 3.10 3.9 3.8 3.7 3.6)
     find_package(PythonInterp) # Deprecated since version 3.12
     if(PYTHONINTERP_FOUND)
-        set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
+      set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
     endif()
   else()
+    # backward compatibility
+    if(PYTHON_EXECUTABLE)
+      set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
+    endif()
     find_package(Python COMPONENTS Interpreter)
   endif()
   if(Python_EXECUTABLE)

cmake/Modules/LAMMPSInterfacePlugin.cmake

@@ -6,6 +6,9 @@ if(${CMAKE_SOURCE_DIR} STREQUAL ${CMAKE_BINARY_DIR})
     "Please remove CMakeCache.txt and CMakeFiles first.")
 endif()
 
+set(LAMMPS_THIRDPARTY_URL "https://download.lammps.org/thirdparty"
+  CACHE STRING "URL for thirdparty package downloads")
+
 # global LAMMPS/plugin build settings
 set(LAMMPS_SOURCE_DIR ""  CACHE PATH "Location of LAMMPS sources folder")
 if(NOT LAMMPS_SOURCE_DIR)
@@ -78,6 +81,13 @@ function(get_newest_file path variable)
   set(${variable} ${_bestfile} PARENT_SCOPE)
 endfunction()
 
+# get LAMMPS version date
+function(get_lammps_version version_header variable)
+    file(STRINGS ${version_header} line REGEX LAMMPS_VERSION)
+    string(REGEX REPLACE "#define LAMMPS_VERSION \"([0-9]+) ([A-Za-z]+) ([0-9]+)\"" "\\1\\2\\3" date "${line}")
+    set(${variable} "${date}" PARENT_SCOPE)
+endfunction()
+
 #################################################################################
 # LAMMPS C++ interface. We only need the header related parts except on windows.
 add_library(lammps INTERFACE)
@@ -89,6 +99,7 @@ endif()
 ################################################################################
 # MPI configuration
 if(NOT CMAKE_CROSSCOMPILING)
+  set(MPI_CXX_SKIP_MPICXX TRUE)
   find_package(MPI QUIET)
   option(BUILD_MPI "Build MPI version" ${MPI_FOUND})
 else()

cmake/Modules/Packages/GPU.cmake

@@ -233,7 +233,8 @@ elseif(GPU_API STREQUAL "OPENCL")
 elseif(GPU_API STREQUAL "HIP")
   if(NOT DEFINED HIP_PATH)
       if(NOT DEFINED ENV{HIP_PATH})
-          set(HIP_PATH "/opt/rocm/hip" CACHE PATH "Path to HIP installation")
+          message(FATAL_ERROR "GPU_API=HIP requires HIP_PATH to be defined.\n"
+          "Either pass the HIP_PATH as a CMake option via -DHIP_PATH=... or set the HIP_PATH environment variable.")
       else()
           set(HIP_PATH $ENV{HIP_PATH} CACHE PATH "Path to HIP installation")
       endif()
@@ -261,6 +262,8 @@ elseif(GPU_API STREQUAL "HIP")
 
   if(HIP_PLATFORM STREQUAL "hcc" OR HIP_PLATFORM STREQUAL "amd")
     set(HIP_ARCH "gfx906" CACHE STRING "HIP target architecture")
+  elseif(HIP_PLATFORM STREQUAL "spirv")
+    set(HIP_ARCH "spirv" CACHE STRING "HIP target architecture")
   elseif(HIP_PLATFORM STREQUAL "nvcc")
     find_package(CUDA REQUIRED)
     set(HIP_ARCH "sm_50" CACHE STRING "HIP primary CUDA architecture (e.g. sm_60)")
@@ -340,7 +343,14 @@ elseif(GPU_API STREQUAL "HIP")
           VERBATIM COMMAND ${HIP_HIPCC_EXECUTABLE} --fatbin --use_fast_math -DUSE_HIP -D_${GPU_PREC_SETTING} -DLAMMPS_${LAMMPS_SIZES} ${HIP_CUDA_GENCODE} -I${LAMMPS_LIB_SOURCE_DIR}/gpu -o ${CUBIN_FILE} ${CU_FILE}
           DEPENDS ${CU_FILE}
           COMMENT "Generating ${CU_NAME}.cubin")
-    endif()
+    elseif(HIP_PLATFORM STREQUAL "spirv")
+      configure_file(${CU_FILE} ${CU_CPP_FILE} COPYONLY)
+
+      add_custom_command(OUTPUT ${CUBIN_FILE}
+        VERBATIM COMMAND ${HIP_HIPCC_EXECUTABLE} -c -O3 -DUSE_HIP -D_${GPU_PREC_SETTING} -DLAMMPS_${LAMMPS_SIZES} -I${LAMMPS_LIB_SOURCE_DIR}/gpu -o ${CUBIN_FILE} ${CU_CPP_FILE}
+        DEPENDS ${CU_CPP_FILE}
+        COMMENT "Gerating ${CU_NAME}.cubin")
+      endif()
 
     add_custom_command(OUTPUT ${CUBIN_H_FILE}
       COMMAND ${CMAKE_COMMAND} -D SOURCE_DIR=${CMAKE_CURRENT_SOURCE_DIR} -D VARNAME=${CU_NAME} -D HEADER_FILE=${CUBIN_H_FILE} -D SOURCE_FILE=${CUBIN_FILE} -P ${CMAKE_CURRENT_SOURCE_DIR}/Modules/GenerateBinaryHeader.cmake

cmake/Modules/Packages/KOKKOS.cmake

@@ -47,8 +47,8 @@ if(DOWNLOAD_KOKKOS)
   list(APPEND KOKKOS_LIB_BUILD_ARGS "-DCMAKE_CXX_EXTENSIONS=${CMAKE_CXX_EXTENSIONS}")
   list(APPEND KOKKOS_LIB_BUILD_ARGS "-DCMAKE_TOOLCHAIN_FILE=${CMAKE_TOOLCHAIN_FILE}")
   include(ExternalProject)
-  set(KOKKOS_URL "https://github.com/kokkos/kokkos/archive/3.6.00.tar.gz" CACHE STRING "URL for KOKKOS tarball")
-  set(KOKKOS_MD5 "b5c44ea961031795f434002cd7b31c20" CACHE STRING "MD5 checksum of KOKKOS tarball")
+  set(KOKKOS_URL "https://github.com/kokkos/kokkos/archive/3.6.01.tar.gz" CACHE STRING "URL for KOKKOS tarball")
+  set(KOKKOS_MD5 "0ec97fc0c356dd65bd2487defe81a7bf" CACHE STRING "MD5 checksum of KOKKOS tarball")
   mark_as_advanced(KOKKOS_URL)
   mark_as_advanced(KOKKOS_MD5)
   ExternalProject_Add(kokkos_build
@@ -72,7 +72,7 @@ if(DOWNLOAD_KOKKOS)
   add_dependencies(LAMMPS::KOKKOSCORE kokkos_build)
   add_dependencies(LAMMPS::KOKKOSCONTAINERS kokkos_build)
 elseif(EXTERNAL_KOKKOS)
-  find_package(Kokkos 3.6.00 REQUIRED CONFIG)
+  find_package(Kokkos 3.6.01 REQUIRED CONFIG)
   target_link_libraries(lammps PRIVATE Kokkos::kokkos)
   target_link_libraries(lmp PRIVATE Kokkos::kokkos)
 else()

cmake/Modules/Packages/LATTE.cmake

@@ -23,8 +23,9 @@ if(DOWNLOAD_LATTE)
   # CMake cannot pass BLAS or LAPACK library variable to external project if they are a list
   list(LENGTH BLAS_LIBRARIES} NUM_BLAS)
   list(LENGTH LAPACK_LIBRARIES NUM_LAPACK)
-  if((NUM_BLAS GREATER 1) OR (NUM_LAPACK GREATER 1))
-    message(FATAL_ERROR "Cannot compile downloaded LATTE library due to a technical limitation")
+  if((NUM_BLAS GREATER 1) OR (NUM_LAPACK GREATER 1) AND NOT USE_INTERNAL_LINALG)
+    message(FATAL_ERROR "Cannot compile downloaded LATTE library due to a technical limitation. "
+                        "Try to configure LAMMPS with '-D USE_INTERNAL_LINALG=on' added as a workaround.")
   endif()
 
   include(ExternalProject)

cmake/Modules/Packages/MDI.cmake

@@ -8,8 +8,8 @@ option(DOWNLOAD_MDI "Download and compile the MDI library instead of using an al
 
 if(DOWNLOAD_MDI)
   message(STATUS "MDI download requested - we will build our own")
-  set(MDI_URL "https://github.com/MolSSI-MDI/MDI_Library/archive/v1.3.2.tar.gz" CACHE STRING "URL for MDI tarball")
-  set(MDI_MD5 "836f5da400d8cff0f0e4435640f9454f" CACHE STRING "MD5 checksum for MDI tarball")
+  set(MDI_URL "https://github.com/MolSSI-MDI/MDI_Library/archive/v1.4.1.tar.gz" CACHE STRING "URL for MDI tarball")
+  set(MDI_MD5 "f9505fccd4c79301a619f6452dad4ad9" CACHE STRING "MD5 checksum for MDI tarball")
   mark_as_advanced(MDI_URL)
   mark_as_advanced(MDI_MD5)
   enable_language(C)
@@ -26,8 +26,21 @@ if(DOWNLOAD_MDI)
   # detect if we have python development support and thus can enable python plugins
   set(MDI_USE_PYTHON_PLUGINS OFF)
   if(CMAKE_VERSION VERSION_LESS 3.12)
+    if(NOT PYTHON_VERSION_STRING)
+      set(Python_ADDITIONAL_VERSIONS 3.12 3.11 3.10 3.9 3.8 3.7 3.6)
+      # search for interpreter first, so we have a consistent library
+      find_package(PythonInterp) # Deprecated since version 3.12
+      if(PYTHONINTERP_FOUND)
+        set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
+      endif()
+    endif()
+    # search for the library matching the selected interpreter
+    set(Python_ADDITIONAL_VERSIONS ${PYTHON_VERSION_MAJOR}.${PYTHON_VERSION_MINOR})
     find_package(PythonLibs QUIET) # Deprecated since version 3.12
     if(PYTHONLIBS_FOUND)
+      if(NOT (PYTHON_VERSION_STRING STREQUAL PYTHONLIBS_VERSION_STRING))
+        message(FATAL_ERROR "Python Library version ${PYTHONLIBS_VERSION_STRING} does not match Interpreter version ${PYTHON_VERSION_STRING}")
+      endif()
       set(MDI_USE_PYTHON_PLUGINS ON)
     endif()
   else()
@@ -44,7 +57,7 @@ if(DOWNLOAD_MDI)
   ExternalProject_Add(mdi_build
     URL     ${MDI_URL}
     URL_MD5 ${MDI_MD5}
-    CMAKE_ARGS ${CMAKE_REQUEST_PIC}
+    CMAKE_ARGS
     -DCMAKE_INSTALL_PREFIX=<INSTALL_DIR>
     -DCMAKE_C_COMPILER=${CMAKE_C_COMPILER}
     -DCMAKE_BUILD_TYPE=${CMAKE_BUILD_TYPE}
@@ -54,6 +67,7 @@ if(DOWNLOAD_MDI)
     -Dlanguage=C
     -Dlibtype=STATIC
     -Dmpi=${MDI_USE_MPI}
+    -Dplugins=ON
     -Dpython_plugins=${MDI_USE_PYTHON_PLUGINS}
     UPDATE_COMMAND ""
     INSTALL_COMMAND ""

cmake/Modules/Packages/MISC.cmake

@@ -0,0 +1,13 @@
+# pair style and fix srp/react depend on the fixes bond/break and bond/create from the MC package
+if(NOT PKG_MC)
+  get_property(LAMMPS_FIX_HEADERS GLOBAL PROPERTY FIX)
+  list(REMOVE_ITEM LAMMPS_FIX_HEADERS ${LAMMPS_SOURCE_DIR}/MISC/fix_srp_react.h)
+  set_property(GLOBAL PROPERTY FIX "${LAMMPS_FIX_HEADERS}")
+  get_property(LAMMPS_PAIR_HEADERS GLOBAL PROPERTY PAIR)
+  list(REMOVE_ITEM LAMMPS_PAIR_HEADERS ${LAMMPS_SOURCE_DIR}/MISC/pair_srp_react.h)
+  set_property(GLOBAL PROPERTY PAIR "${LAMMPS_PAIR_HEADERS}")
+  get_target_property(LAMMPS_SOURCES lammps SOURCES)
+  list(REMOVE_ITEM LAMMPS_SOURCES ${LAMMPS_SOURCE_DIR}/MISC/fix_srp_react.cpp)
+  list(REMOVE_ITEM LAMMPS_SOURCES ${LAMMPS_SOURCE_DIR}/MISC/pair_srp_react.cpp)
+  set_property(TARGET lammps PROPERTY SOURCES "${LAMMPS_SOURCES}")
+endif()

cmake/Modules/Packages/ML-HDNNP.cmake

@@ -44,7 +44,9 @@ if(DOWNLOAD_N2P2)
   else()
     # get path to MPI include directory
     get_target_property(N2P2_MPI_INCLUDE MPI::MPI_CXX INTERFACE_INCLUDE_DIRECTORIES)
-    set(N2P2_PROJECT_OPTIONS "-I${N2P2_MPI_INCLUDE}")
+    foreach (_INCL ${N2P2_MPI_INCLUDE})
+      set(N2P2_PROJECT_OPTIONS "${N2P2_PROJECT_OPTIONS} -I${_INCL}")
+    endforeach()
   endif()
 
   # prefer GNU make, if available. N2P2 lib seems to need it.
@@ -75,7 +77,7 @@ if(DOWNLOAD_N2P2)
     UPDATE_COMMAND ""
     CONFIGURE_COMMAND ""
     PATCH_COMMAND sed -i -e "s/\\(MPI_\\(P\\|Unp\\)ack(\\)/\\1(void *) /" src/libnnpif/LAMMPS/InterfaceLammps.cpp
-    BUILD_COMMAND ${N2P2_MAKE} -f makefile libnnpif ${N2P2_BUILD_OPTIONS}
+    BUILD_COMMAND ${N2P2_MAKE} -C <SOURCE_DIR>/src -f makefile libnnpif ${N2P2_BUILD_OPTIONS}
     BUILD_ALWAYS YES
     INSTALL_COMMAND ""
     BUILD_IN_SOURCE 1

cmake/Modules/Packages/PYTHON.cmake

@@ -1,8 +1,28 @@
 if(CMAKE_VERSION VERSION_LESS 3.12)
+  if(NOT PYTHON_VERSION_STRING)
+    set(Python_ADDITIONAL_VERSIONS 3.12 3.11 3.10 3.9 3.8 3.7 3.6)
+    # search for interpreter first, so we have a consistent library
+    find_package(PythonInterp) # Deprecated since version 3.12
+    if(PYTHONINTERP_FOUND)
+      set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
+    endif()
+  endif()
+  # search for the library matching the selected interpreter
+  set(Python_ADDITIONAL_VERSIONS ${PYTHON_VERSION_MAJOR}.${PYTHON_VERSION_MINOR})
   find_package(PythonLibs REQUIRED) # Deprecated since version 3.12
+  if(NOT (PYTHON_VERSION_STRING STREQUAL PYTHONLIBS_VERSION_STRING))
+    message(FATAL_ERROR "Python Library version ${PYTHONLIBS_VERSION_STRING} does not match Interpreter version ${PYTHON_VERSION_STRING}")
+  endif()
   target_include_directories(lammps PRIVATE ${PYTHON_INCLUDE_DIRS})
   target_link_libraries(lammps PRIVATE ${PYTHON_LIBRARIES})
 else()
+  if(NOT Python_INTERPRETER)
+    # backward compatibility
+    if(PYTHON_EXECUTABLE)
+      set(Python_EXECUTABLE ${PYTHON_EXECUTABLE})
+    endif()
+    find_package(Python COMPONENTS Interpreter)
+  endif()
   find_package(Python REQUIRED COMPONENTS Interpreter Development)
   target_link_libraries(lammps PRIVATE Python::Python)
 endif()

cmake/presets/all_off.cmake

@@ -3,6 +3,7 @@
 
 set(ALL_PACKAGES
   ADIOS
+  AMOEBA
   ASPHERE
   ATC
   AWPMD
@@ -11,7 +12,7 @@ set(ALL_PACKAGES
   BPM
   BROWNIAN
   CG-DNA
-  CG-SDK
+  CG-SPICA
   CLASS2
   COLLOID
   COLVARS

cmake/presets/all_on.cmake

@@ -5,6 +5,7 @@
 
 set(ALL_PACKAGES
   ADIOS
+  AMOEBA
   ASPHERE
   ATC
   AWPMD
@@ -13,7 +14,7 @@ set(ALL_PACKAGES
   BPM
   BROWNIAN
   CG-DNA
-  CG-SDK
+  CG-SPICA
   CLASS2
   COLLOID
   COLVARS

cmake/presets/mingw-cross.cmake

@@ -1,4 +1,5 @@
 set(WIN_PACKAGES
+  AMOEBA
   ASPHERE
   ATC
   AWPMD
@@ -7,7 +8,7 @@ set(WIN_PACKAGES
   BPM
   BROWNIAN
   CG-DNA
-  CG-SDK
+  CG-SPICA
   CLASS2
   COLLOID
   COLVARS

cmake/presets/most.cmake

@@ -3,13 +3,14 @@
 # are removed. The resulting binary should be able to run most inputs.
 
 set(ALL_PACKAGES
+  AMOEBA
   ASPHERE
   BOCS
   BODY
   BPM
   BROWNIAN
   CG-DNA
-  CG-SDK
+  CG-SPICA
   CLASS2
   COLLOID
   COLVARS

cmake/presets/windows.cmake

@@ -1,11 +1,12 @@
 set(WIN_PACKAGES
+  AMOEBA
   ASPHERE
   BOCS
   BODY
   BPM
   BROWNIAN
   CG-DNA
-  CG-SDK
+  CG-SPICA
   CLASS2
   COLLOID
   COLVARS

doc/Makefile

@@ -242,7 +242,6 @@ $(MATHJAX):
 $(ANCHORCHECK): $(VENV)
 	@( \
 		. $(VENV)/bin/activate; \
-		(cd utils/converters;\
-		python setup.py develop);\
+		pip $(PIP_OPTIONS) install -e utils/converters;\
 		deactivate;\
 	)

doc/lammps.1

@@ -1,7 +1,7 @@
-.TH LAMMPS "1" "23 June 2022" "2022-6-23"
+.TH LAMMPS "1" "3 August 2022" "2022-8-3"
 .SH NAME
 .B LAMMPS
-\- Molecular Dynamics Simulator.  Version 23 June 2022
+\- Molecular Dynamics Simulator.  Version 3 August 2022
 
 .SH SYNOPSIS
 .B lmp
@@ -161,7 +161,7 @@ list references for specific cite-able features used during a
 run.
 .TP
 \fB\-pk <style> [options]\fR or \fB\-package <style> [options]\fR
-Invoke the \fBpackage\R command with <style> and optional arguments.
+Invoke the \fBpackage\fR command with <style> and optional arguments.
 The syntax is the same as if the command appeared in an input script.
 For example "-pk gpu 2" is the same as "package gpu 2" in the input
 script. The possible styles and options are discussed in the

doc/src/Build_extras.rst

@@ -123,6 +123,7 @@ CMake build
    -D GPU_API=value             # value = opencl (default) or cuda or hip
    -D GPU_PREC=value            # precision setting
                                 # value = double or mixed (default) or single
+   -D HIP_PATH                  # path to HIP installation. Must be set if GPU_API=HIP
    -D GPU_ARCH=value            # primary GPU hardware choice for GPU_API=cuda
                                 # value = sm_XX, see below
                                 # default is sm_50
@@ -179,10 +180,17 @@ set appropriate environment variables. Some variables such as
 :code:`HCC_AMDGPU_TARGET` (for ROCm <= 4.0) or :code:`CUDA_PATH` are necessary for :code:`hipcc`
 and the linker to work correctly.
 
+Using CHIP-SPV implementation of HIP is now supported. It allows one to run HIP
+code on Intel GPUs via the OpenCL or Level Zero backends. To use CHIP-SPV, you must
+set :code:`-DHIP_USE_DEVICE_SORT=OFF` in your CMake command line as CHIP-SPV does not
+yet support hipCUB. The use of HIP for Intel GPUs is still experimental so you
+should only use this option in preparations to run on Aurora system at ANL.
+
 .. code:: bash
 
    # AMDGPU target (ROCm <= 4.0)
    export HIP_PLATFORM=hcc
+   export HIP_PATH=/path/to/HIP/install
    export HCC_AMDGPU_TARGET=gfx906
    cmake -D PKG_GPU=on -D GPU_API=HIP -D HIP_ARCH=gfx906 -D CMAKE_CXX_COMPILER=hipcc ..
    make -j 4
@@ -191,6 +199,7 @@ and the linker to work correctly.
 
    # AMDGPU target (ROCm >= 4.1)
    export HIP_PLATFORM=amd
+   export HIP_PATH=/path/to/HIP/install
    cmake -D PKG_GPU=on -D GPU_API=HIP -D HIP_ARCH=gfx906 -D CMAKE_CXX_COMPILER=hipcc ..
    make -j 4
 
@@ -199,10 +208,20 @@ and the linker to work correctly.
    # CUDA target (not recommended, use GPU_ARCH=cuda)
    # !!! DO NOT set CMAKE_CXX_COMPILER !!!
    export HIP_PLATFORM=nvcc
+   export HIP_PATH=/path/to/HIP/install
    export CUDA_PATH=/usr/local/cuda
    cmake -D PKG_GPU=on -D GPU_API=HIP -D HIP_ARCH=sm_70 ..
    make -j 4
 
+.. code:: bash
+
+   # SPIR-V target (Intel GPUs)
+   export HIP_PLATFORM=spirv
+   export HIP_PATH=/path/to/HIP/install
+   export CMAKE_CXX_COMPILER=<hipcc/clang++>
+   cmake -D PKG_GPU=on -D GPU_API=HIP ..
+   make -j 4
+
 Traditional make
 ^^^^^^^^^^^^^^^^
 
@@ -788,8 +807,10 @@ library.
 
       .. code-block:: bash
 
-         -D DOWNLOAD_LATTE=value    # download LATTE for build, value = no (default) or yes
-         -D LATTE_LIBRARY=path      # LATTE library file (only needed if a custom location)
+         -D DOWNLOAD_LATTE=value      # download LATTE for build, value = no (default) or yes
+         -D LATTE_LIBRARY=path        # LATTE library file (only needed if a custom location)
+         -D USE_INTERNAL_LINALG=value # Use the internal linear algebra library instead of LAPACK
+                                      #   value = no (default) or yes
 
       If ``DOWNLOAD_LATTE`` is set, the LATTE library will be downloaded
       and built inside the CMake build directory.  If the LATTE library
@@ -797,6 +818,13 @@ library.
       ``LATTE_LIBRARY`` is the filename (plus path) of the LATTE library
       file, not the directory the library file is in.
 
+      The LATTE library requires LAPACK (and BLAS) and CMake can identify
+      their locations and pass that info to the LATTE build script. But
+      on some systems this triggers a (current) limitation of CMake and
+      the configuration will fail. Try enabling ``USE_INTERNAL_LINALG`` in
+      those cases to use the bundled linear algebra library and work around
+      the limitation.
+
    .. tab:: Traditional make
 
       You can download and build the LATTE library manually if you
@@ -1913,14 +1941,25 @@ within CMake will download the non-commercial use version.
 
       .. code-block:: bash
 
-         -D DOWNLOAD_QUIP=value   # download OpenKIM API v2 for build, value = no (default) or yes
-         -D QUIP_LIBRARY=path     # path to libquip.a (only needed if a custom location)
+         -D DOWNLOAD_QUIP=value       # download QUIP library for build, value = no (default) or yes
+         -D QUIP_LIBRARY=path         # path to libquip.a (only needed if a custom location)
+         -D USE_INTERNAL_LINALG=value # Use the internal linear algebra library instead of LAPACK
+                                      #   value = no (default) or yes
 
-      CMake will try to download and build the QUIP library from GitHub, if it is not
-      found on the local machine. This requires to have git installed. It will use the same compilers
-      and flags as used for compiling LAMMPS.  Currently this is only supported for the GNU and the
-      Intel compilers. Set the ``QUIP_LIBRARY`` variable if you want to use a previously compiled
-      and installed QUIP library and CMake cannot find it.
+      CMake will try to download and build the QUIP library from GitHub,
+      if it is not found on the local machine. This requires to have git
+      installed. It will use the same compilers and flags as used for
+      compiling LAMMPS.  Currently this is only supported for the GNU
+      and the Intel compilers. Set the ``QUIP_LIBRARY`` variable if you
+      want to use a previously compiled and installed QUIP library and
+      CMake cannot find it.
+
+      The QUIP library requires LAPACK (and BLAS) and CMake can identify
+      their locations and pass that info to the QUIP build script. But
+      on some systems this triggers a (current) limitation of CMake and
+      the configuration will fail. Try enabling ``USE_INTERNAL_LINALG`` in
+      those cases to use the bundled linear algebra library and work around
+      the limitation.
 
    .. tab:: Traditional make
 

doc/src/Build_manual.rst

@@ -48,18 +48,15 @@ Build using GNU make
 
 The LAMMPS manual is written in `reStructuredText <rst_>`_ format which
 can be translated to different output format using the `Sphinx
-<sphinx_>`_ document generator tool.  It also incorporates programmer
-documentation extracted from the LAMMPS C++ sources through the `Doxygen
-<https://doxygen.nl>`_ program.  Currently the translation to HTML, PDF
-(via LaTeX), ePUB (for many e-book readers) and MOBI (for Amazon Kindle
-readers) are supported.  For that to work a Python 3 interpreter, the
-``doxygen`` tools and internet access to download additional files and
-tools are required.  This download is usually only required once or
-after the documentation folder is returned to a pristine state with
-``make clean-all``.
-
-.. _rst: https://docutils.readthedocs.io/en/sphinx-docs/user/rst/quickstart.html
-.. _sphinx: https://www.sphinx-doc.org
+<https://sphinx-doc.org>`_ document generator tool.  It also
+incorporates programmer documentation extracted from the LAMMPS C++
+sources through the `Doxygen <https://doxygen.nl>`_ program.  Currently
+the translation to HTML, PDF (via LaTeX), ePUB (for many e-book readers)
+and MOBI (for Amazon Kindle readers) are supported.  For that to work a
+Python 3 interpreter, the ``doxygen`` tools and internet access to
+download additional files and tools are required.  This download is
+usually only required once or after the documentation folder is returned
+to a pristine state with ``make clean-all``.
 
 For the documentation build a python virtual environment is set up in
 the folder ``doc/docenv`` and various python packages are installed into
@@ -252,6 +249,5 @@ manual with ``make spelling``.  This requires `a library called enchant
 positives* (e.g. keywords, names, abbreviations) those can be added to
 the file ``lammps/doc/utils/sphinx-config/false_positives.txt``.
 
-.. _rst: https://docutils.readthedocs.io/en/sphinx-docs/user/rst/quickstart.html
-
 .. _lws: https://www.lammps.org
+.. _rst: https://www.sphinx-doc.org/en/master/usage/restructuredtext/index.html

doc/src/Commands.rst

@@ -21,6 +21,7 @@ commands in it are used to define a LAMMPS simulation.
    Commands_pair
    Commands_bond
    Commands_kspace
+   Commands_dump
 
 .. toctree::
    :maxdepth: 1

doc/src/Commands_all.rst

@@ -10,6 +10,7 @@
    * :ref:`Dihedral styles <dihedral>`
    * :ref:`Improper styles <improper>`
    * :doc:`KSpace styles <Commands_kspace>`
+   * :doc:`Dump styles <Commands_dump>`
 
 General commands
 ================

doc/src/Commands_bond.rst

@@ -10,6 +10,7 @@
    * :ref:`Dihedral styles <dihedral>`
    * :ref:`Improper styles <improper>`
    * :doc:`KSpace styles <Commands_kspace>`
+   * :doc:`Dump styles <Commands_dump>`
 
 .. _bond:
 
@@ -74,6 +75,7 @@ OPT.
    *
    *
    *
+   * :doc:`amoeba <angle_amoeba>`
    * :doc:`charmm (iko) <angle_charmm>`
    * :doc:`class2 (ko) <angle_class2>`
    * :doc:`class2/p6 <angle_class2>`
@@ -92,7 +94,7 @@ OPT.
    * :doc:`harmonic (iko) <angle_harmonic>`
    * :doc:`mm3 <angle_mm3>`
    * :doc:`quartic (o) <angle_quartic>`
-   * :doc:`sdk (o) <angle_sdk>`
+   * :doc:`spica (o) <angle_spica>`
    * :doc:`table (o) <angle_table>`
 
 .. _dihedral:
@@ -152,6 +154,7 @@ OPT.
    *
    *
    *
+   * :doc:`amoeba <improper_amoeba>`
    * :doc:`class2 (ko) <improper_class2>`
    * :doc:`cossq (o) <improper_cossq>`
    * :doc:`cvff (io) <improper_cvff>`

doc/src/Commands_compute.rst

@@ -10,6 +10,7 @@
    * :ref:`Dihedral styles <dihedral>`
    * :ref:`Improper styles <improper>`
    * :doc:`KSpace styles <Commands_kspace>`
+   * :doc:`Dump styles <Commands_dump>`
 
 Compute commands
 ================
@@ -138,6 +139,8 @@ KOKKOS, o = OPENMP, t = OPT.
    * :doc:`smd/vol <compute_smd_vol>`
    * :doc:`snap <compute_sna_atom>`
    * :doc:`sna/atom <compute_sna_atom>`
+   * :doc:`sna/grid <compute_sna_atom>`
+   * :doc:`sna/grid/local <compute_sna_atom>`
    * :doc:`snad/atom <compute_sna_atom>`
    * :doc:`snav/atom <compute_sna_atom>`
    * :doc:`sph/e/atom <compute_sph_e_atom>`

doc/src/Commands_dump.rst

@@ -0,0 +1,56 @@
+.. table_from_list::
+   :columns: 3
+
+   * :doc:`General commands <Commands_all>`
+   * :doc:`Fix styles <Commands_fix>`
+   * :doc:`Compute styles <Commands_compute>`
+   * :doc:`Pair styles <Commands_pair>`
+   * :ref:`Bond styles <bond>`
+   * :ref:`Angle styles <angle>`
+   * :ref:`Dihedral styles <dihedral>`
+   * :ref:`Improper styles <improper>`
+   * :doc:`KSpace styles <Commands_kspace>`
+   * :doc:`Dump styles <Commands_dump>`
+
+Dump commands
+=============
+
+An alphabetic list of all LAMMPS :doc:`dump <dump>` commands.
+
+.. table_from_list::
+   :columns: 5
+
+   * :doc:`atom <dump>`
+   * :doc:`atom/adios <dump_adios>`
+   * :doc:`atom/gz <dump>`
+   * :doc:`atom/mpiio <dump>`
+   * :doc:`atom/zstd <dump>`
+   * :doc:`cfg <dump>`
+   * :doc:`cfg/gz <dump>`
+   * :doc:`cfg/mpiio <dump>`
+   * :doc:`cfg/uef <dump_cfg_uef>`
+   * :doc:`cfg/zstd <dump>`
+   * :doc:`custom <dump>`
+   * :doc:`custom/adios <dump_adios>`
+   * :doc:`custom/gz <dump>`
+   * :doc:`custom/mpiio <dump>`
+   * :doc:`custom/zstd <dump>`
+   * :doc:`dcd <dump>`
+   * :doc:`deprecated <dump>`
+   * :doc:`h5md <dump_h5md>`
+   * :doc:`image <dump_image>`
+   * :doc:`local <dump>`
+   * :doc:`local/gz <dump>`
+   * :doc:`local/zstd <dump>`
+   * :doc:`molfile <dump_molfile>`
+   * :doc:`movie <dump_image>`
+   * :doc:`netcdf <dump_netcdf>`
+   * :doc:`netcdf/mpiio <dump>`
+   * :doc:`vtk <dump_vtk>`
+   * :doc:`xtc <dump>`
+   * :doc:`xyz <dump>`
+   * :doc:`xyz/gz <dump>`
+   * :doc:`xyz/mpiio <dump>`
+   * :doc:`xyz/zstd <dump>`
+   * :doc:`yaml <dump>`
+

doc/src/Commands_fix.rst

@@ -10,6 +10,7 @@
    * :ref:`Dihedral styles <dihedral>`
    * :ref:`Improper styles <improper>`
    * :doc:`KSpace styles <Commands_kspace>`
+   * :doc:`Dump styles <Commands_dump>`
 
 Fix commands
 ============
@@ -28,6 +29,8 @@ OPT.
    * :doc:`adapt/fep <fix_adapt_fep>`
    * :doc:`addforce <fix_addforce>`
    * :doc:`addtorque <fix_addtorque>`
+   * :doc:`amoeba/bitorsion <fix_amoeba_bitorsion>`
+   * :doc:`amoeba/pitorsion <fix_amoeba_pitorsion>`
    * :doc:`append/atoms <fix_append_atoms>`
    * :doc:`atc <fix_atc>`
    * :doc:`atom/swap <fix_atom_swap>`
@@ -103,7 +106,7 @@ OPT.
    * :doc:`lb/viscous <fix_lb_viscous>`
    * :doc:`lineforce <fix_lineforce>`
    * :doc:`manifoldforce <fix_manifoldforce>`
-   * :doc:`mdi/aimd <fix_mdi_aimd>`
+   * :doc:`mdi/qm <fix_mdi_qm>`
    * :doc:`meso/move <fix_meso_move>`
    * :doc:`mol/swap <fix_mol_swap>`
    * :doc:`momentum (k) <fix_momentum>`

doc/src/Commands_kspace.rst

@@ -10,6 +10,7 @@
    * :ref:`Dihedral styles <dihedral>`
    * :ref:`Improper styles <improper>`
    * :doc:`KSpace styles <Commands_kspace>`
+   * :doc:`Dump styles <Commands_dump>`
 
 KSpace solvers
 ==============

doc/src/Commands_pair.rst

@@ -10,6 +10,7 @@
    * :ref:`Dihedral styles <dihedral>`
    * :ref:`Improper styles <improper>`
    * :doc:`KSpace styles <Commands_kspace>`
+   * :doc:`Dump styles <Commands_dump>`
 
 Pair_style potentials
 ======================
@@ -38,6 +39,7 @@ OPT.
    * :doc:`agni (o) <pair_agni>`
    * :doc:`airebo (io) <pair_airebo>`
    * :doc:`airebo/morse (io) <pair_airebo>`
+   * :doc:`amoeba <pair_amoeba>`
    * :doc:`atm <pair_atm>`
    * :doc:`awpmd/cut <pair_awpmd>`
    * :doc:`beck (go) <pair_beck>`
@@ -124,6 +126,7 @@ OPT.
    * :doc:`hbond/dreiding/lj (o) <pair_hbond_dreiding>`
    * :doc:`hbond/dreiding/morse (o) <pair_hbond_dreiding>`
    * :doc:`hdnnp <pair_hdnnp>`
+   * :doc:`hippo <pair_amoeba>`
    * :doc:`ilp/graphene/hbn (t) <pair_ilp_graphene_hbn>`
    * :doc:`ilp/tmd (t) <pair_ilp_tmd>`
    * :doc:`kolmogorov/crespi/full <pair_kolmogorov_crespi_full>`
@@ -179,9 +182,9 @@ OPT.
    * :doc:`lj/long/tip4p/long (o) <pair_lj_long>`
    * :doc:`lj/mdf <pair_mdf>`
    * :doc:`lj/relres (o) <pair_lj_relres>`
-   * :doc:`lj/sdk (gko) <pair_sdk>`
-   * :doc:`lj/sdk/coul/long (go) <pair_sdk>`
-   * :doc:`lj/sdk/coul/msm (o) <pair_sdk>`
+   * :doc:`lj/spica (gko) <pair_spica>`
+   * :doc:`lj/spica/coul/long (go) <pair_spica>`
+   * :doc:`lj/spica/coul/msm (o) <pair_spica>`
    * :doc:`lj/sf/dipole/sf (go) <pair_dipole>`
    * :doc:`lj/smooth (go) <pair_lj_smooth>`
    * :doc:`lj/smooth/linear (o) <pair_lj_smooth_linear>`
@@ -194,7 +197,7 @@ OPT.
    * :doc:`lubricateU/poly <pair_lubricateU>`
    * :doc:`mdpd <pair_mesodpd>`
    * :doc:`mdpd/rhosum <pair_mesodpd>`
-   * :doc:`meam <pair_meam>`
+   * :doc:`meam (k) <pair_meam>`
    * :doc:`meam/spline (o) <pair_meam_spline>`
    * :doc:`meam/sw/spline <pair_meam_sw_spline>`
    * :doc:`mesocnt <pair_mesocnt>`
@@ -268,6 +271,7 @@ OPT.
    * :doc:`spin/magelec <pair_spin_magelec>`
    * :doc:`spin/neel <pair_spin_neel>`
    * :doc:`srp <pair_srp>`
+   * :doc:`srp/react <pair_srp>`
    * :doc:`sw (giko) <pair_sw>`
    * :doc:`sw/angle/table <pair_sw_angle_table>`
    * :doc:`sw/mod (o) <pair_sw>`

doc/src/Commands_parse.rst

@@ -123,14 +123,15 @@ LAMMPS:
 .. _six:
 
 6. If you want text with spaces to be treated as a single argument, it
-   can be enclosed in either single or double or triple quotes.  A long
-   single argument enclosed in single or double quotes can span multiple
-   lines if the "&" character is used, as described above.  When the
-   lines are concatenated together (and the "&" characters and line
-   breaks removed), the text will become a single line.  If you want
-   multiple lines of an argument to retain their line breaks, the text
-   can be enclosed in triple quotes, in which case "&" characters are
-   not needed.  For example:
+   can be enclosed in either single (') or double (") or triple (""")
+   quotes.  A long single argument enclosed in single or double quotes
+   can span multiple lines if the "&" character is used, as described
+   in :ref:`1 <one>` above.  When the lines are concatenated together
+   by LAMMPS (and the "&" characters and line breaks removed), the
+   combined text will become a single line.  If you want multiple lines
+   of an argument to retain their line breaks, the text can be enclosed
+   in triple quotes, in which case "&" characters are not needed and do
+   not function as line continuation character.  For example:
 
    .. code-block:: LAMMPS
 
@@ -144,8 +145,9 @@ LAMMPS:
       System temperature = $t
       """
 
-   In each case, the single, double, or triple quotes are removed when
-   the single argument they enclose is stored internally.
+   In each of these cases, the single, double, or triple quotes are
+   removed and the enclosed text stored internally as a single
+   argument.
 
    See the :doc:`dump modify format <dump_modify>`, :doc:`print
    <print>`, :doc:`if <if>`, and :doc:`python <python>` commands for

doc/src/Developer.rst

@@ -17,6 +17,7 @@ of time and requests from the LAMMPS user community.
    Developer_flow
    Developer_write
    Developer_notes
+   Developer_updating
    Developer_plugins
    Developer_unittest
    Classes

doc/src/Developer_updating.rst

@@ -0,0 +1,324 @@
+Notes for updating code written for older LAMMPS versions
+---------------------------------------------------------
+
+This section documents how C++ source files that are available *outside
+of the LAMMPS source distribution* (e.g. in external USER packages or as
+source files provided as a supplement to a publication) that are written
+for an older version of LAMMPS and thus need to be updated to be
+compatible with the current version of LAMMPS.  Due to the active
+development of LAMMPS it is likely to always be incomplete.  Please
+contact developer@lammps.org in case you run across an issue that is not
+(yet) listed here.  Please also review the latest information about the
+LAMMPS :doc:`programming style conventions <Modify_style>`, especially
+if you are considering to submit the updated version for inclusion into
+the LAMMPS distribution.
+
+Available topics in mostly chronological order are:
+
+- `Setting flags in the constructor`_
+- `Rename of pack/unpack_comm() to pack/unpack_forward_comm()`_
+- `Use ev_init() to initialize variables derived from eflag and vflag`_
+- `Use utils::numeric() functions instead of force->numeric()`_
+- `Use utils::open_potential() function to open potential files`_
+- `Simplify customized error messages`_
+- `Use of "override" instead of "virtual"`_
+- `Simplified and more compact neighbor list requests`_
+
+----
+
+Setting flags in the constructor
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+As LAMMPS gains additional functionality, new flags may need to be set
+in the constructor or a class to signal compatibility with such features.
+Most of the time the defaults are chosen conservatively, but sometimes
+the conservative choice is the uncommon choice, and then those settings
+need to be made when updating code.
+
+Pair styles:
+
+  - ``manybody_flag``: set to 1 if your pair style is not pair-wise additive
+  - ``restartinfo``: set to 0 if your pair style does not store data in restart files
+
+
+Rename of pack/unpack_comm() to pack/unpack_forward_comm()
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 8Aug2014
+
+In this change set the functions to pack data into communication buffers
+and to unpack data from communication buffers for :doc:`forward
+communications <Developer_comm_ops>` were renamed from ``pack_comm()``
+and ``unpack_comm()`` to ``pack_forward_comm()`` and
+``unpack_forward_comm()``, respectively.  Also the meaning of the return
+value of these functions was changed: rather than returning the number
+of items per atom stored in the buffer, now the total number of items
+added (or unpacked) needs to be returned.  Here is an example from the
+`PairEAM` class.  Of course the member function declaration in corresponding
+header file needs to be updated accordingly.
+
+Old:
+
+.. code-block:: C++
+
+   int PairEAM::pack_comm(int n, int *list, double *buf, int pbc_flag, int *pbc)
+   {
+     int m = 0;
+     for (int i = 0; i < n; i++) {
+       int j = list[i];
+       buf[m++] = fp[j];
+     }
+     return 1;
+   }
+
+New:
+
+.. code-block:: C++
+
+   int PairEAM::pack_forward_comm(int n, int *list, double *buf, int pbc_flag, int *pbc)
+   {
+     int m = 0;
+     for (int i = 0; i < n; i++) {
+       int j = list[i];
+       buf[m++] = fp[j];
+     }
+     return m;
+   }
+
+.. note::
+
+   Because the various "pack" and "unpack" functions are defined in the
+   respective base classes as dummy functions doing nothing, and because
+   of the the name mismatch the custom versions in the derived class
+   will no longer be called, there will be no compilation error when
+   this change is not applied.  Only calculations will suddenly produce
+   incorrect results because the required forward communication calls
+   will cease to function correctly.
+
+Use ev_init() to initialize variables derived from eflag and vflag
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 29Mar2019
+
+There are several variables that need to be initialized based on
+the values of the "eflag" and "vflag" variables and since sometimes
+there are new bits added and new variables need to be set to 1 or 0.
+To make this consistent, across all styles, there is now an inline
+function ``ev_init(eflag, vflag)`` that makes those settings
+consistently and calls either ``ev_setup()`` or ``ev_unset()``.
+Example from a pair style:
+
+Old:
+
+.. code-block:: C++
+
+   if (eflag || vflag) ev_setup(eflag, vflag);
+   else evflag = vflag_fdotr = eflag_global = eflag_atom = 0;
+
+New:
+
+.. code-block:: C++
+
+   ev_init(eflag, vflag);
+
+Not applying this change will not cause a compilation error, but
+can lead to inconsistent behavior and incorrect tallying of
+energy or virial.
+
+Use utils::numeric() functions instead of force->numeric()
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 18Sep2020
+
+The "numeric()" conversion functions (including "inumeric()",
+"bnumeric()", and "tnumeric()") have been moved from the Force class to
+the utils namespace.  Also they take an additional argument that selects
+whether the ``Error::all()`` or ``Error::one()`` function should be
+called in case of an error.  The former should be used when *all* MPI
+processes call the conversion function and the latter *must* be used
+when they are called from only one or a subset of the MPI processes.
+
+Old:
+
+.. code-block:: C++
+
+    val = force->numeric(FLERR, arg[1]);
+    num = force->inumeric(FLERR, arg[2]);
+
+New:
+
+.. code-block:: C++
+
+    val = utils::numeric(FLERR, true, arg[1], lmp);
+    num = utils::inumeric(FLERR, false, arg[2], lmp);
+
+.. seealso::
+
+   :cpp:func:`utils::numeric() <LAMMPS_NS::utils::numeric>`,
+   :cpp:func:`utils::inumeric() <LAMMPS_NS::utils::inumeric>`,
+   :cpp:func:`utils::bnumeric() <LAMMPS_NS::utils::bnumeric>`,
+   :cpp:func:`utils::tnumeric() <LAMMPS_NS::utils::tnumeric>`
+
+Use utils::open_potential() function to open potential files
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 18Sep2020
+
+The :cpp:func:`utils::open_potential()
+<LAMMPS_NS::utils::open_potential>` function must be used to replace
+calls to ``force->open_potential()`` and should be used to replace
+``fopen()`` for opening potential files for reading.  The custom
+function does three additional steps compared to ``fopen()``: 1) it will
+try to parse the ``UNITS:`` and ``DATE:`` metadata will stop with an
+error on a units mismatch and will print the date info, if present, in
+the log file; 2) for pair styles that support it, it will set up
+possible automatic unit conversions based on the embedded unit
+information and LAMMPS' current units setting; 3) it will not only try
+to open a potential file at the given path, but will also search in the
+folders listed in the ``LAMMPS_POTENTIALS`` environment variable.  This
+allows to keep potential files in a common location instead of having to
+copy them around for simulations.
+
+Old:
+
+.. code-block:: C++
+
+   fp = force->open_potential(filename);
+   fp = fopen(filename, "r");
+
+New:
+
+.. code-block:: C++
+
+   fp = utils::open_potential(filename, lmp);
+
+Simplify customized error messages
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 14May2021
+
+Aided by features of the bundled {fmt} library, error messages now
+can have a variable number of arguments and the string will be interpreted
+as a {fmt} style format string so that custom error messages can be
+easily customized without having to use temporary buffers and ``sprintf()``.
+Example:
+
+Old:
+
+.. code-block:: C++
+
+   if (fptr == NULL) {
+     char str[128];
+     sprintf(str,"Cannot open AEAM potential file %s",filename);
+     error->one(FLERR,str);
+   }
+
+New:
+
+.. code-block:: C++
+
+   if (fptr == nullptr)
+     error->one(FLERR, "Cannot open AEAM potential file {}: {}", filename, utils::getsyserror());
+
+Use of "override" instead of "virtual"
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 17Feb2022
+
+Since LAMMPS requires C++11 we switched to use the "override" keyword
+instead of "virtual" to indicate polymorphism in derived classes.  This
+allows the C++ compiler to better detect inconsistencies when an
+override is intended or not.  Please note that "override" has to be
+added to **all** polymorph functions in derived classes and "virtual"
+*only* to the function in the base class (or the destructor).  Here is
+an example from the ``FixWallReflect`` class:
+
+Old:
+
+.. code-block:: C++
+
+   FixWallReflect(class LAMMPS *, int, char **);
+   virtual ~FixWallReflect();
+   int setmask();
+   void init();
+   void post_integrate();
+
+New:
+
+.. code-block:: C++
+
+   FixWallReflect(class LAMMPS *, int, char **);
+   ~FixWallReflect() override;
+   int setmask() override;
+   void init() override;
+   void post_integrate() override;
+
+This change set will neither cause a compilation failure, nor will it
+change functionality, but if you plan to submit the updated code for
+inclusion into the LAMMPS distribution, it will be requested for achieve
+a consistent :doc:`programming style <Modify_style>`.
+
+Simplified function names for forward and reverse communication
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 24Mar2022
+
+Rather then using the function name to distinguish between the different
+forward and reverse communication functions for styles, LAMMPS now uses
+the type of the "this" pointer argument.
+
+Old:
+
+.. code-block:: C++
+
+   comm->forward_comm_pair(this);
+   comm->forward_comm_fix(this);
+   comm->forward_comm_compute(this);
+   comm->forward_comm_dump(this);
+   comm->reverse_comm_pair(this);
+   comm->reverse_comm_fix(this);
+   comm->reverse_comm_compute(this);
+   comm->reverse_comm_dump(this);
+
+New:
+
+.. code-block:: C++
+
+   comm->forward_comm(this);
+   comm->reverse_comm(this);
+
+This change is required or else the code will not compile.
+
+Simplified and more compact neighbor list requests
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+.. versionchanged:: 24Mar2022
+
+This change set reduces the amount of code required to request a
+neighbor list.  It enforces consistency and no longer requires to change
+internal data of the request.  More information on neighbor list
+requests can be :doc:`found here <Developer_notes>`. Example from the
+``ComputeRDF`` class:
+
+Old:
+
+.. code-block:: C++
+
+   int irequest = neighbor->request(this,instance_me);
+   neighbor->requests[irequest]->pair = 0;
+   neighbor->requests[irequest]->compute = 1;
+   neighbor->requests[irequest]->occasional = 1;
+   if (cutflag) {
+     neighbor->requests[irequest]->cut = 1;
+     neighbor->requests[irequest]->cutoff = mycutneigh;
+   }
+
+New:
+
+.. code-block:: C++
+
+   auto req = neighbor->add_request(this, NeighConst::REQ_OCCASIONAL);
+   if (cutflag) req->set_cutoff(mycutneigh);
+
+Public access to the ``NeighRequest`` class data members has been
+removed so this update is *required* to avoid compilation failure.

doc/src/Developer_utils.rst

@@ -154,6 +154,9 @@ and parsing files or arguments.
 .. doxygenfunction:: trim_and_count_words
    :project: progguide
 
+.. doxygenfunction:: join_words
+   :project: progguide
+
 .. doxygenfunction:: split_words
    :project: progguide
 

doc/src/Errors_messages.rst

@@ -476,65 +476,6 @@ Doc page with :doc:`WARNING messages <Errors_warnings>`
 *Bonds defined but no bond types*
    The data file header lists bonds but no bond types.
 
-*Bond/react: Cannot use fix bond/react with non-molecular systems*
-   Only systems with bonds that can be changed can be used. Atom_style
-   template does not qualify.
-
-*Bond/react: Invalid template atom ID in map file*
-   Atom IDs in molecule templates range from 1 to the number of atoms in the template.
-
-*Bond/react: Rmax cutoff is longer than pairwise cutoff*
-   This is not allowed because bond creation is done using the pairwise
-   neighbor list.
-
-*Bond/react: Molecule template ID for fix bond/react does not exist*
-   A valid molecule template must have been created with the molecule
-   command.
-
-*Bond/react: Reaction templates must contain the same number of atoms*
-   There should be a one-to-one correspondence between atoms in the
-   pre-reacted and post-reacted templates, as specified by the map file.
-
-*Bond/react: Unknown section in map file*
-   Please ensure reaction map files are properly formatted.
-
-*Bond/react: Atom/Bond type affected by reaction too close to template edge*
-   This means an atom which changes type or connectivity during the
-   reaction is too close to an 'edge' atom defined in the map
-   file.  This could cause incorrect assignment of bonds, angle, etc.
-   Generally, this means you must include more atoms in your templates,
-   such that there are at least two atoms between each atom involved in
-   the reaction and an edge atom.
-
-*Bond/react: Fix bond/react needs ghost atoms from farther away*
-   This is because a processor needs to map the entire unreacted
-   molecule template onto simulation atoms it knows about. The
-   comm_modify cutoff command can be used to extend the communication
-   range.
-
-*Bond/react: A deleted atom cannot remain bonded to an atom that is not deleted*
-   Self-explanatory.
-
-*Bond/react: First neighbors of chiral atoms must be of mutually different types*
-   Self-explanatory.
-
-*Bond/react: Chiral atoms must have exactly four first neighbors*
-   Self-explanatory.
-
-*Bond/react: Molecule template 'Coords' section required for chiralIDs keyword*
-   The coordinates of atoms in the pre-reacted template are used to determine
-   chirality.
-
-*Bond/react special bond generation overflow*
-   The number of special bonds per-atom created by a reaction exceeds the
-   system setting. See the read_data or create_box command for how to
-   specify this value.
-
-*Bond/react topology/atom exceed system topology/atom*
-   The number of bonds, angles etc per-atom created by a reaction exceeds
-   the system setting. See the read_data or create_box command for how to
-   specify this value.
-
 *Both restart files must use % or neither*
    Self-explanatory.
 
@@ -1291,7 +1232,7 @@ Doc page with :doc:`WARNING messages <Errors_warnings>`
 *Cannot use chosen neighbor list style with lj/gromacs/kk*
    Self-explanatory.
 
-*Cannot use chosen neighbor list style with lj/sdk/kk*
+*Cannot use chosen neighbor list style with lj/spica/kk*
    That style is not supported by Kokkos.
 
 *Cannot use chosen neighbor list style with pair eam/kk*
@@ -1659,10 +1600,10 @@ Doc page with :doc:`WARNING messages <Errors_warnings>`
 *Cannot use newton pair with lj/gromacs/gpu pair style*
    Self-explanatory.
 
-*Cannot use newton pair with lj/sdk/coul/long/gpu pair style*
+*Cannot use newton pair with lj/spica/coul/long/gpu pair style*
    Self-explanatory.
 
-*Cannot use newton pair with lj/sdk/gpu pair style*
+*Cannot use newton pair with lj/spica/gpu pair style*
    Self-explanatory.
 
 *Cannot use newton pair with lj96/cut/gpu pair style*
@@ -3521,6 +3462,65 @@ Doc page with :doc:`WARNING messages <Errors_warnings>`
    acquire needed info, The comm_modify cutoff command can be used to
    extend the communication range.
 
+*Fix bond/react: Cannot use fix bond/react with non-molecular systems*
+   Only systems with bonds that can be changed can be used. Atom_style
+   template does not qualify.
+
+*Fix bond/react: Invalid template atom ID in map file*
+   Atom IDs in molecule templates range from 1 to the number of atoms in the template.
+
+*Fix bond/react: Rmax cutoff is longer than pairwise cutoff*
+   This is not allowed because bond creation is done using the pairwise
+   neighbor list.
+
+*Fix bond/react: Molecule template ID for fix bond/react does not exist*
+   A valid molecule template must have been created with the molecule
+   command.
+
+*Fix bond/react: Reaction templates must contain the same number of atoms*
+   There should be a one-to-one correspondence between atoms in the
+   pre-reacted and post-reacted templates, as specified by the map file.
+
+*Fix bond/react: Unknown section in map file*
+   Please ensure reaction map files are properly formatted.
+
+*Fix bond/react: Atom/Bond type affected by reaction too close to template edge*
+   This means an atom which changes type or connectivity during the
+   reaction is too close to an 'edge' atom defined in the map
+   file.  This could cause incorrect assignment of bonds, angle, etc.
+   Generally, this means you must include more atoms in your templates,
+   such that there are at least two atoms between each atom involved in
+   the reaction and an edge atom.
+
+*Fix bond/react: Fix bond/react needs ghost atoms from farther away*
+   This is because a processor needs to map the entire unreacted
+   molecule template onto simulation atoms it knows about. The
+   comm_modify cutoff command can be used to extend the communication
+   range.
+
+*Fix bond/react: A deleted atom cannot remain bonded to an atom that is not deleted*
+   Self-explanatory.
+
+*Fix bond/react: First neighbors of chiral atoms must be of mutually different types*
+   Self-explanatory.
+
+*Fix bond/react: Chiral atoms must have exactly four first neighbors*
+   Self-explanatory.
+
+*Fix bond/react: Molecule template 'Coords' section required for chiralIDs keyword*
+   The coordinates of atoms in the pre-reacted template are used to determine
+   chirality.
+
+*Fix bond/react special bond generation overflow*
+   The number of special bonds per-atom created by a reaction exceeds the
+   system setting. See the read_data or create_box command for how to
+   specify this value.
+
+*Fix bond/react topology/atom exceed system topology/atom*
+   The number of bonds, angles etc per-atom created by a reaction exceeds
+   the system setting. See the read_data or create_box command for how to
+   specify this value.
+
 *Fix bond/swap cannot use dihedral or improper styles*
    These styles cannot be defined when using this fix.
 
@@ -6782,7 +6782,7 @@ keyword to allow for additional bonds to be formed
    This is because the computation of constraint forces within a water
    molecule adds forces to atoms owned by other processors.
 
-*Pair style lj/sdk/coul/long/gpu requires atom attribute q*
+*Pair style lj/spica/coul/long/gpu requires atom attribute q*
    The atom style defined does not have this attribute.
 
 *Pair style nb3b/harmonic requires atom IDs*

doc/src/Errors_warnings.rst

@@ -68,14 +68,6 @@ Doc page with :doc:`ERROR messages <Errors_messages>`
    length, multiplying by the number of bonds in the interaction (e.g. 3
    for a dihedral) and adding a small amount of stretch.
 
-*Bond/react: Atom affected by reaction too close to template edge*
-   This means an atom which changes type or connectivity during the
-   reaction is too close to an 'edge' atom defined in the superimpose
-   file. This could cause incorrect assignment of bonds, angle, etc.
-   Generally, this means you must include more atoms in your templates,
-   such that there are at least two atoms between each atom involved in
-   the reaction and an edge atom.
-
 *Both groups in compute group/group have a net charge; the Kspace boundary correction to energy will be non-zero*
    Self-explanatory.
 
@@ -206,12 +198,20 @@ Doc page with :doc:`ERROR messages <Errors_messages>`
 *Fix SRD walls overlap but fix srd overlap not set*
    You likely want to set this in your input script.
 
-* Fix bond/create is used multiple times or with fix bond/break - may not work as expected*
+*Fix bond/create is used multiple times or with fix bond/break - may not work as expected*
    When using fix bond/create multiple times or in combination with
    fix bond/break, the individual fix instances do not share information
    about changes they made at the same time step and thus it may result
    in unexpected behavior.
 
+*Fix bond/react: Atom affected by reaction too close to template edge*
+   This means an atom which changes type or connectivity during the
+   reaction is too close to an 'edge' atom defined in the superimpose
+   file. This could cause incorrect assignment of bonds, angle, etc.
+   Generally, this means you must include more atoms in your templates,
+   such that there are at least two atoms between each atom involved in
+   the reaction and an edge atom.
+
 *Fix bond/swap will ignore defined angles*
    See the page for fix bond/swap for more info on this
    restriction.
@@ -470,6 +470,12 @@ This will most likely cause errors in kinetic fluctuations.
 *More than one compute sna/atom*
    Self-explanatory.
 
+*More than one compute sna/grid*
+   Self-explanatory.
+
+*More than one compute sna/grid/local*
+   Self-explanatory.
+
 *More than one compute snad/atom*
    Self-explanatory.
 
@@ -804,4 +810,3 @@ This will most likely cause errors in kinetic fluctuations.
 
 *Using pair tail corrections with pair_modify compute no*
    The tail corrections will thus not be computed.
-

doc/src/Howto.rst

@@ -65,6 +65,7 @@ Force fields howto
    :maxdepth: 1
 
    Howto_bioFF
+   Howto_amoeba
    Howto_tip3p
    Howto_tip4p
    Howto_spc

doc/src/Howto_amoeba.rst

@@ -0,0 +1,324 @@
+AMOEBA and HIPPO force fields
+=============================
+
+The AMOEBA and HIPPO polarizable force fields were developed by Jay
+Ponder's group at the U Washington at St Louis.  Their implementation
+in LAMMPS was done using F90 code provided by the Ponder group from
+their `Tinker MD code <https://dasher.wustl.edu/tinker/>`_.
+
+The current implementation (July 2022) of AMOEBA in LAMMPS matches the
+version discussed in :ref:`(Ponder) <amoeba-Ponder>`, :ref:`(Ren)
+<amoeba-Ren>`, and :ref:`(Shi) <amoeba-Shi>`.  Likewise the current
+implementation of HIPPO in LAMMPS matches the version discussed in
+:ref:`(Rackers) <amoeba-Rackers>`.
+
+These force fields can be used when polarization effects are desired
+in simulations of water, organic molecules, and biomolecules including
+proteins, provided that parameterizations (Tinker PRM force field
+files) are available for the systems you are interested in.  Files in
+the LAMMPS potentials directory with a "amoeba" or "hippo" suffix can
+be used.  The Tinker distribution and website have additional force
+field files as well:
+`https://github.com/TinkerTools/tinker/tree/release/params
+<https://github.com/TinkerTools/tinker/tree/release/params>`_.
+
+Note that currently, HIPPO can only be used for water systems, but
+HIPPO files for a variety of small organic and biomolecules are in
+preparation by the Ponder group.  Those force field files will be
+included in the LAMMPS distribution when available.
+
+To use the AMOEBA or HIPPO force fields, a simulation must be 3d, and
+fully periodic or fully non-periodic, and use an orthogonal (not
+triclinic) simulation box.
+
+----------
+
+The AMOEBA and HIPPO force fields contain the following terms in their
+energy (U) computation.  Further details for AMOEBA equations are in
+:ref:`(Ponder) <amoeba-Ponder>`, further details for the HIPPO
+equations are in :ref:`(Rackers) <amoeba-Rackers>`.
+
+.. math::
+
+  U & = U_{intermolecular} + U_{intramolecular} \\
+  U_{intermolecular} & = U_{hal} + U_{repulsion} + U_{dispersion} + U_{multipole} + U_{polar} + U_{qxfer} \\
+  U_{intramolecular} & = U_{bond} + U_{angle} + U_{torsion} + U_{oop} + U_{b\theta} + U_{UB} + U_{pitorsion} + U_{bitorsion}
+
+For intermolecular terms, the AMOEBA force field includes only the
+:math:`U_{hal}`, :math:`U_{multipole}`, :math:`U_{polar}` terms.  The
+HIPPO force field includes all but the :math:`U_{hal}` term.  In
+LAMMPS, these are all computed by the :doc:`pair_style amoeba or hippo
+<pair_style>` command.  Note that the :math:`U_{multipole}` and
+:math:`U_{polar}` terms in this formula are not the same for the
+AMOEBA and HIPPO force fields.
+
+For intramolecular terms, the :math:`U_{bond}`, :math:`U_{angle}`,
+:math:`U_{torsion}`, :math:`U_{oop}` terms are computed by the
+:doc:`bond_style class2 <bond_class2>` :doc:`angle_style amoeba
+<angle_amoeba>`, :doc:`dihedral_style fourier <dihedral_fourier>`, and
+:doc:`improper_style amoeba <improper_amoeba>` commands respectively.
+The :doc:`angle_style amoeba <angle_amoeba>` command includes the
+:math:`U_{b\theta}` bond-angle cross term, and the :math:`U_{UB}` term
+for a Urey-Bradley bond contribution between the I,K atoms in the IJK
+angle.
+
+The :math:`U_{pitorsion}` term is computed by the :doc:`fix
+amoeba/pitorsion <fix_amoeba_pitorsion>` command.  It computes 6-body
+interaction between a pair of bonded atoms which each have 2
+additional bond partners.
+
+The :math:`U_{bitorsion}` term is computed by the :doc:`fix
+amoeba/bitorsion <fix_amoeba_bitorsion>` command.  It computes 5-body
+interaction between two 4-body torsions (dihedrals) which overlap,
+having 3 atoms in common.
+
+These command doc pages have additional details on the terms they
+compute:
+
+* :doc:`pair_style amoeba or hippo <pair_amoeba>`
+* :doc:`bond_style class2 <bond_class2>`
+* :doc:`angle_style amoeba <angle_amoeba>`
+* :doc:`dihedral_style fourier <dihedral_fourier>`
+* :doc:`improper_style amoeba <improper_amoeba>`
+* :doc:`fix amoeba/pitorsion <fix_amoeba_pitorsion>`
+* :doc:`fix amoeba/bitorsion <fix_amoeba_bitorsion>`
+
+----------
+
+To use the AMOEBA or HIPPO force fields in LAMMPS, use commands like
+the following appropriately in your input script.  The only change
+needed for AMOEBA vs HIPPO simulation is for the :doc:`pair_style
+<pair_style>` and :doc:`pair_coeff <pair_coeff>` commands, as shown
+below.  See examples/amoeba for example input scripts for both AMOEBA
+and HIPPO.
+
+.. code-block:: LAMMPS
+
+   units              real                           # required
+   atom_style         amoeba
+   bond_style         class2                         # CLASS2 package
+   angle_style        amoeba
+   dihedral_style     fourier                        # EXTRA-MOLECULE package
+   improper_style     amoeba
+                                                     # required per-atom data
+   fix                amtype all property/atom i_amtype ghost yes
+   fix                extra all property/atom &
+                      i_amgroup i_ired i_xaxis i_yaxis i_zaxis d_pval ghost yes
+   fix                polaxe all property/atom i_polaxe
+
+   fix                pit all amoeba/pitorsion       # PiTorsion terms in FF
+   fix_modify         pit energy yes
+                                                     # Bitorsion terms in FF
+   fix                bit all amoeba/bitorsion bitorsion.ubiquitin.data
+   fix_modify         bit energy yes
+
+   read_data          data.ubiquitin fix amtype NULL "Tinker Types" &
+                      fix pit "pitorsion types" "PiTorsion Coeffs" &
+                      fix pit pitorsions PiTorsions &
+                      fix bit bitorsions BiTorsions
+
+   pair_style         amoeba                          # AMOEBA FF
+   pair_coeff         * * amoeba_ubiquitin.prm amoeba_ubiquitin.key
+
+   pair_style         hippo                           # HIPPO FF
+   pair_coeff         * * hippo_water.prm hippo_water.key
+
+   special_bonds      lj/coul 0.5 0.5 0.5 one/five yes     # 1-5 neighbors
+
+The data file read by the :doc:`read_data <read_data>` command should
+be created by the tools/tinker/tinker2lmp.py conversion program
+described below.  It will create a section in the data file with the
+header "Tinker Types".  A :doc:`fix property/atom <fix_property_atom>`
+command for the data must be specified before the read_data command.
+In the example above the fix ID is *amtype*.
+
+Similarly, if the system you are simulating defines AMOEBA/HIPPO
+pitorsion or bitorsion interactions, there will be entries in the data
+file for those interactions.  They require a :doc:`fix
+amoeba/pitortion <fix_amoeba_pitorsion>` and :doc:`fix
+amoeba/bitorsion <fix_amoeba_bitorsion>` command be defined.  In the
+example above, the IDs for these two fixes are *pit* and *bit*.
+
+Of course, if the system being modeled does not have one or more of
+the following -- bond, angle, dihedral, improper, pitorsion,
+bitorsion interactions -- then the corresponding style and fix
+commands above do not need to be used.  See the example scripts in
+examples/amoeba for water systems as examples; they are simpler than
+what is listed above.
+
+The two :doc:`fix property/atom <fix_property_atom>` commands with IDs
+(in the example above) *extra* and *polaxe* are also needed to define
+internal per-atom quantities used by the AMOEBA and HIPPO force
+fields.
+
+The :doc:`pair_coeff <pair_coeff>` command used for either the AMOEBA
+or HIPPO force field takes two arguments for Tinker force field files,
+namely a PRM and KEY file.  The keyfile can be specified as NULL and
+default values for a various settings will be used.  Note that these 2
+files are meant to allow use of native Tinker files as-is.  However
+LAMMPS does not support all the options which can be included
+in a Tinker PRM or KEY file.  See specifics below.
+
+A :doc:`special_bonds <special_bonds>` command with the *one/five*
+option is required, since the AMOEBA/HIPPO force fields define
+weighting factors for not only 1-2, 1-3, 1-4 interactions, but also
+1-5 interactions.  This command will trigger a per-atom list of 1-5
+neighbors to be generated.  The AMOEBA and HIPPO force fields define
+their own custom weighting factors for all the 1-2, 1-3, 1-4, 1-5
+terms which in the Tinker PRM and KEY files; they can be different for
+different terms in the force field.
+
+In addition to the list above, these command doc pages have additional
+details:
+
+* :doc:`atom_style amoeba <atom_style>`
+* :doc:`fix property/atom <fix_property_atom>`
+* :doc:`special_bonds <special_bonds>`
+
+----------
+
+Tinker PRM and KEY files
+
+A Tinker PRM file is composed of sections, each of which has multiple
+lines.  This is the list of PRM sections LAMMPS knows how to parse and
+use.  Any other sections are skipped:
+
+* Angle Bending Parameters
+* Atom Type Definitions
+* Atomic Multipole Parameters
+* Bond Stretching Parameters
+* Charge Penetration Parameters
+* Charge Transfer Parameters
+* Dipole Polarizability Parameters
+* Dispersion Parameters
+* Force Field Definition
+* Literature References
+* Out-of-Plane Bend Parameters
+* Pauli Repulsion Parameters
+* Pi-Torsion Parameters
+* Stretch-Bend Parameters
+* Torsion-Torsion Parameters
+* Torsional Parameters
+* Urey-Bradley Parameters
+* Van der Waals Pair Parameters
+* Van der Waals Parameters
+
+A Tinker KEY file is composed of lines, each of which has a keyword
+followed by zero or more parameters.  This is the list of keywords
+LAMMPS knows how to parse and use in the same manner Tinker does.  Any
+other keywords are skipped.  The value in parenthesis is the default
+value for the keyword if it is not specified, or if the keyfile in the
+:doc:`pair_coeff <pair_coeff>` command is specified as NULL:
+
+* a-axis (0.0)
+* b-axis (0.0)
+* c-axis (0.0)
+* ctrn-cutoff (6.0)
+* ctrn-taper (0.9 * ctrn-cutoff)
+* cutoff
+* delta-halgren (0.07)
+* dewald (no long-range dispersion unless specified)
+* dewald-alpha (0.4)
+* dewald-cutoff (7.0)
+* dispersion-cutoff (9.0)
+* dispersion-taper (9.0 * dispersion-cutoff)
+* dpme-grid
+* dpme-order (4)
+* ewald (no long-range electrostatics unless specified)
+* ewald-alpha (0.4)
+* ewald-cutoff (7.0)
+* gamma-halgren (0.12)
+* mpole-cutoff (9.0)
+* mpole-taper (0.65 * mpole-cutoff)
+* pcg-guess (enabled by default)
+* pcg-noguess (disable pcg-guess if specified)
+* pcg-noprecond (disable pcg-precond if specified)
+* pcg-peek (1.0)
+* pcg-precond (enabled by default)
+* pewald-alpha (0.4)
+* pme-grid
+* pme-order (5)
+* polar-eps (1.0e-6)
+* polar-iter (100)
+* polar-predict (no prediction operation unless specified)
+* ppme-order (5)
+* repulsion-cutoff (6.0)
+* repulsion-taper (0.9 * repulsion-cutoff)
+* taper
+* usolve-cutoff (4.5)
+* usolve-diag (2.0)
+* vdw-cutoff (9.0)
+* vdw-taper (0.9 * vdw-cutoff)
+
+----------
+
+Tinker2lmp.py tool
+
+This conversion tool is found in the tools/tinker directory.
+As shown in examples/amoeba/README, these commands produce
+the data files found in examples/amoeba, and also illustrate
+all the options available to use with the tinker2lmp.py script:
+
+.. code-block:: bash
+
+   % python tinker2lmp.py -xyz water_dimer.xyz -amoeba amoeba_water.prm -data data.water_dimer.amoeba                # AMOEBA non-periodic system
+   % python tinker2lmp.py -xyz water_dimer.xyz -hippo hippo_water.prm -data data.water_dimer.hippo                   # HIPPO non-periodic system
+   % python tinker2lmp.py -xyz water_box.xyz -amoeba amoeba_water.prm -data data.water_box.amoeba -pbc 18.643 18.643 18.643    # AMOEBA periodic system
+   % python tinker2lmp.py -xyz water_box.xyz -hippo hippo_water.prm -data data.water_box.hippo -pbc 18.643 18.643 18.643       # HIPPO periodic system
+   % python tinker2lmp.py -xyz ubiquitin.xyz -amoeba amoeba_ubiquitin.prm -data data.ubiquitin.new -pbc 54.99 41.91 41.91 -bitorsion bitorsion.ubiquitin.data.new   # system with bitorsions
+
+Switches and their arguments may be specified in any order.
+
+The -xyz switch is required and specifies an input XYZ file as an
+argument.  The format of this file is an extended XYZ format defined
+and used by Tinker for its input.  Example \*.xyz files are in the
+examples/amoeba directory.  The file lists the atoms in the system.
+Each atom has the following information: Tinker species name (ignored
+by LAMMPS), xyz coordinates, Tinker numeric type, and a list of atom
+IDs the atom is bonded to.
+
+Here is more information about the extended XYZ format defined and
+used by Tinker, and links to programs that convert standard PDB files
+to the extended XYZ format:
+
+* `http://openbabel.org/docs/current/FileFormats/Tinker_XYZ_format.html <http://openbabel.org/docs/current/FileFormats/Tinker_XYZ_format.html>`_
+* `https://github.com/emleddin/pdbxyz-xyzpdb <https://github.com/emleddin/pdbxyz-xyzpdb>`_
+* `https://github.com/TinkerTools/tinker/blob/release/source/pdbxyz.f <https://github.com/TinkerTools/tinker/blob/release/source/pdbxyz.f>`_
+
+The -amoeba or -hippo switch is required.  It specifies an input
+AMOEBA or HIPPO PRM force field file as an argument.  This should be
+the same file used by the :doc:`pair_style <pair_style>` command in
+the input script.
+
+The -data switch is required.  It specifies an output file name for
+the LAMMPS data file that will be produced.
+
+For periodic systems, the -pbc switch is required.  It specifies the
+periodic box size for each dimension (x,y,z).  For a Tinker simulation
+these are specified in the KEY file.
+
+The -bitorsion switch is only needed if the system contains Tinker
+bitorsion interactions.  The data for each type of bitorsion
+interaction will be written to the specified file, and read by the
+:doc:`fix amoeba/bitorsion <fix_amoeba_bitorsion>` command.  The data
+includes 2d arrays of values to which splines are fit, and thus is not
+compatible with the LAMMPS data file format.
+
+----------
+
+.. _howto-Ponder:
+
+**(Ponder)** Ponder, Wu, Ren, Pande, Chodera, Schnieders, Haque,  Mobley, Lambrecht, DiStasio Jr, M. Head-Gordon, Clark,  Johnson, T. Head-Gordon, J Phys Chem B, 114, 2549-2564 (2010).
+
+.. _howto-Rackers:
+
+**(Rackers)** Rackers, Silva, Wang, Ponder, J Chem Theory Comput, 17, 7056-7084 (2021).
+
+.. _howto-Ren:
+
+**(Ren)** Ren and Ponder, J Phys Chem B, 107, 5933 (2003).
+
+.. _howto-Shi:
+
+**(Shi)** Shi, Xia, Zhang, Best, Wu, Ponder, Ren, J Chem Theory Comp, 9, 4046, 2013.
+

doc/src/Howto_mdi.rst

@@ -5,9 +5,9 @@ Client/server coupling of two (or more) codes is where one code is the
 "client" and sends request messages (data) to one (or more) "server"
 code(s).  A server responds to each request with a reply message
 (data).  This enables two (or more) codes to work in tandem to perform
-a simulation.  LAMMPS can act as either a client or server code; it
-does this by using the `MolSSI Driver Interface (MDI) library
-<https://molssi-mdi.github.io/MDI_Library/html/index.html>`_,
+a simulation.  In this context, LAMMPS can act as either a client or
+server code.  It does this by using the `MolSSI Driver Interface (MDI)
+library <https://molssi-mdi.github.io/MDI_Library/html/index.html>`_,
 developed by the `Molecular Sciences Software Institute (MolSSI)
 <https://molssi.org>`_, which is supported by the :ref:`MDI <PKG-MDI>`
 package.
@@ -63,22 +63,39 @@ The package also provides a :doc:`mdi plugin <mdi>` command which
 enables LAMMPS to operate as an MDI driver and load an MDI engine as a
 plugin library.
 
-The package also has a `fix mdi/aimd <fix_mdi_aimd>` command in which
-LAMMPS operates as an MDI driver to perform *ab initio* MD simulations
-in conjunction with a quantum mechanics code.  Its post_force() method
-illustrates how a driver issues MDI commands to another code.  This
-command can be used to couple to an MDI engine which is either a
-stand-alone code or a plugin library.
+The package also has a `fix mdi/qm <fix_mdi_qm>` command in which
+LAMMPS operates as an MDI driver in conjunction with a quantum
+mechanics code as an MDI engine.  The post_force() method of the
+fix_mdi_qm.cpp file shows how a driver issues MDI commands to another
+code.  This command can be used to couple to an MDI engine which is
+either a stand-alone code or a plugin library.
+
+As explained on the `fix mdi/qm <fix_mdi_qm>` command doc page, it can
+be used to perform *ab initio* MD simulations or energy minimizations,
+or to evaluate the quantum energy and forces for a series of
+independent systems.  The examples/mdi directory has example input
+scripts for all of these use cases.
 
 ----------
 
 The examples/mdi directory contains Python scripts and LAMMPS input
 script which use LAMMPS as either an MDI driver or engine or both.
-Three example use cases are provided:
+Currently, 5 example use cases are provided:
 
-* Run ab initio MD (AIMD) using 2 instances of LAMMPS, one as driver
-  and one as an engine.  As an engine, LAMMPS is a surrogate for a
-  quantum code.
+* Run ab initio MD (AIMD) using 2 instances of LAMMPS.  As a driver
+  LAMMPS performs the timestepping in either NVE or NPT mode.  As an
+  engine, LAMMPS computes forces and is a surrogate for a quantum
+  code.
+
+* As a driver, LAMMPS runs an MD simulation.  Every N steps it passes
+  the current snapshot to an MDI engine to evaluate the energy,
+  virial, and peratom forces.  As the engine LAMMPS is a surrogate for
+  a quantum code.
+
+* As a driver, LAMMPS loops over a series of data files and passes the
+  configuration to an MDI engine to evaluate the energy, virial, and
+  peratom forces.  As the engine LAMMPS is a surrogate for a quantum
+  code.
 
 * A Python script driver invokes a sequence of unrelated LAMMPS
   calculations.  Calculations can be single-point energy/force
@@ -91,20 +108,22 @@ Three example use cases are provided:
 
 Note that in any of these example where LAMMPS is used as an engine,
 an actual QM code (which supports MDI) could be used in its place,
-without modifying other code or scripts, except to specify the name of
-the QM code.
+without modifying the input scripts or launch commands, except to
+specify the name of the QM code.
 
-The examples/mdi/README file explains how to launch both driver and
+The examples/mdi/Run.sh file illustrates how to launch both driver and
 engine codes so that they communicate using the MDI library via either
-MPI or sockets.
+MPI or sockets.  Or using the engine as a stand-alone code or plugin
+library.
 
 -------------
 
-Currently there are two quantum DFT codes which have direct MDI
-support, `Quantum ESPRESSO (QE) <https://www.quantum-espresso.org/>`_
-and `INQ <https://qsg.llnl.gov/node/101.html>`_.  There are also
-several QM codes which have indirect support through QCEngine or i-PI.
-The former means they require a wrapper program (QCEngine) with MDI
+Currently there are at least two quantum DFT codes which have direct
+MDI support, `Quantum ESPRESSO (QE)
+<https://www.quantum-espresso.org/>`_ and `INQ
+<https://qsg.llnl.gov/node/101.html>`_.  There are also several QM
+codes which have indirect support through QCEngine or i-PI.  The
+former means they require a wrapper program (QCEngine) with MDI
 support which writes/read files to pass data to the quantum code
 itself.  The list of QCEngine-supported and i-PI-supported quantum
 codes is on the `MDI webpage

doc/src/Howto_viscosity.rst

@@ -68,7 +68,8 @@ liquid Ar via the GK formalism:
    # Sample LAMMPS input script for viscosity of liquid Ar
 
    units       real
-   variable    T equal 86.4956
+   variable    T equal 200.0       # run temperature
+   variable    Tinit equal 250.0   # equilibration temperature
    variable    V equal vol
    variable    dt equal 4.0
    variable    p equal 400     # correlation length
@@ -99,12 +100,14 @@ liquid Ar via the GK formalism:
 
    # equilibration and thermalization
 
-   velocity     all create $T 102486 mom yes rot yes dist gaussian
-   fix          NVT all nvt temp $T $T 10 drag 0.2
+   velocity     all create ${Tinit} 102486 mom yes rot yes dist gaussian
+   fix          NVT all nvt temp ${Tinit} ${Tinit} 10 drag 0.2
    run          8000
 
    # viscosity calculation, switch to NVE if desired
 
+   velocity     all create $T 102486 mom yes rot yes dist gaussian
+   fix          NVT all nvt temp $T $T 10 drag 0.2
    #unfix       NVT
    #fix         NVE all nve
 
@@ -122,7 +125,7 @@ liquid Ar via the GK formalism:
    run          100000
    variable     v equal (v_v11+v_v22+v_v33)/3.0
    variable     ndens equal count(all)/vol
-   print        "average viscosity: $v [Pa.s] @ $T K, ${ndens} /A^3"
+   print        "average viscosity: $v [Pa.s] @ $T K, ${ndens} atoms/A^3"
 
 The fifth method is related to the above Green-Kubo method,
 but uses the Einstein formulation, analogous to the Einstein
@@ -131,9 +134,9 @@ time-integrated momentum fluxes play the role of Cartesian
 coordinates, whose mean-square displacement increases linearly
 with time at sufficiently long times.
 
-The sixth is periodic perturbation method. It is also a non-equilibrium MD method.
-However, instead of measure the momentum flux in response of applied velocity gradient,
-it measures the velocity profile in response of applied stress.
+The sixth is the periodic perturbation method, which is also a non-equilibrium MD method.
+However, instead of measuring the momentum flux in response to an applied velocity gradient,
+it measures the velocity profile in response to applied stress.
 A cosine-shaped periodic acceleration is added to the system via the
 :doc:`fix accelerate/cos <fix_accelerate_cos>` command,
 and the :doc:`compute viscosity/cos<compute_viscosity_cos>` command is used to monitor the

doc/src/Install.rst

@@ -3,10 +3,20 @@ Install LAMMPS
 
 You can download LAMMPS as an executable or as source code.
 
-With source code, you also have to :doc:`build LAMMPS <Build>`.  But you
-have more flexibility as to what features to include or exclude in the
-build.  If you plan to :doc:`modify or extend LAMMPS <Modify>`, then you
-need the source code.
+When downloading the LAMMPS source code, you also have to :doc:`build
+LAMMPS <Build>`.  But you have more flexibility as to what features to
+include or exclude in the build.  When you download and install
+pre-compiled LAMMPS executables, you are limited to install which
+version of LAMMPS is available and which features are included of these
+builds.  If you plan to :doc:`modify or extend LAMMPS <Modify>`, then
+you **must** build LAMMPS from the source code.
+
+.. note::
+
+   If you have questions about the pre-compiled LAMMPS executables, you
+   need to contact the people preparing those executables.  The LAMMPS
+   developers have no control over their choices of how they configure
+   and build their packages and when they update them.
 
 .. toctree::
    :maxdepth: 1

doc/src/Install_conda.rst

@@ -38,3 +38,10 @@ up the Conda capability.
 .. _openkim: https://openkim.org
 .. _conda: https://docs.conda.io/en/latest/index.html
 .. _mini_conda_install: https://docs.conda.io/en/latest/miniconda.html
+
+.. note::
+
+   If you have questions about these pre-compiled LAMMPS executables,
+   you need to contact the people preparing those packages.  The LAMMPS
+   developers have no control over their choices of how they configure
+   and build their packages and when they update them.

doc/src/Install_linux.rst

@@ -3,13 +3,19 @@ Download an executable for Linux
 
 Binaries are available for different versions of Linux:
 
-| :ref:`Pre-built Ubuntu Linux executables <ubuntu>`
-| :ref:`Pre-built Fedora Linux executables <fedora>`
-| :ref:`Pre-built EPEL Linux executables (RHEL, CentOS) <epel>`
-| :ref:`Pre-built OpenSuse Linux executables <opensuse>`
-| :ref:`Gentoo Linux executable <gentoo>`
-| :ref:`Arch Linux build-script <arch>`
-|
+- :ref:`Pre-built Ubuntu Linux executables <ubuntu>`
+- :ref:`Pre-built Fedora Linux executables <fedora>`
+- :ref:`Pre-built EPEL Linux executables (RHEL, CentOS) <epel>`
+- :ref:`Pre-built OpenSuse Linux executables <opensuse>`
+- :ref:`Gentoo Linux executable <gentoo>`
+- :ref:`Arch Linux build-script <arch>`
+
+.. note::
+
+   If you have questions about these pre-compiled LAMMPS executables,
+   you need to contact the people preparing those packages.  The LAMMPS
+   developers have no control over their choices of how they configure
+   and build their packages and when they update them.
 
 ----------
 
@@ -18,41 +24,28 @@ Binaries are available for different versions of Linux:
 Pre-built Ubuntu Linux executables
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
-A pre-built LAMMPS executable suitable for running on the latest
-Ubuntu Linux versions, can be downloaded as a Debian package.  This
-allows you to install LAMMPS with a single command, and stay
-up-to-date with the current stable version of LAMMPS by simply updating
-your operating system.  Please note, that the repository below offers
-two LAMMPS packages, ``lammps-daily`` and ``lammps-stable``.  The
-LAMMPS developers recommend to use the ``lammps-stable`` package for
-any production simulations.  The ``lammps-daily`` package is built
-from the LAMMPS development sources, and those versions may have known
-issues and bugs when new features are added and the software has not
-undergone full release testing.
-
-To install the appropriate personal-package archives (PPAs), do the
-following once:
-
-.. code-block:: bash
-
-   $ sudo add-apt-repository ppa:gladky-anton/lammps
-   $ sudo add-apt-repository ppa:openkim/latest
-   $ sudo apt-get update
+A pre-built LAMMPS executable suitable for running on the latest Ubuntu
+Linux versions, can be downloaded as a Debian package.  This allows you
+to install LAMMPS with a single command, and stay (mostly) up-to-date
+with the current stable version of LAMMPS by simply updating your
+operating system.
 
 To install LAMMPS do the following once:
 
 .. code-block:: bash
 
-   $ sudo apt-get install lammps-stable
+   $ sudo apt-get install lammps
 
-This downloads an executable named ``lmp_stable`` to your box, which
-can then be used in the usual way to run input scripts:
+This downloads an executable named ``lmp`` to your box and multiple
+packages with supporting data, examples and libraries as well as any
+missing dependencies.  This executable can then be used in the usual way
+to run input scripts:
 
 .. code-block:: bash
 
-   $ lmp_stable -in in.lj
+   $ lmp -in in.lj
 
-To update LAMMPS to the most current stable version, do the following:
+To update LAMMPS to the latest packaged version, do the following:
 
 .. code-block:: bash
 
@@ -60,44 +53,24 @@ To update LAMMPS to the most current stable version, do the following:
 
 which will also update other packages on your system.
 
-To get a copy of the current documentation and examples:
-
-.. code-block:: bash
-
-   $ sudo apt-get install lammps-stable-doc
-
-which will download the doc files in
-``/usr/share/doc/lammps-stable-doc/doc`` and example problems in
-``/usr/share/doc/lammps-doc/examples``.
-
-To get a copy of the current potentials files:
-
-.. code-block:: bash
-
-   $ sudo apt-get install lammps-stable-data
-
-which will download the potentials files to
-``/usr/share/lammps-stable/potentials``.  The ``lmp_stable`` binary is
-hard-coded to look for potential files in this directory (it does not
-use the ``LAMMPS_POTENTIALS`` environment variable, as described
-in :doc:`pair_coeff <pair_coeff>` command).
-
-The ``lmp_stable`` binary is built with the :ref:`KIM package <kim>` which
-results in the above command also installing the ``kim-api`` binaries when LAMMPS
-is installed.  In order to use potentials from `openkim.org <openkim_>`_, you
-can install the ``openkim-models`` package
+The ``lmp`` binary is built with the :ref:`KIM package <kim>` included,
+which results in the above command also installing the ``kim-api``
+binaries when LAMMPS is installed.  In order to use potentials from
+`openkim.org <openkim_>`_, you can also install the ``openkim-models``
+package
 
 .. code-block:: bash
 
    $ sudo apt-get install openkim-models
 
+Or use the KIM-API commands to download and install individual models.
 To un-install LAMMPS, do the following:
 
 .. code-block:: bash
 
-   $ sudo apt-get remove lammps-stable
+   $ sudo apt-get remove lammps
 
-Please use ``lmp_stable -help`` to see which compilation options, packages,
+Please use ``lmp -help`` to see which compilation options, packages,
 and styles are included in the binary.
 
 Thanks to Anton Gladky (gladky.anton at gmail.com) for setting up this
@@ -110,21 +83,21 @@ Ubuntu package capability.
 Pre-built Fedora Linux executables
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
-Pre-built LAMMPS packages for stable releases are available
-in the Fedora Linux distribution as of version 28. The packages
-can be installed via the dnf package manager. There are 3 basic
-varieties (lammps = no MPI, lammps-mpich = MPICH MPI library,
-lammps-openmpi = OpenMPI MPI library) and for each support for
-linking to the C library interface (lammps-devel, lammps-mpich-devel,
-lammps-openmpi-devel), the header for compiling programs using
-the C library interface (lammps-headers), and the LAMMPS python
-module for Python 3. All packages can be installed at the same
-time and the name of the LAMMPS executable is ``lmp`` and ``lmp_openmpi``
-or ``lmp_mpich`` respectively.  By default, ``lmp`` will refer to the
-serial executable, unless one of the MPI environment modules is loaded
-(``module load mpi/mpich-x86_64`` or ``module load mpi/openmpi-x86_64``).
-Then the corresponding parallel LAMMPS executable can be used.
-The same mechanism applies when loading the LAMMPS python module.
+Pre-built LAMMPS packages for stable releases are available in the
+Fedora Linux distribution as of Fedora version 28. The packages can be
+installed via the dnf package manager. There are 3 basic varieties
+(lammps = no MPI, lammps-mpich = MPICH MPI library, lammps-openmpi =
+OpenMPI MPI library) and for each support for linking to the C library
+interface (lammps-devel, lammps-mpich-devel, lammps-openmpi-devel), the
+header for compiling programs using the C library interface
+(lammps-headers), and the LAMMPS python module for Python 3. All
+packages can be installed at the same time and the name of the LAMMPS
+executable is ``lmp`` and ``lmp_openmpi`` or ``lmp_mpich`` respectively.
+By default, ``lmp`` will refer to the serial executable, unless one of
+the MPI environment modules is loaded (``module load mpi/mpich-x86_64``
+or ``module load mpi/openmpi-x86_64``).  Then the corresponding parallel
+LAMMPS executable can be used.  The same mechanism applies when loading
+the LAMMPS python module.
 
 To install LAMMPS with OpenMPI and run an input ``in.lj`` with 2 CPUs do:
 
@@ -273,3 +246,10 @@ Alternatively, you may use an AUR helper to install these packages.
 
 Note that the AUR provides build-scripts that download the source and
 the build the package on your machine.
+
+.. note::
+
+   It looks like the Arch Linux AUR repository build scripts for LAMMPS
+   have not been updated since the 29 October 2020 version.  You may want
+   to consider installing a more current version of LAMMPS from source
+   directly.

doc/src/Intro_citing.rst

@@ -33,9 +33,9 @@ initial versions of LAMMPS is:
 DOI for the LAMMPS source code
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
-LAMMPS developers use the `Zenodo service at CERN <https://zenodo.org/>`_
-to create digital object identifies (DOI) for stable releases of the
-LAMMPS source code. There are two types of DOIs for the LAMMPS source code.
+The LAMMPS developers use the `Zenodo service at CERN <https://zenodo.org/>`_
+to create digital object identifiers (DOI) for stable releases of the
+LAMMPS source code.  There are two types of DOIs for the LAMMPS source code.
 
 The canonical DOI for **all** versions of LAMMPS, which will always
 point to the **latest** stable release version is:

doc/src/Manual.rst

@@ -110,6 +110,7 @@ Command Reference
    angles
    dihedrals
    impropers
+   dumps
    fix_modify_atc_commands
    Bibliography
 

doc/src/Modify_fix.rst

@@ -23,6 +23,8 @@ derived class.  See fix.h for details.
 +---------------------------+--------------------------------------------------------------------------------------------+
 | init                      | initialization before a run (optional)                                                     |
 +---------------------------+--------------------------------------------------------------------------------------------+
+| init_list                 | store pointer to neighbor list; called by neighbor list code (optional)                    |
++---------------------------+--------------------------------------------------------------------------------------------+
 | setup_pre_exchange        | called before atom exchange in setup (optional)                                            |
 +---------------------------+--------------------------------------------------------------------------------------------+
 | setup_pre_force           | called before force computation in setup (optional)                                        |

doc/src/Modify_style.rst

@@ -100,13 +100,14 @@ Documentation (strict)
 
 Contributions that add new styles or commands or augment existing ones
 must include the corresponding new or modified documentation in
-`ReStructuredText format <rst>`_ (.rst files in the ``doc/src/`` folder). The
-documentation shall be written in American English and the .rst file
-must use only ASCII characters so it can be cleanly translated to PDF
-files (via `sphinx <sphinx>`_ and PDFLaTeX).  Special characters may be included via
-embedded math expression typeset in a LaTeX subset.
+`ReStructuredText format <rst_>`_ (.rst files in the ``doc/src/``
+folder). The documentation shall be written in American English and the
+.rst file must use only ASCII characters so it can be cleanly translated
+to PDF files (via `sphinx <https://www.sphinx-doc.org>`_ and PDFLaTeX).
+Special characters may be included via embedded math expression typeset
+in a LaTeX subset.
 
-.. _rst: https://docutils.readthedocs.io/en/sphinx-docs/user/rst/quickstart.html
+.. _rst: https://www.sphinx-doc.org/en/master/usage/restructuredtext/index.html
 
 When adding new commands, they need to be integrated into the sphinx
 documentation system, and the corresponding command tables and lists
@@ -133,7 +134,7 @@ error free completion of the HTML and PDF build will be performed and
 also a spell check, a check for correct anchors and labels, and a check
 for completeness of references all styles in their corresponding tables
 and lists is run.  In case the spell check reports false positives they
-can be added to the file doc/utils/sphinx-config/false_positives.txt
+can be added to the file ``doc/utils/sphinx-config/false_positives.txt``
 
 Contributions that add or modify the library interface or "public" APIs
 from the C++ code or the Fortran module must include suitable doxygen

doc/src/Packages_details.rst

@@ -27,6 +27,7 @@ page gives those details.
    :columns: 6
 
    * :ref:`ADIOS <PKG-ADIOS>`
+   * :ref:`AMOEBA <PKG-AMOEBA>`
    * :ref:`ASPHERE <PKG-ASPHERE>`
    * :ref:`ATC <PKG-ATC>`
    * :ref:`AWPMD <PKG-AWPMD>`
@@ -35,7 +36,7 @@ page gives those details.
    * :ref:`BPM <PKG-BPM>`
    * :ref:`BROWNIAN <PKG-BROWNIAN>`
    * :ref:`CG-DNA <PKG-CG-DNA>`
-   * :ref:`CG-SDK <PKG-CG-SDK>`
+   * :ref:`CG-SPICA <PKG-CG-SPICA>`
    * :ref:`CLASS2 <PKG-CLASS2>`
    * :ref:`COLLOID <PKG-COLLOID>`
    * :ref:`COLVARS <PKG-COLVARS>`
@@ -149,6 +150,24 @@ This package has :ref:`specific installation instructions <adios>` on the :doc:`
 
 ----------
 
+.. _PKG-AMOEBA:
+
+AMOEBA package
+---------------
+
+**Contents:**
+
+TODO
+
+**Supporting info:**
+
+* src/AMOEBA: filenames -> commands
+* :doc:`AMOEBA and HIPPO howto <Howto_amoeba>`
+* examples/amoeba
+* TODO
+
+----------
+
 .. _PKG-ASPHERE:
 
 ASPHERE package
@@ -298,6 +317,8 @@ models for mesoscale simulations of solids and fracture.  See the
 
 **Authors:** Joel T. Clemmer (Sandia National Labs)
 
+.. versionadded:: 4May2022
+
 **Supporting info:**
 
 * src/BPM filenames -> commands
@@ -365,28 +386,30 @@ The CG-DNA package requires that also the `MOLECULE <PKG-MOLECULE>`_ and
 
 ----------
 
-.. _PKG-CG-SDK:
+.. _PKG-CG-SPICA:
 
-CG-SDK package
+CG-SPICA package
 ------------------
 
 **Contents:**
 
 Several pair styles and an angle style which implement the
-coarse-grained SDK model of Shinoda, DeVane, and Klein which enables
-simulation of ionic liquids, electrolytes, lipids and charged amino
-acids.
+coarse-grained SPICA (formerly called SDK) model which enables
+simulation of biological or soft material systems.
 
-**Author:** Axel Kohlmeyer (Temple U).
+**Original Author:** Axel Kohlmeyer (Temple U).
+
+**Maintainers:** Yusuke Miyazaki and Wataru Shinoda (Okayama U).
 
 **Supporting info:**
 
-* src/CG-SDK: filenames -> commands
-* src/CG-SDK/README
-* :doc:`pair_style lj/sdk/\* <pair_sdk>`
-* :doc:`angle_style sdk <angle_sdk>`
-* examples/PACKAGES/cgsdk
+* src/CG-SPICA: filenames -> commands
+* src/CG-SPICA/README
+* :doc:`pair_style lj/spica/\* <pair_spica>`
+* :doc:`angle_style spica <angle_spica>`
+* examples/PACKAGES/cgspica
 * https://www.lammps.org/pictures.html#cg
+* https://www.spica-ff.org/
 
 ----------
 
@@ -824,6 +847,8 @@ groups of atoms that interact with the remaining atoms as electrolyte.
 Ahrens-Iwers (TUHH, Hamburg, Germany), Shern Tee (UQ, Brisbane, Australia) and
 Robert Meissner (TUHH, Hamburg, Germany).
 
+.. versionadded:: 4May2022
+
 **Install:**
 
 This package has :ref:`specific installation instructions <electrode>` on the
@@ -1479,7 +1504,7 @@ the :doc:`Build extras <Build_extras>` page.
 * lib/mdi/README
 * :doc:`Howto MDI <Howto_mdi>`
 * :doc:`mdi <mdi>`
-* :doc:`fix mdi/aimd <fix_mdi_aimd>`
+* :doc:`fix mdi/qm <fix_mdi_qm>`
 * examples/PACKAGES/mdi
 
 ----------
@@ -1801,6 +1826,8 @@ computes which analyze attributes of the potential.
 * src/ML-SNAP: filenames -> commands
 * :doc:`pair_style snap <pair_snap>`
 * :doc:`compute sna/atom <compute_sna_atom>`
+* :doc:`compute sna/grid <compute_sna_atom>`
+* :doc:`compute sna/grid/local <compute_sna_atom>`
 * :doc:`compute snad/atom <compute_sna_atom>`
 * :doc:`compute snav/atom <compute_sna_atom>`
 * examples/snap

doc/src/Packages_list.rst

@@ -33,6 +33,11 @@ whether an extra library is needed to build and use the package:
      - :doc:`dump adios <dump_adios>`
      - PACKAGES/adios
      - ext
+   * - :ref:`AMOEBA <PKG-AMOEBA>`
+     - AMOEBA and HIPPO force fields
+     - :doc:`AMOEBA and HIPPO howto <Howto_amoeba>`
+     - amoeba
+     - no
    * - :ref:`ASPHERE <PKG-ASPHERE>`
      - aspherical particle models
      - :doc:`Howto spherical <Howto_spherical>`
@@ -73,10 +78,10 @@ whether an extra library is needed to build and use the package:
      - src/CG-DNA/README
      - PACKAGES/cgdna
      - no
-   * - :ref:`CG-SDK <PKG-CG-SDK>`
-     - SDK coarse-graining model
-     - :doc:`pair_style lj/sdk <pair_sdk>`
-     - PACKAGES/cgsdk
+   * - :ref:`CG-SPICA <PKG-CG-SPICA>`
+     - SPICA (SDK) coarse-graining model
+     - :doc:`pair_style lj/spica <pair_spica>`
+     - PACKAGES/cgspica
      - no
    * - :ref:`CLASS2 <PKG-CLASS2>`
      - class 2 force fields

doc/src/angle_amoeba.rst

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

doc/src/angle_class2.rst

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

doc/src/angle_spica.rst

@@ -1,32 +1,32 @@
-.. index:: angle_style sdk
-.. index:: angle_style sdk/omp
+.. index:: angle_style spica
+.. index:: angle_style spica/omp
 
-angle_style sdk command
-=======================
+angle_style spica command
+=========================
 
-Accelerator Variants: *sdk/omp*
+Accelerator Variants: *spica/omp*
 
 Syntax
 """"""
 
 .. code-block:: LAMMPS
 
-   angle_style sdk
+   angle_style spica
 
-   angle_style sdk/omp
+   angle_style spica/omp
 
 Examples
 """"""""
 
 .. code-block:: LAMMPS
 
-   angle_style sdk
+   angle_style spica
    angle_coeff 1 300.0 107.0
 
 Description
 """""""""""
 
-The *sdk* angle style is a combination of the harmonic angle potential,
+The *spica* angle style is a combination of the harmonic angle potential,
 
 .. math::
 
@@ -34,10 +34,10 @@ The *sdk* angle style is a combination of the harmonic angle potential,
 
 where :math:`\theta_0` is the equilibrium value of the angle and
 :math:`K` a prefactor, with the *repulsive* part of the non-bonded
-*lj/sdk* pair style between the atoms 1 and 3.  This angle potential is
-intended for coarse grained MD simulations with the CMM parameterization
-using the :doc:`pair_style lj/sdk <pair_sdk>`.  Relative to the
-pair_style *lj/sdk*, however, the energy is shifted by
+*lj/spica* pair style between the atoms 1 and 3.  This angle potential is
+intended for coarse grained MD simulations with the SPICA (formerly called SDK) parameterization
+using the :doc:`pair_style lj/spica <pair_spica>`.  Relative to the
+pair_style *lj/spica*, however, the energy is shifted by
 :math:`\epsilon`, to avoid sudden jumps.  Note that the usual 1/2 factor
 is included in :math:`K`.
 
@@ -51,9 +51,12 @@ The following coefficients must be defined for each angle type via the
 radians internally; hence :math:`K` is effectively energy per
 radian\^2.
 
-The required *lj/sdk* parameters are extracted automatically from the
+The required *lj/spica* parameters are extracted automatically from the
 pair_style.
 
+Style *sdk*, the original implementation of style *spica*, is available
+for backward compatibility.
+
 ----------
 
 .. include:: accel_styles.rst
@@ -64,14 +67,14 @@ Restrictions
 """"""""""""
 
 This angle style can only be used if LAMMPS was built with the
-CG-SDK package.  See the :doc:`Build package <Build_package>` doc
+CG-SPICA package.  See the :doc:`Build package <Build_package>` doc
 page for more info.
 
 Related commands
 """"""""""""""""
 
-:doc:`angle_coeff <angle_coeff>`, :doc:`angle_style harmonic <angle_harmonic>`, :doc:`pair_style lj/sdk <pair_sdk>`,
-:doc:`pair_style lj/sdk/coul/long <pair_sdk>`
+:doc:`angle_coeff <angle_coeff>`, :doc:`angle_style harmonic <angle_harmonic>`, :doc:`pair_style lj/spica <pair_spica>`,
+:doc:`pair_style lj/spica/coul/long <pair_spica>`
 
 Default
 """""""

doc/src/angle_style.rst

@@ -73,6 +73,7 @@ of (g,i,k,o,t) to indicate which accelerated styles exist.
 * :doc:`zero <angle_zero>` - topology but no interactions
 * :doc:`hybrid <angle_hybrid>` - define multiple styles of angle interactions
 
+* :doc:`amoeba <angle_amoeba>` - AMOEBA angle
 * :doc:`charmm <angle_charmm>` - CHARMM angle
 * :doc:`class2 <angle_class2>` - COMPASS (class 2) angle
 * :doc:`class2/p6 <angle_class2>` - COMPASS (class 2) angle expanded to 6th order
@@ -91,7 +92,7 @@ of (g,i,k,o,t) to indicate which accelerated styles exist.
 * :doc:`harmonic <angle_harmonic>` - harmonic angle
 * :doc:`mm3 <angle_mm3>` - anharmonic angle
 * :doc:`quartic <angle_quartic>` - angle with cubic and quartic terms
-* :doc:`sdk <angle_sdk>` - harmonic angle with repulsive SDK pair style between 1-3 atoms
+* :doc:`spica <angle_spica>` - harmonic angle with repulsive SPICA pair style between 1-3 atoms
 * :doc:`table <angle_table>` - tabulated by angle
 
 ----------

doc/src/atom_style.rst

@@ -10,7 +10,7 @@ Syntax
 
    atom_style style args
 
-* style = *angle* or *atomic* or *body* or *bond* or *charge* or *dipole* or  *dpd* or *edpd* or *electron* or *ellipsoid* or *full* or *line* or *mdpd* or *molecular* or *oxdna* or *peri* or *smd* or *sph* or *sphere* or *bpm/sphere* or *spin* or *tdpd* or *tri* or *template* or *hybrid*
+* style = *amoeba* or *angle* or *atomic* or *body* or *bond* or *charge* or *dipole* or *dpd* or *edpd* or *electron* or *ellipsoid* or *full* or *line* or *mdpd* or *molecular* or *oxdna* or *peri* or *smd* or *sph* or *sphere* or or *bpm/sphere* or *spin* or *tdpd* or *tri* or *template* or *hybrid*
 
   .. parsed-literal::
 
@@ -78,6 +78,8 @@ coordinates, velocities, atom IDs and types.  See the
 :doc:`set <set>` commands for info on how to set these various
 quantities.
 
++--------------+-----------------------------------------------------+--------------------------------------+
+| *amoeba*     | molecular + charge + 1/5 neighbors                  | AMOEBA/HIPPO polarized force fields  |
 +--------------+-----------------------------------------------------+--------------------------------------+
 | *angle*      | bonds and angles                                    | bead-spring polymers with stiffness  |
 +--------------+-----------------------------------------------------+--------------------------------------+
@@ -137,11 +139,13 @@ quantities.
 .. note::
 
    It is possible to add some attributes, such as a molecule ID, to
-   atom styles that do not have them via the :doc:`fix property/atom <fix_property_atom>` command.  This command also
-   allows new custom attributes consisting of extra integer or
-   floating-point values to be added to atoms.  See the :doc:`fix property/atom <fix_property_atom>` page for examples of cases
-   where this is useful and details on how to initialize, access, and
-   output the custom values.
+   atom styles that do not have them via the :doc:`fix property/atom
+   <fix_property_atom>` command.  This command also allows new custom
+   attributes consisting of extra integer or floating-point values to
+   be added to atoms.  See the :doc:`fix property/atom
+   <fix_property_atom>` page for examples of cases where this is
+   useful and details on how to initialize, access, and output the
+   custom values.
 
 All of the above styles define point particles, except the *sphere*,
 *bpm/sphere*, *ellipsoid*, *electron*, *peri*, *wavepacket*, *line*,
@@ -154,19 +158,20 @@ per-type basis, using the :doc:`mass <mass>` command, The finite-size
 particle styles assign mass to individual particles on a per-particle
 basis.
 
-For the *sphere* and *bpm/sphere* styles, the particles are spheres and each stores a
-per-particle diameter and mass.  If the diameter > 0.0, the particle
-is a finite-size sphere.  If the diameter = 0.0, it is a point
-particle.  Note that by use of the *disc* keyword with the :doc:`fix
-nve/sphere <fix_nve_sphere>`, :doc:`fix nvt/sphere <fix_nvt_sphere>`,
-:doc:`fix nph/sphere <fix_nph_sphere>`, :doc:`fix npt/sphere
-<fix_npt_sphere>` commands for the *sphere* style, spheres can be effectively treated as 2d
-discs for a 2d simulation if desired.  See also the :doc:`set
-density/disc <set>` command.  The *sphere* and *bpm/sphere* styles take an optional 0
-or 1 argument.  A value of 0 means the radius of each sphere is
-constant for the duration of the simulation.  A value of 1 means the
-radii may vary dynamically during the simulation, e.g. due to use of
-the :doc:`fix adapt <fix_adapt>` command.
+For the *sphere* and *bpm/sphere* styles, the particles are spheres
+and each stores a per-particle diameter and mass.  If the diameter >
+0.0, the particle is a finite-size sphere.  If the diameter = 0.0, it
+is a point particle.  Note that by use of the *disc* keyword with the
+:doc:`fix nve/sphere <fix_nve_sphere>`, :doc:`fix nvt/sphere
+<fix_nvt_sphere>`, :doc:`fix nph/sphere <fix_nph_sphere>`, :doc:`fix
+npt/sphere <fix_npt_sphere>` commands for the *sphere* style, spheres
+can be effectively treated as 2d discs for a 2d simulation if desired.
+See also the :doc:`set density/disc <set>` command.  The *sphere* and
+*bpm/sphere* styles take an optional 0 or 1 argument.  A value of 0
+means the radius of each sphere is constant for the duration of the
+simulation.  A value of 1 means the radii may vary dynamically during
+the simulation, e.g. due to use of the :doc:`fix adapt <fix_adapt>`
+command.
 
 For the *ellipsoid* style, the particles are ellipsoids and each
 stores a flag which indicates whether it is a finite-size ellipsoid or
@@ -175,15 +180,16 @@ vector with the 3 diameters of the ellipsoid and a quaternion 4-vector
 with its orientation.
 
 For the *dielectric* style, each particle can be either a physical
-particle (e.g. an ion), or an interface particle representing a boundary
-element. For physical particles, the per-particle properties are
-the same as atom_style full.  For interface particles, in addition to
-these properties, each particle also has an area, a normal unit vector,
-a mean local curvature, the mean and difference of the dielectric constants
-of two sides of the interface, and the local dielectric constant at the
-boundary element.  The distinction between the physical and interface
-particles is only meaningful when :doc:`fix polarize <fix_polarize>`
-commands are applied to the interface particles.
+particle (e.g. an ion), or an interface particle representing a
+boundary element. For physical particles, the per-particle properties
+are the same as atom_style full.  For interface particles, in addition
+to these properties, each particle also has an area, a normal unit
+vector, a mean local curvature, the mean and difference of the
+dielectric constants of two sides of the interface, and the local
+dielectric constant at the boundary element.  The distinction between
+the physical and interface particles is only meaningful when :doc:`fix
+polarize <fix_polarize>` commands are applied to the interface
+particles.
 
 For the *dipole* style, a point dipole is defined for each point
 particle.  Note that if you wish the particles to be finite-size
@@ -272,16 +278,17 @@ showing the use of the *template* atom style versus *molecular*.
 
 .. note::
 
-   When using the *template* style with a :doc:`molecule template <molecule>` that contains multiple molecules, you should
-   insure the atom types, bond types, angle_types, etc in all the
-   molecules are consistent.  E.g. if one molecule represents H2O and
-   another CO2, then you probably do not want each molecule file to
-   define 2 atom types and a single bond type, because they will conflict
-   with each other when a mixture system of H2O and CO2 molecules is
-   defined, e.g. by the :doc:`read_data <read_data>` command.  Rather the
-   H2O molecule should define atom types 1 and 2, and bond type 1.  And
-   the CO2 molecule should define atom types 3 and 4 (or atom types 3 and
-   2 if a single oxygen type is desired), and bond type 2.
+   When using the *template* style with a :doc:`molecule template
+   <molecule>` that contains multiple molecules, you should insure the
+   atom types, bond types, angle_types, etc in all the molecules are
+   consistent.  E.g. if one molecule represents H2O and another CO2,
+   then you probably do not want each molecule file to define 2 atom
+   types and a single bond type, because they will conflict with each
+   other when a mixture system of H2O and CO2 molecules is defined,
+   e.g. by the :doc:`read_data <read_data>` command.  Rather the H2O
+   molecule should define atom types 1 and 2, and bond type 1.  And
+   the CO2 molecule should define atom types 3 and 4 (or atom types 3
+   and 2 if a single oxygen type is desired), and bond type 2.
 
 For the *body* style, the particles are arbitrary bodies with internal
 attributes defined by the "style" of the bodies, which is specified by
@@ -339,6 +346,8 @@ Many of the styles listed above are only enabled if LAMMPS was built
 with a specific package, as listed below.  See the :doc:`Build package
 <Build_package>` page for more info.
 
+The *amoeba* style is part of the AMOEBA package.
+
 The *angle*, *bond*, *full*, *molecular*, and *template* styles are
 part of the MOLECULE package.
 
@@ -350,9 +359,11 @@ The *dipole* style is part of the DIPOLE package.
 
 The *peri* style is part of the PERI package for Peridynamics.
 
-The *oxdna* style is part of the CG-DNA package for coarse-grained simulation of DNA and RNA.
+The *oxdna* style is part of the CG-DNA package for coarse-grained
+simulation of DNA and RNA.
 
-The *electron* style is part of the EFF package for :doc:`electronic force fields <pair_eff>`.
+The *electron* style is part of the EFF package for :doc:`electronic
+force fields <pair_eff>`.
 
 The *dpd* style is part of the DPD-REACT package for dissipative
 particle dynamics (DPD).
@@ -363,7 +374,8 @@ dissipative particle dynamics (mDPD), and transport dissipative particle
 dynamics (tDPD), respectively.
 
 The *sph* style is part of the SPH package for smoothed particle
-hydrodynamics (SPH).  See `this PDF guide <PDF/SPH_LAMMPS_userguide.pdf>`_ to using SPH in LAMMPS.
+hydrodynamics (SPH).  See `this PDF guide
+<PDF/SPH_LAMMPS_userguide.pdf>`_ to using SPH in LAMMPS.
 
 The *mesont* style is part of the MESONT package.
 

doc/src/compute.rst

@@ -284,6 +284,8 @@ The individual style names on the :doc:`Commands compute <Commands_compute>` pag
 * :doc:`smd/vol <compute_smd_vol>` - per-particle volumes and their sum in Smooth Mach Dynamics
 * :doc:`snap <compute_sna_atom>` - gradients of SNAP energy and forces w.r.t. linear coefficients and related quantities for fitting SNAP potentials
 * :doc:`sna/atom <compute_sna_atom>` - bispectrum components for each atom
+* :doc:`sna/grid <compute_sna_atom>` - global array of bispectrum components on a regular grid
+* :doc:`sna/grid/local <compute_sna_atom>` - local array of bispectrum components on a regular grid
 * :doc:`snad/atom <compute_sna_atom>` - derivative of bispectrum components for each atom
 * :doc:`snav/atom <compute_sna_atom>` - virial contribution from bispectrum components for each atom
 * :doc:`sph/e/atom <compute_sph_e_atom>` - per-atom internal energy of Smooth-Particle Hydrodynamics atoms

doc/src/compute_born_matrix.rst

@@ -33,6 +33,8 @@ Examples
 Description
 """""""""""
 
+.. versionadded:: 4May2022
+
 Define a compute that calculates
 :math:`\frac{\partial{}^2U}{\partial\varepsilon_{i}\partial\varepsilon_{j}}` the
 second derivatives of the potential energy :math:`U` w.r.t. strain

doc/src/compute_contact_atom.rst

@@ -8,10 +8,11 @@ Syntax
 
 .. parsed-literal::
 
-   compute ID group-ID contact/atom
+   compute ID group-ID contact/atom group2-ID
 
 * ID, group-ID are documented in :doc:`compute <compute>` command
 * contact/atom = style name of this compute command
+* group2-ID = optional argument to restrict which atoms to consider for contacts (see below)
 
 Examples
 """"""""
@@ -19,6 +20,7 @@ Examples
 .. code-block:: LAMMPS
 
    compute 1 all contact/atom
+   compute 1 all contact/atom mygroup
 
 Description
 """""""""""
@@ -45,6 +47,9 @@ overview of LAMMPS output options.
 The per-atom vector values will be a number >= 0.0, as explained
 above.
 
+The optional *group2-ID* argument allows to specify from which group atoms
+contribute to the coordination number. Default setting is group 'all'.
+
 Restrictions
 """"""""""""
 
@@ -63,4 +68,7 @@ Related commands
 Default
 """""""
 
+*group2-ID* = all
+
+
 none

doc/src/compute_sna_atom.rst

@@ -2,6 +2,8 @@
 .. index:: compute snad/atom
 .. index:: compute snav/atom
 .. index:: compute snap
+.. index:: compute sna/grid
+.. index:: compute sna/grid/local
 
 compute sna/atom command
 ========================
@@ -15,6 +17,12 @@ compute snav/atom command
 compute snap command
 ====================
 
+compute sna/grid command
+========================
+
+compute sna/grid/local command
+==============================
+
 Syntax
 """"""
 
@@ -24,6 +32,9 @@ Syntax
    compute ID group-ID snad/atom rcutfac rfac0 twojmax R_1 R_2 ... w_1 w_2 ... keyword values ...
    compute ID group-ID snav/atom rcutfac rfac0 twojmax R_1 R_2 ... w_1 w_2 ... keyword values ...
    compute ID group-ID snap rcutfac rfac0 twojmax R_1 R_2 ... w_1 w_2 ... keyword values ...
+   compute ID group-ID snap rcutfac rfac0 twojmax R_1 R_2 ... w_1 w_2 ... keyword values ...
+   compute ID group-ID sna/grid nx ny nz rcutfac rfac0 twojmax R_1 R_2 ... w_1 w_2 ... keyword values ...
+   compute ID group-ID sna/grid/local nx ny nz rcutfac rfac0 twojmax R_1 R_2 ... w_1 w_2 ... keyword values ...
 
 * ID, group-ID are documented in :doc:`compute <compute>` command
 * sna/atom = style name of this compute command
@@ -32,8 +43,9 @@ Syntax
 * twojmax = band limit for bispectrum components (non-negative integer)
 * R_1, R_2,... = list of cutoff radii, one for each type (distance units)
 * w_1, w_2,... = list of neighbor weights, one for each type
+* nx, ny, nz = number of grid points in x, y, and z directions (positive integer)
 * zero or more keyword/value pairs may be appended
-* keyword = *rmin0* or *switchflag* or *bzeroflag* or *quadraticflag* or *chem* or *bnormflag* or *wselfallflag* or *bikflag* or *switchinnerflag* or *sinner* or *dinner*
+* keyword = *rmin0* or *switchflag* or *bzeroflag* or *quadraticflag* or *chem* or *bnormflag* or *wselfallflag* or *bikflag* or *switchinnerflag* or *sinner* or *dinner* or *dgradflag*
 
   .. parsed-literal::
 
@@ -56,9 +68,6 @@ Syntax
        *wselfallflag* value = *0* or *1*
           *0* = self-contribution only for element of central atom
           *1* = self-contribution for all elements
-       *bikflag* value = *0* or *1* (only implemented for compute snap)
-          *0* = per-atom bispectrum descriptors are summed over atoms
-          *1* = per-atom bispectrum descriptors are not summed over atoms
        *switchinnerflag* value = *0* or *1*
           *0* = do not use inner switching function
           *1* = use inner switching function
@@ -66,6 +75,12 @@ Syntax
           *sinnerlist* = *ntypes* values of *Sinner* (distance units)
        *dinner* values = *dinnerlist*
           *dinnerlist* = *ntypes* values of *Dinner* (distance units)
+       *bikflag* value = *0* or *1* (only implemented for compute snap)
+          *0* = descriptors are summed over atoms of each type
+          *1* = descriptors are listed separately for each atom
+       *dgradflag* value = *0* or *1* (only implemented for compute snap)
+          *0* = descriptor gradients are summed over atoms of each type
+          *1* = descriptor gradients are listed separately for each atom pair
 
 Examples
 """"""""
@@ -78,6 +93,7 @@ Examples
    compute snap all snap 1.4 0.95 6 2.0 1.0
    compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 chem 2 0 1
    compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 switchinnerflag 1 sinner 1.35 1.6 dinner 0.25 0.3
+   compute bgrid all sna/grid/local 200 200 200 1.4 0.95 6 2.0 1.0
 
 Description
 """""""""""
@@ -212,6 +228,46 @@ command:
 See section below on output for a detailed explanation of the data
 layout in the global array.
 
+The compute *sna/grid* and *sna/grid/local* commands calculate
+bispectrum components for a regular grid of points.
+These are calculated from the local density of nearby atoms *i'*
+around each grid point, as if there was a central atom *i*
+at the grid point. This is useful for characterizing fine-scale
+structure in a configuration of atoms, and it is used
+in the `MALA package <https://github.com/casus/mala>`_
+to build machine-learning surrogates for finite-temperature Kohn-Sham
+density functional theory (:ref:`Ellis et al. <Ellis2021>`)
+Neighbor atoms not in the group do not contribute to the
+bispectrum components of the grid points. The distance cutoff :math:`R_{ii'}`
+assumes that *i* has the same type as the neighbor atom *i'*.
+
+Compute *sna/grid* calculates a global array containing bispectrum
+components for a regular grid of points.
+The grid is aligned with the current box dimensions, with the
+first point at the box origin, and forming a regular 3d array with
+*nx*, *ny*, and *nz* points in the x, y, and z directions. For triclinic
+boxes, the array is congruent with the periodic lattice vectors
+a, b, and c. The array contains one row for each of the
+:math:`nx \times ny \times nz` grid points, looping over the index for *ix* fastest,
+then *iy*, and *iz* slowest.  Each row of the array contains the *x*, *y*,
+and *z* coordinates of the grid point, followed by the bispectrum
+components. See section below on output for a detailed explanation of the data
+layout in the global array.
+
+Compute *sna/grid/local* calculates bispectrum components of a regular
+grid of points similarly to compute *sna/grid* described above.
+However, because the array is local, it contains only rows for grid points
+that are local to the processor sub-domain. The global grid
+of :math:`nx \times ny \times nz` points is still laid out in space the same as for *sna/grid*,
+but grid points are strictly partitioned, so that every grid point appears in
+one and only one local array.  The array contains one row for each of the
+local grid points, looping over the global index *ix* fastest,
+then *iy*, and *iz* slowest.  Each row of the array contains
+the global indexes *ix*, *iy*, and *iz* first, followed by the *x*, *y*,
+and *z* coordinates of the grid point, followed by the bispectrum
+components. See section below on output for a detailed explanation of the data
+layout in the global array.
+
 The value of all bispectrum components will be zero for atoms not in
 the group. Neighbor atoms not in the group do not contribute to the
 bispectrum of atoms in the group.
@@ -307,15 +363,6 @@ This option is typically used in conjunction with the *chem* keyword,
 and LAMMPS will generate a warning if both *chem* and *bnormflag*
 are not both set or not both unset.
 
-The keyword *bikflag* determines whether or not to expand the bispectrum
-rows of the global array returned by compute snap. If *bikflag* is set
-to *1* then the bispectrum row, which is typically the per-atom bispectrum
-descriptors :math:`B_{i,k}` summed over all atoms *i* to produce
-:math:`B_k`, becomes bispectrum rows equal to the number of atoms. Thus,
-the resulting bispectrum rows are :math:`B_{i,k}` instead of just
-:math:`B_k`. In this case, the entries in the final column for these rows
-are set to zero.
-
 The keyword *switchinnerflag* with value 1
 activates an additional radial switching
 function similar to :math:`f_c(r)` above, but acting to switch off
@@ -340,6 +387,36 @@ When the central atom and the neighbor atom have different types,
 the values of :math:`S_{inner}` and :math:`D_{inner}` are
 the arithmetic means of the values for both types.
 
+The keywords *bikflag* and *dgradflag* are only used by compute *snap*.
+The keyword *bikflag* determines whether or not to list the descriptors
+of each atom separately, or sum them together and list in a single row.
+If *bikflag* is set
+to *0* then a single bispectrum row is used, which contains the per-atom bispectrum
+descriptors :math:`B_{i,k}` summed over all atoms *i* to produce
+:math:`B_k`.  If *bikflag* is set
+to *1* this is replaced by a separate per-atom bispectrum row for each atom.
+In this case, the entries in the final column for these rows
+are set to zero.
+
+The keyword *dgradflag* determines whether to sum atom gradients or list
+them separately. If *dgradflag* is set to 0, the bispectrum
+descriptor gradients w.r.t. atom *j* are summed over all atoms *i'*
+of type *I* (similar to *snad/atom* above).
+If *dgradflag* is set to 1, gradients are listed separately for each pair of atoms.
+Each row corresponds
+to a single term :math:`\frac{\partial {B_{i,k}  }}{\partial {r}^a_j}`
+where :math:`{r}^a_j` is the *a-th* position coordinate of the atom with global
+index *j*. This also changes
+the number of columns to be equal to the number of bispectrum components, with 3
+additional columns representing the indices :math:`i`, :math:`j`, and :math:`a`,
+as explained more in the Output info section below. The option *dgradflag=1*
+requires that *bikflag=1*.
+
+.. note::
+
+   Using *dgradflag* = 1 produces a global array with :math:`N + 3N^2 + 1` rows
+   which becomes expensive for systems with more than 1000 atoms.
+
 .. note::
 
    If you have a bonded system, then the settings of :doc:`special_bonds
@@ -414,6 +491,21 @@ number of columns in the global array generated by *snap* are 31, and
 931, respectively, while the number of rows is 1+3\*\ *N*\ +6, where *N*
 is the total number of atoms.
 
+Compute *sna/grid* evaluates a global array.
+The array contains one row for each of the
+:math:`nx \times ny \times nz` grid points, looping over the index for *ix* fastest,
+then *iy*, and *iz* slowest.  Each row of the array contains the *x*, *y*,
+and *z* coordinates of the grid point, followed by the bispectrum
+components.
+
+Compute *sna/grid/local* evaluates a local array.
+The array contains one row for each of the
+local grid points, looping over the global index *ix* fastest,
+then *iy*, and *iz* slowest.  Each row of the array contains
+the global indexes *ix*, *iy*, and *iz* first, followed by the *x*, *y*,
+and *z* coordinates of the grid point, followed by the bispectrum
+components.
+
 If the *quadratic* keyword value is set to 1, then additional columns
 are generated, corresponding to the products of all distinct pairs of
 bispectrum components. If the number of bispectrum components is *K*,
@@ -435,6 +527,42 @@ components. For the purposes of handling contributions to force, virial,
 and quadratic combinations, these :math:`N_{elem}^3` sub-blocks are
 treated as a single block of :math:`K N_{elem}^3` columns.
 
+If the *bik* keyword is set to 1, the structure of the snap array is expanded.
+The first :math:`N` rows of the snap array
+correspond to :math:`B_{i,k}` instead of a single row summed over atoms :math:`i`.
+In this case, the entries in the final column for these rows
+are set to zero. Also, each row contains only non-zero entries for the
+columns corresponding to the type of that atom. This is not true in the case
+of *dgradflag* keyword = 1 (see below).
+
+If the *dgradflag* keyword is set to 1, this changes the structure of the
+global array completely.
+Here the *snad/atom* quantities are replaced with rows corresponding to
+descriptor gradient components on single atoms:
+
+.. math::
+
+  \frac{\partial {B_{i,k}  }}{\partial {r}^a_j}
+
+where :math:`{r}^a_j` is the *a-th* position coordinate of the atom with global
+index *j*. The rows are
+organized in chunks, where each chunk corresponds to an atom with global index
+:math:`j`. The rows in an atom :math:`j` chunk correspond to
+atoms with global index :math:`i`. The total number of rows for
+these descriptor gradients is therefore :math:`3N^2`.
+The number of columns is equal to the number of bispectrum components,
+plus 3 additional left-most columns representing the global atom indices
+:math:`i`, :math:`j`,
+and Cartesian direction :math:`a`  (0, 1, 2, for x, y, z).
+The first 3 columns of the first :math:`N` rows belong to the reference
+potential force components. The remaining K columns contain the
+:math:`B_{i,k}` per-atom descriptors corresponding to the non-zero entries
+obtained when *bikflag* = 1.
+The first column of the last row, after the first
+:math:`N + 3N^2` rows, contains the reference potential
+energy. The virial components are not used with this option. The total number of
+rows is therefore :math:`N + 3N^2 + 1` and the number of columns is :math:`K + 3`.
+
 These values can be accessed by any command that uses per-atom values
 from a compute as input.  See the :doc:`Howto output <Howto_output>` doc
 page for an overview of LAMMPS output options. To see how this command
@@ -464,8 +592,7 @@ The optional keyword defaults are *rmin0* = 0,
 
 .. _Thompson20141:
 
-**(Thompson)** Thompson, Swiler, Trott, Foiles, Tucker, under review, preprint
-available at `arXiv:1409.3880 <http://arxiv.org/abs/1409.3880>`_
+**(Thompson)** Thompson, Swiler, Trott, Foiles, Tucker, J Comp Phys, 285, 316, (2015).
 
 .. _Bartok20101:
 
@@ -486,4 +613,8 @@ of Angular Momentum, World Scientific, Singapore (1987).
 
 .. _Cusentino2020:
 
-**(Cusentino)** Cusentino, Wood, and Thompson, J Phys Chem A, xxx, xxxxx, (2020)
+**(Cusentino)** Cusentino, Wood, Thompson, J Phys Chem A, 124, 5456, (2020)
+
+.. _Ellis2021:
+
+**(Ellis)** Ellis, Fiedler, Popoola, Modine, Stephens, Thompson, Cangi, Rajamanickam,  Phys Rev B, 104, 035120, (2021)

doc/src/create_atoms.rst

@@ -136,6 +136,8 @@ positions.
             :align: right
             :target: _images/marble_race.jpg
 
+.. versionadded:: 2Jun2022
+
 For the *mesh* style, a file with a triangle mesh in `STL format
 <https://en.wikipedia.org/wiki/STL_(file_format)>`_ is read and one or
 more particles are placed into the area of each triangle.  The reader
@@ -391,6 +393,8 @@ the atom style.  Its value is a prefactor (must be > 0.0, default is
 individual triangles in the triangle mesh that the particle corresponds
 to.
 
+.. versionadded:: 2Jun2022
+
 The *overlap* keyword only applies to the *random* style.  It prevents
 newly created particles from being created closer than the specified
 *Doverlap* distance from any other particle.  When the particles being

doc/src/dump.rst

@@ -1,4 +1,26 @@
 .. index:: dump
+.. index:: dump atom
+.. index:: dump cfg
+.. index:: dump custom
+.. index:: dump dcd
+.. index:: dump local
+.. index:: dump xtc
+.. index:: dump yaml
+.. index:: dump xyz
+.. index:: dump atom/gz
+.. index:: dump cfg/gz
+.. index:: dump custom/gz
+.. index:: dump local/gz
+.. index:: dump xyz/gz
+.. index:: dump atom/mpiio
+.. index:: dump cfg/mpiio
+.. index:: dump custom/mpiio
+.. index:: dump xyz/mpiio
+.. index:: dump atom/zstd
+.. index:: dump cfg/zstd
+.. index:: dump custom/zstd
+.. index:: dump xyz/zstd
+.. index:: dump local/zstd
 
 dump command
 ============
@@ -27,6 +49,9 @@ dump command
 :doc:`dump custom/adios <dump_adios>` command
 =============================================
 
+:doc:`dump cfg/uef <dump_cfg_uef>` command
+==========================================
+
 Syntax
 """"""
 
@@ -36,7 +61,7 @@ Syntax
 
 * ID = user-assigned name for the dump
 * group-ID = ID of the group of atoms to be dumped
-* style = *atom* or *atom/gz* or *atom/zstd or *atom/mpiio* or *cfg* or *cfg/gz* or *cfg/zstd* or *cfg/mpiio* or *custom* or *custom/gz* or *custom/zstd* or *custom/mpiio* or *dcd* or *h5md* or *image* or *local* or *local/gz* or *local/zstd* or *molfile* or *movie* or *netcdf* or *netcdf/mpiio* or *vtk* or *xtc* or *xyz* or *xyz/gz* or *xyz/zstd* or *xyz/mpiio* or *yaml*
+* style = *atom* or *atom/gz* or *atom/zstd or *atom/mpiio* or *cfg* or *cfg/gz* or *cfg/zstd* or *cfg/mpiio* or *cfg/uef* or *custom* or *custom/gz* or *custom/zstd* or *custom/mpiio* or *dcd* or *h5md* or *image* or *local* or *local/gz* or *local/zstd* or *molfile* or *movie* or *netcdf* or *netcdf/mpiio* or *vtk* or *xtc* or *xyz* or *xyz/gz* or *xyz/zstd* or *xyz/mpiio* or *yaml*
 * N = dump every this many timesteps
 * file = name of file to write dump info to
 * args = list of arguments for a particular style
@@ -47,22 +72,23 @@ Syntax
        *atom/gz* args = none
        *atom/zstd* args = none
        *atom/mpiio* args = none
-       *atom/adios* args = none,  discussed on :doc:`dump atom/adios <dump_adios>` doc page
+       *atom/adios* args = none,  discussed on :doc:`dump atom/adios <dump_adios>` page
        *cfg* args = same as *custom* args, see below
        *cfg/gz* args = same as *custom* args, see below
        *cfg/zstd* args = same as *custom* args, see below
        *cfg/mpiio* args = same as *custom* args, see below
+       *cfg/uef* args = same as *custom* args, discussed on :doc:`dump cfg/uef <dump_cfg_uef>` page
        *custom*, *custom/gz*, *custom/zstd*, *custom/mpiio* args = see below
-       *custom/adios* args = same as *custom* args, discussed on :doc:`dump custom/adios <dump_adios>` doc page
+       *custom/adios* args = same as *custom* args, discussed on :doc:`dump custom/adios <dump_adios>` page
        *dcd* args = none
-       *h5md* args = discussed on :doc:`dump h5md <dump_h5md>` doc page
-       *image* args = discussed on :doc:`dump image <dump_image>` doc page
+       *h5md* args = discussed on :doc:`dump h5md <dump_h5md>` page
+       *image* args = discussed on :doc:`dump image <dump_image>` page
        *local*, *local/gz*, *local/zstd* args = see below
-       *molfile* args = discussed on :doc:`dump molfile <dump_molfile>` doc page
-       *movie* args = discussed on :doc:`dump image <dump_image>` doc page
-       *netcdf* args = discussed on :doc:`dump netcdf <dump_netcdf>` doc page
-       *netcdf/mpiio* args = discussed on :doc:`dump netcdf <dump_netcdf>` doc page
-       *vtk* args = same as *custom* args, see below, also :doc:`dump vtk <dump_vtk>` doc page
+       *molfile* args = discussed on :doc:`dump molfile <dump_molfile>` page
+       *movie* args = discussed on :doc:`dump image <dump_image>` page
+       *netcdf* args = discussed on :doc:`dump netcdf <dump_netcdf>` page
+       *netcdf/mpiio* args = discussed on :doc:`dump netcdf <dump_netcdf>` page
+       *vtk* args = same as *custom* args, see below, also :doc:`dump vtk <dump_vtk>` page
        *xtc* args = none
        *xyz* args = none
        *xyz/gz* args = none
@@ -155,7 +181,7 @@ timesteps in one of several styles.  The *image* and *movie* styles are
 the exception: the *image* style renders a JPG, PNG, or PPM image file
 of the atom configuration every N timesteps while the *movie* style
 combines and compresses them into a movie file; both are discussed in
-detail on the :doc:`dump image <dump_image>` doc page.  The timesteps on
+detail on the :doc:`dump image <dump_image>` page.  The timesteps on
 which dump output is written can also be controlled by a variable.
 See the :doc:`dump_modify every <dump_modify>` command.
 
@@ -194,7 +220,7 @@ or multiple smaller files).
 For the *atom*, *custom*, *cfg*, and *local* styles, sorting is off by
 default.  For the *dcd*, *xtc*, *xyz*, and *molfile* styles, sorting
 by atom ID is on by default. See the :doc:`dump_modify <dump_modify>`
-doc page for details.
+page for details.
 
 The *atom/gz*, *cfg/gz*, *custom/gz*, *local/gz*, and *xyz/gz* styles
 are identical in command syntax to the corresponding styles without
@@ -204,7 +230,7 @@ alternative approach to writing compressed files via a pipe, as done
 by the regular dump styles, which may be required on clusters where
 the interface to the high-speed network disallows using the fork()
 library call (which is needed for a pipe).  For the remainder of this
-doc page, you should thus consider the *atom* and *atom/gz* styles
+page, you should thus consider the *atom* and *atom/gz* styles
 (etc) to be inter-changeable, with the exception of the required
 filename suffix.
 
@@ -218,7 +244,7 @@ As explained below, the *atom/mpiio*, *cfg/mpiio*, *custom/mpiio*, and
 *xyz/mpiio* styles are identical in command syntax and in the format
 of the dump files they create, to the corresponding styles without
 "mpiio", except the single dump file they produce is written in
-parallel via the MPI-IO library.  For the remainder of this doc page,
+parallel via the MPI-IO library.  For the remainder of this page,
 you should thus consider the *atom* and *atom/mpiio* styles (etc) to
 be inter-changeable.  The one exception is how the filename is
 specified for the MPI-IO styles, as explained below.
@@ -388,6 +414,8 @@ from using the (numerical) atom type to an element name (or some
 other label). This will help many visualization programs to guess
 bonds and colors.
 
+.. versionadded:: 4May2022
+
 Dump style *yaml* has the same command syntax as style *custom* and
 writes YAML format files that can be easily parsed by a variety of data
 processing tools and programming languages.  Each timestep will be
@@ -664,7 +692,7 @@ so that each value is 0.0 to 1.0.  If the simulation box is triclinic
 (tilted), then all atom coords will still be between 0.0 and 1.0.
 I.e. actual unscaled (x,y,z) = xs\*A + ys\*B + zs\*C, where (A,B,C) are
 the non-orthogonal vectors of the simulation box edges, as discussed
-on the :doc:`Howto triclinic <Howto_triclinic>` doc page.
+on the :doc:`Howto triclinic <Howto_triclinic>` page.
 
 Use *xu*, *yu*, *zu* if you want the coordinates "unwrapped" by the
 image flags for each atom.  Unwrapped means that if the atom has
@@ -779,6 +807,11 @@ To write gzipped dump files, you must either compile LAMMPS with the
 -DLAMMPS_GZIP option or use the styles from the COMPRESS package.
 See the :doc:`Build settings <Build_settings>` page for details.
 
+While a dump command is active (i.e. has not been stopped by using
+the undump command), no commands may be used that will change the
+timestep (e.g. :doc:`reset_timestep <reset_timestep>`).  LAMMPS
+will terminate with an error otherwise.
+
 The *atom/gz*, *cfg/gz*, *custom/gz*, and *xyz/gz* styles are part of
 the COMPRESS package.  They are only enabled if LAMMPS was built with
 that package.  See the :doc:`Build package <Build_package>` page for
@@ -787,7 +820,7 @@ more info.
 The *atom/mpiio*, *cfg/mpiio*, *custom/mpiio*, and *xyz/mpiio* styles
 are part of the MPIIO package.  They are only enabled if LAMMPS was
 built with that package.  See the :doc:`Build package <Build_package>`
-doc page for more info.
+page for more info.
 
 The *xtc*, *dcd* and *yaml* styles are part of the EXTRA-DUMP package.
 They are only enabled if LAMMPS was built with that package.  See the
@@ -797,12 +830,12 @@ Related commands
 """"""""""""""""
 
 :doc:`dump atom/adios <dump_adios>`, :doc:`dump custom/adios <dump_adios>`,
-:doc:`dump h5md <dump_h5md>`, :doc:`dump image <dump_image>`,
-:doc:`dump molfile <dump_molfile>`, :doc:`dump_modify <dump_modify>`,
-:doc:`undump <undump>`
+:doc:`dump cfg/uef <dump_cfg_uef>`, :doc:`dump h5md <dump_h5md>`,
+:doc:`dump image <dump_image>`, :doc:`dump molfile <dump_molfile>`,
+:doc:`dump_modify <dump_modify>`, :doc:`undump <undump>`, :doc:`write_dump <write_dump>`
 
 Default
 """""""
 
 The defaults for the *image* and *movie* styles are listed on the
-:doc:`dump image <dump_image>` doc page.
+:doc:`dump image <dump_image>` page.

doc/src/dump_image.rst

@@ -1,4 +1,5 @@
 .. index:: dump image
+.. index:: dump movie
 
 dump image command
 ==================

doc/src/dump_modify.rst

@@ -380,6 +380,8 @@ performed with dump style *xtc*\ .
 
 ----------
 
+.. versionadded:: 4May2022
+
 The *colname* keyword can be used to change the default header keyword
 for dump styles: *atom*, *custom*, and *cfg* and their compressed, ADIOS,
 and MPIIO variants.  The setting for *ID string* replaces the default

doc/src/dump_netcdf.rst

@@ -1,4 +1,5 @@
 .. index:: dump netcdf
+.. index:: dump netcdf/mpiio
 
 dump netcdf command
 ===================

doc/src/dumps.rst

@@ -0,0 +1,8 @@
+Dump Styles
+###############
+
+.. toctree::
+   :maxdepth: 1
+   :glob:
+
+   dump*

doc/src/fix.rst

@@ -171,6 +171,8 @@ accelerated styles exist.
 * :doc:`adapt/fep <fix_adapt_fep>` - enhanced version of fix adapt
 * :doc:`addforce <fix_addforce>` - add a force to each atom
 * :doc:`addtorque <fix_addtorque>` - add a torque to a group of atoms
+* :doc:`amoeba/bitorsion <fix_amoeba_bitorsion>` - torsion/torsion terms in AMOEBA force field
+* :doc:`amoeba/pitorsion <fix_amoeba_pitorsion>` - 6-body terms in AMOEBA force field
 * :doc:`append/atoms <fix_append_atoms>` - append atoms to a running simulation
 * :doc:`atc <fix_atc>` - initiates a coupled MD/FE simulation
 * :doc:`atom/swap <fix_atom_swap>` - Monte Carlo atom type swapping
@@ -194,7 +196,7 @@ accelerated styles exist.
 * :doc:`bond/swap <fix_bond_swap>` - Monte Carlo bond swapping
 * :doc:`box/relax <fix_box_relax>` - relax box size during energy minimization
 * :doc:`charge/regulation <fix_charge_regulation>` - Monte Carlo sampling of charge regulation
-* :doc:`cmap <fix_cmap>` - enables CMAP cross-terms of the CHARMM force field
+* :doc:`cmap <fix_cmap>` - CMAP torsion/torsion terms in CHARMM force field
 * :doc:`colvars <fix_colvars>` - interface to the collective variables "Colvars" library
 * :doc:`controller <fix_controller>` - apply control loop feedback mechanism
 * :doc:`damping/cundall <fix_damping_cundall>` - Cundall non-viscous damping for granular simulations
@@ -246,7 +248,7 @@ accelerated styles exist.
 * :doc:`lb/viscous <fix_lb_viscous>` -
 * :doc:`lineforce <fix_lineforce>` - constrain atoms to move in a line
 * :doc:`manifoldforce <fix_manifoldforce>` - restrain atoms to a manifold during minimization
-* :doc:`mdi/aimd <fix_mdi_aimd>` - LAMMPS operates as driver for ab initio MD (AIMD) via the MolSSI Driver Interface (MDI)
+* :doc:`mdi/qm <fix_mdi_qm>` - LAMMPS operates as driver for a quantum code via the MolSSI Driver Interface (MDI)
 * :doc:`meso/move <fix_meso_move>` - move mesoscopic SPH/SDPD particles in a prescribed fashion
 * :doc:`mol/swap <fix_mol_swap>` - Monte Carlo atom type swapping with a molecule
 * :doc:`momentum <fix_momentum>` - zero the linear and/or angular momentum of a group of atoms

doc/src/fix_amoeba_bitorsion.rst

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

doc/src/fix_amoeba_pitorsion.rst

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

doc/src/fix_ave_time.rst

@@ -277,6 +277,8 @@ is the length of the input vectors, and the number of columns is the
 number of values.  Thus the file ends up to be a series of these array
 sections.
 
+.. versionadded:: 4May2022
+
 If the filename ends in '.yaml' or '.yml' then the output format
 conforms to the `YAML standard <https://yaml.org/>`_ which allows
 easy import that data into tools and scripts that support reading YAML
@@ -329,6 +331,8 @@ appropriate fields from the fix ave/time command.
 Restart, fix_modify, output, run start/stop, minimize info
 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
 
+.. versionadded:: 4May2022
+
 No information about this fix is written to :doc:`binary restart files
 <restart>`.  The :doc:`fix_modify colname <fix_modify>` option can be
 used to change the name of the column in the output file.  When writing

doc/src/fix_cmap.rst

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

doc/src/fix_latte.rst

@@ -116,13 +116,6 @@ potential energy of the system as part of :doc:`thermodynamic output
 <thermo_style>`.  The default setting for this fix is :doc:`fix_modify
 energy yes <fix_modify>`.
 
-The :doc:`fix_modify <fix_modify>` *virial* option is supported by
-this fix to add the contribution compute by LATTE to the global
-pressure of the system via the :doc:`compute pressure
-<compute_pressure>` command.  This can be accessed by
-:doc:`thermodynamic output <thermo_style>`.  The default setting for
-this fix is :doc:`fix_modify virial yes <fix_modify>`.
-
 The :doc:`fix_modify <fix_modify>` *virial* option is supported by
 this fix to add the contribution computed by LATTE to the global
 pressure of the system as part of :doc:`thermodynamic output
@@ -137,7 +130,7 @@ energy discussed above.  The scalar value calculated by this fix is
 No parameter of this fix can be used with the *start/stop* keywords of
 the :doc:`run <run>` command.
 
-The DFTB forces computed by LATTE via this fix are imposed during an
+The DFTB forces computed by LATTE via this fix are used during an
 energy minimization, invoked by the :doc:`minimize <minimize>`
 command.
 

doc/src/fix_mdi_aimd.rst

@@ -1,93 +0,0 @@
-.. index:: fix mdi/aimd
-
-fix mdi/aimd command
-======================
-
-Syntax
-""""""
-
-.. parsed-literal::
-
-   fix ID group-ID mdi/aimd keyword
-
-* ID, group-ID are documented in :doc:`fix <fix>` command
-* mdi/aimd = style name of this fix command
-
-Examples
-""""""""
-
-.. code-block:: LAMMPS
-
-   fix 1 all mdi/aimd
-
-Description
-"""""""""""
-
-This command enables LAMMPS to act as a client with another server
-code to couple the two codes together to perform ab initio MD (AIMD)
-simulations.
-
-More specifically, this command causes LAMMPS to begin using the `MDI
-Library <https://molssi-mdi.github.io/MDI_Library/html/index.html>`_
-to run as an MDI driver (client), which sends MDI commands to an
-external MDI engine code (server) which in the case of AIMD is a
-quantum mechanics (QM) code, or could be LAMMPS itself, acting as a
-surrogate for a QM code.  See the :doc:`Howto mdi <Howto_mdi>` page
-for more information about how LAMMPS can operate as either an MDI
-driver or engine.
-
-The examples/mdi directory contains input scripts performing AIMD in
-this manner with LAMMPS acting as both a driver and an engine
-(surrogate for a QM code).  The examples/mdi/README file explains how
-to launch both driver and engine codes so that they communicate using
-the MDI library via either MPI or sockets.  Any QM code that supports
-MDI could be used in place of LAMMPS acting as a QM surrogate.  See
-the :doc:`Howto mdi <Howto_mdi>` page for a current list (March 2022)
-of such QM codes.
-
-The engine code can run either as a stand-alone code, launched at the
-same time as LAMMPS, or as a plugin library.  See the :doc:`mdi plugin
-<mdi>` command for how to trigger LAMMPS to load the plugin library.
-Again, the examples/mdi/README file explains how to launch both driver
-and engine codes so that engine is used in plugin mode.
-
-----------
-
-This fix performs the timestepping portion of an AIMD simulation.
-Both LAMMPS and the engine code (QM or LAMMPS) should define the same
-system (simulation box, atoms and their types) in their respective
-input scripts.  LAMMPS then begins its timestepping.
-
-At the point in each timestep when LAMMPS needs the force on each
-atom, it communicates with the engine code.  It sends the current
-simulation box size and shape (if they change dynamically, e.g. during
-an NPT simulation), and the current atom coordinates.  The engine code
-computes quantum forces on each atom and returns them to LAMMPS.  If
-LAMMPS also needs the system energy and/or virial, it requests those
-values from the engine code as well.
-
-Restrictions
-""""""""""""
-
-This command is part of the MDI package.  It is only enabled if
-LAMMPS was built with that package.  See the :doc:`Build package
-<Build_package>` page for more info.
-
-To use LAMMPS as an MDI driver in conjunction with other MDI-enabled
-atomistic codes, the :doc:`units <units>` command should be used to
-specify *real* or *metal* units.  This will ensure the correct unit
-conversions between LAMMPS and MDI units, which the other codes will
-also perform in their preferred units.
-
-LAMMPS can also be used as an MDI driver in other unit choices it
-supports, e.g. *lj*, but then no unit conversion is performed.
-
-Related commands
-""""""""""""""""
-
-:doc:`mdi engine <mdi>`
-
-Default
-"""""""
-
-none

doc/src/fix_mdi_qm.rst

@@ -0,0 +1,276 @@
+.. index:: fix mdi/qm
+
+fix mdi/qm command
+======================
+
+Syntax
+""""""
+
+.. parsed-literal::
+
+   fix ID group-ID mdi/qm keyword
+
+* ID, group-ID are documented in :doc:`fix <fix>` command
+* mdi/qm = style name of this fix command
+* zero or more keyword/value pairs may be appended
+* keyword = *virial* or *add* or *every* or *connect* or *elements*
+
+  .. parsed-literal::
+
+       *virial* args = *yes* or *no*
+         yes = request virial tensor from server code
+         no = do not request virial tensor from server code
+       *add* args = *yes* or *no*
+         yes = add returned value from server code to LAMMPS quantities
+         no = do not add returned values to LAMMPS quantities
+       *every* args = Nevery
+         Nevery = request values from server code once every Nevery steps
+       *connect* args = *yes* or *no*
+         yes = perform a one-time connection to the MDI engine code
+         no = do not perform the connection operation
+       *elements* args = N_1 N_2 ... N_ntypes
+         N_1,N_2,...N_ntypes = atomic number for each of ntypes LAMMPS atom types
+
+Examples
+""""""""
+
+.. code-block:: LAMMPS
+
+   fix 1 all mdi/qm
+   fix 1 all mdi/qm virial yes
+   fix 1 all mdi/qm add no every 100 elements 13 29
+
+Description
+"""""""""""
+
+This command enables LAMMPS to act as a client with another server
+code that will compute the total energy, per-atom forces, and total
+virial for atom conformations and simulation box size/shapes that
+LAMMPS sends it.
+
+Typically the server code will be a quantum mechanics (QM) code, hence
+the name of the fix.  However this is not required, the server code
+could be another classical molecular dynamics code or LAMMPS itself.
+The server code must support use of the `MDI Library
+<https://molssi-mdi.github.io/MDI_Library/html/index.html>`_ as
+explained below.
+
+These are example use cases for this fix, discussed further below:
+
+* perform an ab initio MD (AIMD) simulation with quantum forces
+* perform an energy minimization with quantum forces
+* perform a nudged elastic band (NEB) calculation with quantum forces
+* perform a QM calculation for a series of independent systems which LAMMPS reads or generates
+
+The code coupling performed by this command is done via the `MDI
+Library <https://molssi-mdi.github.io/MDI_Library/html/index.html>`_.
+LAMMPS runs as an MDI driver (client), and sends MDI commands to an
+external MDI engine code (server), e.g. a QM code which has support
+for MDI.  See the :doc:`Howto mdi <Howto_mdi>` page for more
+information about how LAMMPS can operate as either an MDI driver or
+engine.
+
+The examples/mdi directory contains input scripts using this fix in
+the various use cases discussed below.  In each case, two instances of
+LAMMPS are used, once as an MDI driver, once as an MDI engine
+(surrogate for a QM code).  The examples/mdi/README file explains how
+to launch two codes so that they communicate via the MDI library using
+either MPI or sockets.  Any QM code that supports MDI could be used in
+place of LAMMPS acting as a QM surrogate.  See the :doc:`Howto mdi
+<Howto_mdi>` page for a current list (March 2022) of such QM codes.
+
+Note that an engine code can support MDI in either or both of two
+modes.  It can be used as a stand-alone code, launched at the same
+time as LAMMPS.  Or it can be used as a plugin library, which LAMMPS
+loads.  See the :doc:`mdi plugin <mdi>` command for how to trigger
+LAMMPS to load a plugin library.  The examples/mdi/README file
+explains how to launch the two codes in either mode.
+
+----------
+
+The *virial* keyword setting of yes or no determines whether
+LAMMPS will request the QM code to also compute and return
+a 6-element symmetric virial tensor for the system.
+
+The *add* keyword setting of *yes* or *no* determines whether the
+energy and forces and virial returned by the QM code will be added to
+the LAMMPS internal energy and forces and virial or not.  If the
+setting is *no* then the default :doc:`fix_modify energy <fix_modify>`
+and :doc:`fix_modify virial <fix_modify>` settings are also set to
+*no* and your input scripts should not set them to yes.  See more
+details on these fix_modify settings below.
+
+Whatever the setting for the *add* keyword, the QM energy, forces, and
+virial will be stored by the fix, so they can be accessed by other
+commands.  See details below.
+
+The *every* keyword determines how often the QM code will be invoked
+during a dynamics run with the current LAMMPS simulation box and
+configuration of atoms.  The QM code will be called once every
+*Nevery* timesteps.
+
+The *connect* keyword determines whether this fix performs a one-time
+connection to the QM code.  The default is *yes*.  The only time a
+*no* is needed is if this command is used multiple times in an input
+script.  E.g. if it used inside a loop which also uses the :doc:`clear
+<clear>` command to destroy the system (including any defined fixes).
+See the examples/mdi/in.series.driver script as an example of this,
+where LAMMPS is using the QM code to compute energy and forces for a
+series of system configurations.  In this use case *connect no*
+is used along with the :doc:`mdi connect and exit <mdi>` command
+to one-time initiate/terminate the connection outside the loop.
+
+The *elements* keyword allows specification of what element each
+LAMMPS atom type corresponds to.  This is specified by the atomic
+number of the element, e.g. 13 for Al.  An atomic number must be
+specified for each of the ntypes LAMMPS atom types.  Ntypes is
+typically specified via the create_box command or in the data file
+read by the read_data command.  If this keyword is not specified, then
+this fix will send the LAMMPS atom type for each atom to the MDI
+engine.  If both the LAMMPS driver and the MDI engine are initialized
+so that atom type values are consistent in both codes, then the
+*elements* keyword is not needed.  Otherwise the keyword can be used
+to insure the two codes are consistent in their definition of atomic
+species.
+
+----------
+
+The following 3 example use cases are illustrated in the examples/mdi
+directory.  See its README file for more details.
+
+(1) To run an ab initio MD (AIMD) dynamics simulation, or an energy
+minimization with QM forces, or a multi-replica NEB calculation, use
+*add yes* and *every 1* (the defaults).  This is so that every time
+LAMMPS needs energy and forces, the QM code will be invoked.
+
+Both LAMMPS and the QM code should define the same system (simulation
+box, atoms and their types) in their respective input scripts.  Note
+that on this scenario, it may not be necessary for LAMMPS to define a
+pair style or use a neighbor list.
+
+LAMMPS will then perform the timestepping or minimization iterations
+for the simulation.  At the point in each timestep or iteration when
+LAMMPS needs the force on each atom, it communicates with the engine
+code.  It sends the current simulation box size and shape (if they
+change dynamically, e.g. during an NPT simulation), and the current
+atom coordinates.  The engine code computes quantum forces on each
+atom and the total energy of the system and returns them to LAMMPS.
+
+Note that if the AIMD simulation is an NPT or NPH model, or the energy
+minimization includes :doc:`fix box relax <fix_box_relax>` to
+equilibrate the box size/shape, then LAMMPS computes a pressure.  This
+means the *virial* keyword should be set to *yes* so that the QM
+contribution to the pressure can be included.
+
+(2) To run dynamics with a LAMMPS interatomic potential, and evaluate
+the QM energy and forces once every 1000 steps, use *add no* and
+*every 1000*.  This could be useful for using an MD run to generate
+randomized configurations which are then passed to the QM code to
+produce training data for a machine learning potential.  A :doc:`dump
+custom <dump>` command could be invoked every 1000 steps to dump the
+atom coordinates and QM forces to a file.  Likewise the QM energy and
+virial could be output with the :doc:`thermo_style custom
+<thermo_style>` command.
+
+(3) To do a QM evaluation of energy and forces for a series of *N*
+independent systems (simulation box and atoms), use *add no* and
+*every 1*.  Write a LAMMPS input script which loops over the *N*
+systems.  See the :doc:`Howto multiple <Howto_multiple>` doc page for
+details on looping and removing old systems.  The series of systems
+could be initialized by reading them from data files with
+:doc:`read_data <read_data>` commands.  Or, for example, by using the
+:doc:`lattice <lattice>` , :doc:`create_atoms <create_atoms>`,
+:doc:`delete_atoms <delete_atoms>`, and/or :doc:`displace_atoms
+random <displace_atoms>` commands to generate a series of different
+systems.  At the end of the loop perform :doc:`run 0 <run>` and
+:doc:`write_dump <write_dump>` commands to invoke the QM code and
+output the QM energy and forces.  As in (2) this be useful to produce
+QM data for training a machine learning potential.
+
+----------
+
+Restart, fix_modify, output, run start/stop, minimize info
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+No information about this fix is written to :doc:`binary restart files
+<restart>`.
+
+The :doc:`fix_modify <fix_modify>` *energy* option is supported by
+this fix to add the potential energy computed by the QM code to the
+global potential energy of the system as part of :doc:`thermodynamic
+output <thermo_style>`.  The default setting for this fix is
+:doc:`fix_modify energy yes <fix_modify>`, unless the *add* keyword is
+set to *no*, in which case the default setting is *no*.
+
+The :doc:`fix_modify <fix_modify>` *virial* option is supported by
+this fix to add the contribution computed by the QM code to the global
+pressure of the system as part of :doc:`thermodynamic output
+<thermo_style>`.  The default setting for this fix is :doc:`fix_modify
+virial yes <fix_modify>`, unless the *add* keyword is set to *no*, in
+which case the default setting is *no*.
+
+This fix computes a global scalar which can be accessed by various
+:doc:`output commands <Howto_output>`.  The scalar is the energy
+returned by the QM code.  The scalar value calculated by this fix is
+"extensive".
+
+This fix also computes a global vector with of length 6 which contains
+the symmetric virial tensor values returned by the QM code.  It can
+likewise be accessed by various :doc:`output commands <Howto_output>`.
+
+The ordering of values in the symmetric virial tensor is as follows:
+vxx, vyy, vzz, vxy, vxz, vyz.  The values will be in pressure
+:doc:`units <units>`.
+
+This fix also computes a peratom array with 3 columns which contains
+the peratom forces returned by the QM code.  It can likewise be
+accessed by various :doc:`output commands <Howto_output>`.
+
+No parameter of this fix can be used with the *start/stop* keywords of
+the :doc:`run <run>` command.
+
+Assuming the *add* keyword is set to *yes* (the default), the forces
+computed by the QM code are used during an energy minimization,
+invoked by the :doc:`minimize <minimize>` command.
+
+.. note::
+
+   If you want the potential energy associated with the QM forces to
+   be included in the total potential energy of the system (the
+   quantity being minimized), you MUST not disable the
+   :doc:`fix_modify <fix_modify>` *energy* option for this fix, which
+   means the *add* keyword should also be set to *yes* (the default).
+
+
+Restrictions
+""""""""""""
+
+This command is part of the MDI package.  It is only enabled if
+LAMMPS was built with that package.  See the :doc:`Build package
+<Build_package>` page for more info.
+
+The QM code does not currently compute and return per-atom energy or
+per-atom virial contributions.  So they will not show up as part of
+the calculations performed by the :doc:`compute pe/atom
+<compute_pe_atom>` or :doc:`compute stress/atom <compute_stress_atom>`
+commands.
+
+To use LAMMPS as an MDI driver in conjunction with other MDI-enabled
+codes (MD or QM codes), the :doc:`units <units>` command should be
+used to specify *real* or *metal* units.  This will ensure the correct
+unit conversions between LAMMPS and MDI units.  The other code will
+also perform similar unit conversions into its preferred units.
+
+LAMMPS can also be used as an MDI driver in other unit choices it
+supports, e.g. *lj*, but then no unit conversion is performed.
+
+Related commands
+""""""""""""""""
+
+:doc:`mdi plugin <mdi>`, :doc:`mdi engine <mdi>`
+
+Default
+"""""""
+
+The default for the optional keywords are virial = no, add = yes,
+every = 1, connect = yes.

doc/src/fix_npt_sphere.rst

@@ -22,7 +22,7 @@ Syntax
 
      *disc* value = none = treat particles as 2d discs, not spheres
 
-* additional thermostat and barostat related keyword/value pairs from the :doc:`fix npt <fix_nh>` command can be appended
+* NOTE: additional thermostat and barostat and dipole related keyword/value pairs from the :doc:`fix npt <fix_nh>` command can be appended
 
 Examples
 """"""""
@@ -33,6 +33,7 @@ Examples
    fix 2 all npt/sphere temp 300.0 300.0 100.0 x 5.0 5.0 1000.0
    fix 2 all npt/sphere temp 300.0 300.0 100.0 x 5.0 5.0 1000.0 disc
    fix 2 all npt/sphere temp 300.0 300.0 100.0 x 5.0 5.0 1000.0 drag 0.2
+   fix 2 all npt/sphere temp 300.0 300.0 100.0 x 5.0 5.0 1000.0 update dipole
    fix 2 water npt/sphere temp 300.0 300.0 100.0 aniso 0.0 0.0 1000.0 dilate partial
 
 Description
@@ -61,8 +62,9 @@ The only difference between discs and spheres in this context is their
 moment of inertia, as used in the time integration.
 
 Additional parameters affecting the thermostat and barostat are
-specified by keywords and values documented with the :doc:`fix npt <fix_nh>` command.  See, for example, discussion of the *temp*,
-*iso*, *aniso*, and *dilate* keywords.
+specified by keywords and values documented with the :doc:`fix npt
+<fix_nh>` command.  See, for example, discussion of the *temp*, *iso*,
+*aniso*, and *dilate* keywords.
 
 The particles in the fix group are the only ones whose velocities and
 positions are updated by the velocity/position update portion of the
@@ -87,8 +89,10 @@ this, the fix creates its own computes of style "temp/sphere" and
    compute fix-ID_temp all temp/sphere
    compute fix-ID_press all pressure fix-ID_temp
 
-See the :doc:`compute temp/sphere <compute_temp_sphere>` and :doc:`compute pressure <compute_pressure>` commands for details.  Note that the
-IDs of the new computes are the fix-ID + underscore + "temp" or fix_ID
+See the :doc:`compute temp/sphere <compute_temp_sphere>` and
+:doc:`compute pressure <compute_pressure>` commands for details.  Note
+that the IDs of the new computes are the fix-ID + underscore + "temp"
+or fix_ID
 + underscore + "press", and the group for the new computes is "all"
 since pressure is computed for the entire system.
 
@@ -170,7 +174,9 @@ defined by the :doc:`dimension <dimension>` keyword.
 Related commands
 """"""""""""""""
 
-:doc:`fix npt <fix_nh>`, :doc:`fix nve_sphere <fix_nve_sphere>`, :doc:`fix nvt_sphere <fix_nvt_sphere>`, :doc:`fix npt_asphere <fix_npt_asphere>`, :doc:`fix_modify <fix_modify>`
+:doc:`fix npt <fix_nh>`, :doc:`fix nve_sphere <fix_nve_sphere>`,
+     :doc:`fix nvt_sphere <fix_nvt_sphere>`, :doc:`fix npt_asphere
+     <fix_npt_asphere>`, :doc:`fix_modify <fix_modify>`
 
 Default
 """""""

doc/src/fix_nve_sphere.rst

@@ -51,7 +51,8 @@ If the *update* keyword is used with the *dipole* value, then the
 orientation of the dipole moment of each particle is also updated
 during the time integration.  This option should be used for models
 where a dipole moment is assigned to finite-size particles,
-e.g. spheroids via use of the :doc:`atom_style hybrid sphere dipole <atom_style>` command.
+e.g. spheroids via use of the :doc:`atom_style hybrid sphere dipole
+<atom_style>` command.
 
 The default dipole orientation integrator can be changed to the
 Dullweber-Leimkuhler-McLachlan integration scheme
@@ -75,11 +76,13 @@ moment of inertia, as used in the time integration.
 Restart, fix_modify, output, run start/stop, minimize info
 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
 
-No information about this fix is written to :doc:`binary restart files <restart>`.  None of the :doc:`fix_modify <fix_modify>` options
-are relevant to this fix.  No global or per-atom quantities are stored
-by this fix for access by various :doc:`output commands <Howto_output>`.
+No information about this fix is written to :doc:`binary restart files
+<restart>`.  None of the :doc:`fix_modify <fix_modify>` options are
+relevant to this fix.  No global or per-atom quantities are stored by
+this fix for access by various :doc:`output commands <Howto_output>`.
 No parameter of this fix can be used with the *start/stop* keywords of
-the :doc:`run <run>` command.  This fix is not invoked during :doc:`energy minimization <minimize>`.
+the :doc:`run <run>` command.  This fix is not invoked during
+:doc:`energy minimization <minimize>`.
 
 Restrictions
 """"""""""""

doc/src/fix_nvt_sphere.rst

@@ -22,7 +22,7 @@ Syntax
 
        *disc* value = none = treat particles as 2d discs, not spheres
 
-* additional thermostat related keyword/value pairs from the :doc:`fix nvt <fix_nh>` command can be appended
+* NOTE: additional thermostat and dipole related keyword/value pairs from the :doc:`fix nvt <fix_nh>` command can be appended
 
 Examples
 """"""""
@@ -32,6 +32,7 @@ Examples
    fix 1 all nvt/sphere temp 300.0 300.0 100.0
    fix 1 all nvt/sphere temp 300.0 300.0 100.0 disc
    fix 1 all nvt/sphere temp 300.0 300.0 100.0 drag 0.2
+   fix 1 all nvt/sphere temp 300.0 300.0 100.0 update dipole
 
 Description
 """""""""""
@@ -77,13 +78,13 @@ underscore + "temp", and the group for the new compute is the same as
 the fix group.
 
 Note that this is NOT the compute used by thermodynamic output (see
-the :doc:`thermo_style <thermo_style>` command) with ID = *thermo_temp*.
-This means you can change the attributes of this fix's temperature
-(e.g. its degrees-of-freedom) via the
-:doc:`compute_modify <compute_modify>` command or print this temperature
-during thermodynamic output via the :doc:`thermo_style custom <thermo_style>` command using the appropriate compute-ID.
-It also means that changing attributes of *thermo_temp* will have no
-effect on this fix.
+the :doc:`thermo_style <thermo_style>` command) with ID =
+*thermo_temp*.  This means you can change the attributes of this fix's
+temperature (e.g. its degrees-of-freedom) via the :doc:`compute_modify
+<compute_modify>` command or print this temperature during
+thermodynamic output via the :doc:`thermo_style custom <thermo_style>`
+command using the appropriate compute-ID.  It also means that changing
+attributes of *thermo_temp* will have no effect on this fix.
 
 Like other fixes that perform thermostatting, this fix can be used
 with :doc:`compute commands <compute>` that remove a "bias" from the
@@ -148,7 +149,9 @@ defined by the :doc:`dimension <dimension>` keyword.
 Related commands
 """"""""""""""""
 
-:doc:`fix nvt <fix_nh>`, :doc:`fix nve_sphere <fix_nve_sphere>`, :doc:`fix nvt_asphere <fix_nvt_asphere>`, :doc:`fix npt_sphere <fix_npt_sphere>`, :doc:`fix_modify <fix_modify>`
+:doc:`fix nvt <fix_nh>`, :doc:`fix nve_sphere <fix_nve_sphere>`,
+     :doc:`fix nvt_asphere <fix_nvt_asphere>`, :doc:`fix npt_sphere
+     <fix_npt_sphere>`, :doc:`fix_modify <fix_modify>`
 
 Default
 """""""

doc/src/fix_reaxff_species.rst

@@ -20,7 +20,7 @@ Syntax
 * Nfreq = calculate average bond-order every this many timesteps
 * filename = name of output file
 * zero or more keyword/value pairs may be appended
-* keyword = *cutoff* or *element* or *position*
+* keyword = *cutoff* or *element* or *position* or *delete*
 
   .. parsed-literal::
 
@@ -31,6 +31,14 @@ Syntax
        *position* value = posfreq filepos
          posfreq = write position files every this many timestep
          filepos = name of position output file
+       *delete* value = filedel keyword value
+         filedel = name of delete species output file
+         keyword = *specieslist* or *masslimit*
+           *specieslist* value = Nspecies Species1 Species2 ...
+             Nspecies = number of species in list
+           *masslimit* value = massmin massmax
+             massmin = minimum molecular weight of species to delete
+             massmax = maximum molecular weight of species to delete
 
 Examples
 """"""""
@@ -40,6 +48,7 @@ Examples
    fix 1 all reaxff/species 10 10 100 species.out
    fix 1 all reaxff/species 1 2 20 species.out cutoff 1 1 0.40 cutoff 1 2 0.55
    fix 1 all reaxff/species 1 100 100 species.out element Au O H position 1000 AuOH.pos
+   fix 1 all reaxff/species 1 100 100 species.out delete species.del masslimit 0 50
 
 Description
 """""""""""
@@ -59,13 +68,18 @@ the first line.
 .. warning::
 
    In order to compute averaged data, it is required that there are no
-   neighbor list rebuilds between the *Nfreq* steps. For that reason, fix
-   *reaxff/species* may change your neighbor list settings.  There will
-   be a warning message showing the new settings. Having an *Nfreq*
-   setting that is larger than what is required for correct computation
-   of the ReaxFF force field interactions can thus lead to incorrect
-   results.  For typical ReaxFF calculations a value of 100 is already
-   quite large.
+   neighbor list rebuilds for at least Nrepeat\*Nevery steps preceding
+   each *Nfreq* step.  For that reason, fix *reaxff/species* may
+   change your neighbor list settings.  Reneighboring will occur no
+   more frequently than every Nrepeat\*Nevery timesteps, and will
+   occur less frequently if *Nfreq* is not a multiple of
+   Nrepeat\*Nevery.  There will be a warning message showing the new
+   settings.  Having a *Nfreq* setting that is larger than what is
+   required for correct computation of the ReaxFF force field
+   interactions, in combination with certain *Nrepeat* and *Nevery*
+   settings, can thus lead to incorrect results.  For typical ReaxFF
+   calculations, reneighboring only every 100 steps is already quite a
+   low frequency.
 
 If the filename ends with ".gz", the output file is written in gzipped
 format.  A gzipped dump file will be about 3x smaller than the text version,
@@ -104,6 +118,30 @@ character appears in *filepos*, then one file per snapshot is written
 at *posfreq* and the "\*" character is replaced with the timestep
 value.  For example, AuO.pos.\* becomes AuO.pos.0, AuO.pos.1000, etc.
 
+The optional keyword *delete* enables the periodic removal of
+molecules from the system.  Criteria for deletion can be either a list
+of specific chemical formulae or a range of molecular weights.
+Molecules are deleted every *Nfreq* timesteps, and bond connectivity
+is determined using the *Nevery* and *Nrepeat* keywords.  The
+*filedel* argument is the name of the output file that records the
+species that are removed from the system.  The *specieslist* keyword
+permits specific chemical species to be deleted.  The *Nspecies*
+argument specifies how many species are eligible for deletion and is
+followed by a list of chemical formulae, whose strings are compared to
+species identified by this fix.  For example, "specieslist 2 CO CO2"
+deletes molecules that are identified as "CO" and "CO2" in the species
+output file.  When using the *specieslist* keyword, the *filedel* file
+has the following format: the first line lists the chemical formulae
+eligible for deletion, and each additional line contains the timestep
+on which a molecule deletion occurs and the number of each species
+deleted on that timestep.  The *masslimit* keyword permits deletion of
+molecules with molecular weights between *massmin* and *massmax*.
+When using the *masslimit* keyword, each line of the *filedel* file
+contains the timestep on which deletions occurs, followed by how many
+of each species are deleted (with quantities preceding chemical
+formulae).  The *specieslist* and *masslimit* keywords cannot both be
+used in the same *reaxff/species* fix.
+
 ----------
 
 The *Nevery*, *Nrepeat*, and *Nfreq* arguments specify on what

doc/src/if.rst

@@ -156,27 +156,28 @@ and Boolean operators:
 Each A and B is a number or string or a variable reference like $a or
 ${abc}, or A or B can be another Boolean expression.
 
-If a variable is used it can produce a number when evaluated, like an
-:doc:`equal-style variable <variable>`.  Or it can produce a string,
-like an :doc:`index-style variable <variable>`.  For an individual
-Boolean operator, A and B must both be numbers or must both be
-strings.  You cannot compare a number to a string.
+Note that all variables used will be substituted for before the
+Boolean expression in evaluated.  A variable can produce a number,
+like an :doc:`equal-style variable <variable>`.  Or it can produce a
+string, like an :doc:`index-style variable <variable>`.
+
+The Boolean operators "==" and "!=" can operate on a pair or strings
+or numbers.  They cannot compare a number to a string.  All the other
+Boolean operations can only operate on numbers.
 
 Expressions are evaluated left to right and have the usual C-style
 precedence: the unary logical NOT operator "!" has the highest
 precedence, the 4 relational operators "<", "<=", ">", and ">=" are
 next; the two remaining relational operators "==" and "!=" are next;
 then the logical AND operator "&&"; and finally the logical OR
-operator "\|\|" and logical XOR (exclusive or) operator "\|\^" have the
-lowest precedence.  Parenthesis can be used to group one or more
+operator "\|\|" and logical XOR (exclusive or) operator "\|\^" have
+the lowest precedence.  Parenthesis can be used to group one or more
 portions of an expression and/or enforce a different order of
 evaluation than what would occur with the default precedence.
 
 When the 6 relational operators (first 6 in list above) compare 2
 numbers, they return either a 1.0 or 0.0 depending on whether the
-relationship between A and B is TRUE or FALSE.  When the 6 relational
-operators compare 2 strings, they also return a 1.0 or 0.0 for TRUE or
-FALSE, but the comparison is done by the C function strcmp().
+relationship between A and B is TRUE or FALSE.
 
 When the 3 logical operators (last 3 in list above) compare 2 numbers,
 they also return either a 1.0 or 0.0 depending on whether the
@@ -190,8 +191,16 @@ returns 1.0 if its argument is 0.0, else it returns 0.0.  The 3
 logical operators can only be used to operate on numbers, not on
 strings.
 
-The overall Boolean expression produces a TRUE result if the result is
-non-zero.  If the result is zero, the expression result is FALSE.
+The overall Boolean expression produces a TRUE result if the numeric
+result is non-zero.  If the result is zero, the expression result is
+FALSE.
+
+.. note::
+
+   If the Boolean expression is a single numeric value with no Boolean
+   operators, it will be FALSE if the value = 0.0, otherwise TRUE.  If
+   the Boolean expression is a single string, an error message will be
+   issued.
 
 ----------
 

doc/src/improper_amoeba.rst

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

doc/src/improper_style.rst

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

doc/src/kim_commands.rst

@@ -29,15 +29,16 @@ Examples
 Description
 """""""""""
 
-The *kim command* includes a set of sub-commands that allow LAMMPS users to use
-interatomic models (IM) (potentials and force fields) and their predictions for
-various physical properties archived in the
-`Open Knowledgebase of Interatomic Models (OpenKIM) <https://openkim.org>`_
-repository.
+The *kim command* includes a set of sub-commands that allow LAMMPS
+users to use interatomic models (IM) (potentials and force fields) and
+their predictions for various physical properties archived in the
+`Open Knowledgebase of Interatomic Models (OpenKIM)
+<https://openkim.org>`_ repository.
 
-Using OpenKIM provides LAMMPS users with immediate access to a large number of
-verified IMs and their predictions. OpenKIM IMs have multiple benefits including
-`reliability, reproducibility and convenience <https://openkim.org/doc/overview/kim-features/>`_.
+Using OpenKIM provides LAMMPS users with immediate access to a large
+number of verified IMs and their predictions. OpenKIM IMs have
+multiple benefits including `reliability, reproducibility and
+convenience <https://openkim.org/doc/overview/kim-features/>`_.
 
 .. _IM_types:
 

doc/src/label.rst

@@ -38,7 +38,7 @@ Restrictions
 Related commands
 """"""""""""""""
 
-none
+:doc:`jump <jump>`, :doc:`next <next>`
 
 
 Default

doc/src/mdi.rst

@@ -8,21 +8,26 @@ Syntax
 
 .. parsed-literal::
 
-   mdi mode args
+   mdi option args
 
-* mode = *engine* or *plugin*
+* option = *engine* or *plugin* or *connect* or *exit*
 
   .. parsed-literal::
 
-     *engine* args = none
-     *plugin* args = name keyword value keyword value
+     *engine* args = zero or more keyword arg pairs
+       keywords = *elements*
+         *elements* args = N_1 N_2 ... N_ntypes
+           N_1,N_2,...N_ntypes = atomic number for each of ntypes LAMMPS atom types
+     *plugin* args = name keyword value keyword value ...
        name = name of plugin library, e.g. lammps means a liblammps.so library will be loaded
+       keyword/value pairs in any order, some are required, some are optional
        keywords = *mdi* or *infile* or *extra* or *command*
-         *mdi* value = args passed to MDI for driver to operate with plugins
-         *infile* value = filename the engine will read at start-up
+         *mdi* value = args passed to MDI for driver to operate with plugins (required)
+         *infile* value = filename the engine will read at start-up (optional)
          *extra* value = aditional command-line args to pass to engine library when loaded
-         *command* value = a LAMMPS input script command to execute
-
+         *command* value = a LAMMPS input script command to execute (required)
+     *connect* args = none
+     *exit* args = none
 
 Examples
 """"""""
@@ -30,26 +35,19 @@ Examples
 .. code-block:: LAMMPS
 
    mdi engine
+   mdi engine elements 13 29
    mdi plugin lammps mdi "-role ENGINE -name lammps -method LINK" &
               infile in.aimd.engine extra "-log log.aimd.engine.plugin" &
               command "run 5"
+   mdi connect
+   mdi exit
 
 Description
 """""""""""
 
-This command implements two high-level operations within LAMMPS to use
-the `MDI Library
-<https://molssi-mdi.github.io/MDI_Library/html/index.html>` for
-coupling to other codes in a client/server protocol.
-
-The *engine* mode enables LAMMPS to act as an MDI engine (server),
-responding to requests from an MDI driver (client) code.
-
-The *plugin* mode enables LAMMPS to act as an MDI driver (client), and
-load the MDI engine (server) code as a library plugin.  In this case
-the MDI engine is a library plugin.  It can also be a stand-alone
-code, launched separately from LAMMPS, in which case the mdi plugin
-command is not used.
+This command implements operations within LAMMPS to use the `MDI
+Library <https://molssi-mdi.github.io/MDI_Library/html/index.html>`
+for coupling to other codes in a client/server protocol.
 
 See the Howto MDI doc page for a discussion of all the different ways
 2 or more codes can interact via MDI.
@@ -61,6 +59,22 @@ stand-alone code or as a plugin.  The README file in that directory
 shows how to launch and couple codes for all the 4 usage modes, and so
 they communicate via the MDI library using either MPI or sockets.
 
+The scripts in that directory illustrate the use of all the options
+for this command.
+
+The *engine* option enables LAMMPS to act as an MDI engine (server),
+responding to requests from an MDI driver (client) code.
+
+The *plugin* option enables LAMMPS to act as an MDI driver (client),
+and load the MDI engine (server) code as a library plugin.  In this
+case the MDI engine is a library plugin.  An MDI engine can also be a
+stand-alone code, launched separately from LAMMPS, in which case the
+mdi plugin command is not used.
+
+The *connect* and *exit* options are only used when LAMMPS is acting
+as an MDI driver.  As explained below, these options are normally not
+needed, except for a specific kind of use case.
+
 ----------
 
 The *mdi engine* command is used to make LAMMPS operate as an MDI
@@ -100,6 +114,8 @@ commands, which are described further below.
      - Send/request charge on each atom (N values)
    * - >COORDS or <COORDS
      - Send/request coordinates of each atom (3N values)
+   * - >ELEMENTS
+     - Send elements (atomic numbers) for each atom (N values)
    * - <ENERGY
      - Request total energy (potential + kinetic) of the system (1 value)
    * - >FORCES or <FORCES
@@ -121,11 +137,11 @@ commands, which are described further below.
    * - <PE
      - Request potential energy of the system (1 value)
    * - <STRESS
-     - Request stress tensor (virial) of the system (6 values)
+     - Request symmetric stress tensor (virial) of the system (9 values)
    * - >TOLERANCE
      - Send 4 tolerance parameters for next MD minimization via OPTG command
    * - >TYPES or <TYPES
-     - Send/request the numeric type of each atom (N values)
+     - Send/request the LAMMPS atom type for each atom (N values)
    * - >VELOCITIES or <VELOCITIES
      - Send/request the velocity of each atom (3N values)
    * - @INIT_MD or @INIT_OPTG
@@ -145,9 +161,25 @@ commands, which are described further below.
    builds.  If the change in atom positions is large (since the
    previous >COORDS command), then LAMMPS will do a more expensive
    operation to migrate atoms to new processors as needed and
-   re-neighbor.  If the >NATOMS or >TYPES commands have been sent
-   (since the previous >COORDS command), then LAMMPS assumes the
-   system is new and re-initializes an entirely new simulation.
+   re-neighbor.  If the >NATOMS or >TYPES or >ELEMENTS commands have
+   been sent (since the previous >COORDS command), then LAMMPS assumes
+   the system is new and re-initializes an entirely new simulation.
+
+.. note::
+
+   The >TYPES or >ELEMENTS commands are how the MDI driver tells the
+   LAMMPS engine which LAMMPS atom type to assign to each atom.  If
+   both the MDI driver and the LAMMPS engine are initialized so that
+   atom type values are consistent in both codes, then the >TYPES
+   command can be used.  If not, the optional *elements* keyword can
+   be used to specify what element each LAMMPS atom type corresponds
+   to.  This is specified by the atomic number of the element, e.g. 13
+   for Al.  An atomic number must be specified for each of the ntypes
+   LAMMPS atom types.  Ntypes is typically specified via the
+   create_box command or in the data file read by the read_data
+   command.  In this has been done, the MDI driver can send an
+   >ELEMENTS command to the LAMMPS driver with the atomic number of
+   each atom.
 
 The MD and OPTG commands perform an entire MD simulation or energy
 minimization (to convergence) with no communication from the driver
@@ -270,7 +302,7 @@ The *command* keyword is required.  It specifies a LAMMPS input script
 command (as a single argument in quotes if it is multiple words).
 Once the plugin library is launched, LAMMPS will execute this command.
 Other previously-defined commands in the input script, such as the
-:doc:`fix mdi/aimd <fix_mdi_aimd>` command, should perform MDI
+:doc:`fix mdi/qm <fix_mdi_qm>` command, should perform MDI
 communication with the engine, while the specified *command* executes.
 Note that if *command* is an :doc:`include <include>` command, then it
 could specify a filename with multiple LAMMPS commands.
@@ -284,6 +316,31 @@ could specify a filename with multiple LAMMPS commands.
    "mdi plugin" command could then load the same library plugin or
    a different one if desired.
 
+----------
+
+The *mdi connect* and *mdi exit* commands are only used when LAMMPS is
+operating as an MDI driver.  And when other LAMMPS command(s) which
+send MDI commands and associated data to/from the MDI engine are not
+able to initiate and terminate the connection to the engine code.
+
+The only current MDI driver command in LAMMPS is the :doc:`fix mdi/qm
+<fix_mdi_qm>` command.  If it is only used once in an input script
+then it can initiate and terminate the connection.  But if it is being
+issued multiple times, e.g. in a loop that issues a :doc:`clear
+<clear>` command, then it cannot initiate or terminate the connection
+multiple times.  Instead, the *mdi connect* and *mdi exit* commands
+should be used outside the loop to initiate or terminate the connection.
+
+See the examples/mdi/in.series.driver script for an example of how
+this is done.  The LOOP in that script is reading a series of data
+file configurations and passing them to an MDI engine (e.g. quantum
+code) for energy and force evaluation.  A *clear* command inside the
+loop wipes out the current system so a new one can be defined.  This
+operation also destroys all fixes.  So the :doc:`fix mdi/qm
+<fix_mdi_qm>` command is issued once per loop iteration.  Note that it
+includes a "connect no" option which disables the initiate/terminate
+logic within that fix.
+
 
 Restrictions
 """"""""""""
@@ -304,7 +361,7 @@ supports, e.g. *lj*, but then no unit conversion is performed.
 Related commands
 """"""""""""""""
 
-:doc:`fix mdi/aimd <fix_mdi_aimd>`
+:doc:`fix mdi/qm <fix_mdi_qm>`
 
 Default
 """""""

doc/src/neighbor.rst

@@ -49,29 +49,27 @@ sometimes be faster.  Either style should give the same answers.
 
 The *multi* style is a modified binning algorithm that is useful for
 systems with a wide range of cutoff distances, e.g. due to different
-size particles. For granular pair styles, cutoffs are set to the
-sum of the maximum atomic radii for each atom type.
-For the *bin* style, the bin size is set to 1/2 of
-the largest cutoff distance between any pair of atom types and a
-single set of bins is defined to search over for all atom types.  This
-can be inefficient if one pair of types has a very long cutoff, but
-other type pairs have a much shorter cutoff. The *multi* style uses
-different sized bins for collections of different sized particles, where
-"size" may mean the physical size of the particle or its cutoff
-distance for interacting with other particles. Different
+size particles. For granular pair styles, cutoffs are set to the sum of
+the maximum atomic radii for each atom type.  For the *bin* style, the
+bin size is set to 1/2 of the largest cutoff distance between any pair
+of atom types and a single set of bins is defined to search over for all
+atom types.  This can be inefficient if one pair of types has a very
+long cutoff, but other type pairs have a much shorter cutoff. The
+*multi* style uses different sized bins for collections of different
+sized particles, where "size" may mean the physical size of the particle
+or its cutoff distance for interacting with other particles. Different
 sets of bins are then used to construct the neighbor lists as as further
 described by Shire, Hanley, and Stratford :ref:`(Shire) <bytype-Shire>`.
-This imposes some extra setup overhead, but the searches themselves
-may be much faster. By default, each atom type defines a separate
-collection of particles. For systems where two or more atom types
-have the same size (either physical size or cutoff distance), the
-definition of collections can be customized, which can result in less
-overhead and faster performance. See the :doc:`neigh_modify <neigh_modify>`
-command for how to define custom collections. Whether the collection
-definition is customized or not, also see the
-:doc:`comm_modify mode multi <comm_modify>` command for communication
-options that further improve performance in a manner consistent with
-neighbor style multi.
+This imposes some extra setup overhead, but the searches themselves may
+be much faster. By default, each atom type defines a separate collection
+of particles. For systems where two or more atom types have the same
+size (either physical size or cutoff distance), the definition of
+collections can be customized, which can result in less overhead and
+faster performance. See the :doc:`neigh_modify <neigh_modify>` command
+for how to define custom collections. Whether the collection definition
+is customized or not, also see the :doc:`comm_modify mode multi
+<comm_modify>` command for communication options that further improve
+performance in a manner consistent with neighbor style multi.
 
 An alternate style, *multi/old*, sets the bin size to 1/2 of the shortest
 cutoff distance and multiple sets of bins are defined to search over for
@@ -80,6 +78,16 @@ algorithm in LAMMPS but was found to be significantly slower than the new
 approach. For now we are keeping the old option in case there are use cases
 where multi/old outperforms the new multi style.
 
+.. note::
+
+   If there are multiple sub-styles in a :doc:`hybrid/overlay pair style
+   <pair_hybrid>` that cover the same atom types, but have significantly
+   different cutoffs, the *multi* style does not apply.  Instead, the
+   :doc:`pair_modify neigh/trim <pair_modify>` setting applies (which is
+   *yes* by default).  Please check the neighbor list summary printed at
+   the beginning of a calculation to verify that the desired set of
+   neighbor list builds is performed.
+
 
 The :doc:`neigh_modify <neigh_modify>` command has additional options
 that control how often neighbor lists are built and which pairs are

doc/src/package.rst

@@ -71,7 +71,7 @@ Syntax
            *no_affinity* values = none
        *kokkos* args = keyword value ...
          zero or more keyword/value pairs may be appended
-         keywords = *neigh* or *neigh/qeq* or *neigh/thread* or *newton* or *binsize* or *comm* or *comm/exchange* or *comm/forward* *comm/pair/forward* *comm/fix/forward* or *comm/reverse* or *gpu/aware* or *pair/only*
+         keywords = *neigh* or *neigh/qeq* or *neigh/thread* or *newton* or *binsize* or *comm* or *comm/exchange* or *comm/forward* *comm/pair/forward* *comm/fix/forward* or *comm/reverse* or *comm/pair/reverse* or *gpu/aware* or *pair/only*
            *neigh* value = *full* or *half*
              full = full neighbor list
              half = half neighbor list built in thread-safe manner
@@ -96,6 +96,7 @@ Syntax
            *comm/pair/forward* value = *no* or *device*
            *comm/fix/forward* value = *no* or *device*
            *comm/reverse* value = *no* or *host* or *device*
+           *comm/pair/reverse* value = *no* or *device*
              no = perform communication pack/unpack in non-KOKKOS mode
              host = perform pack/unpack on host (e.g. with OpenMP threading)
              device = perform pack/unpack on device (e.g. on GPU)
@@ -500,7 +501,7 @@ rule of thumb may give too large a binsize and the default should be
 overridden with a smaller value.
 
 The *comm* and *comm/exchange* and *comm/forward* and *comm/pair/forward*
-and *comm/fix/forward* and comm/reverse*
+and *comm/fix/forward* and *comm/reverse* and *comm/pair/reverse*
 keywords determine whether the host or device performs the packing and
 unpacking of data when communicating per-atom data between processors.
 "Exchange" communication happens only on timesteps that neighbor lists
@@ -521,9 +522,16 @@ packing/unpacking data for the communication. A value of *host* means to
 use the host, typically a multi-core CPU, and perform the
 packing/unpacking in parallel with threads. A value of *device* means to
 use the device, typically a GPU, to perform the packing/unpacking
-operation. If a value of *host* is used for the *comm/pair/forward* or
-*comm/fix/forward* keyword, it will be automatically be changed to *no*
-since these keywords don't support *host* mode.
+operation.
+
+For the *comm/pair/forward* or *comm/fix/forward* or *comm/pair/reverse*
+keywords, if a value of *host* is used it will be automatically
+be changed to *no* since these keywords don't support *host* mode. The
+value of *no* will also always be used when running on the CPU, i.e. setting
+the value to *device* will have no effect if the pair/fix style is
+running on the CPU. For the *comm/fix/forward* or *comm/pair/reverse*
+keywords, not all styles support *device* mode and in that case will run
+in *no* mode instead.
 
 The optimal choice for these keywords depends on the input script and
 the hardware used. The *no* value is useful for verifying that the

doc/src/pair_amoeba.rst

@@ -0,0 +1,253 @@
+.. index:: pair_style amoeba
+.. index:: pair_style hippo
+
+pair_style amoeba command
+=========================
+
+pair_style hippo command
+========================
+Syntax
+""""""
+
+.. code-block:: LAMMPS
+
+   pair_style style
+
+* style = *amoeba* or *hippo*
+
+Examples
+""""""""
+
+.. code-block:: LAMMPS
+
+   pair_style amoeba
+   pair_coeff * * protein.prm.amoeba protein.key.amoeba
+
+.. code-block:: LAMMPS
+
+   pair_style hippo
+   pair_coeff * * water.prm.hippo water.key.hippo
+
+
+Additional info
+"""""""""""""""
+
+* :doc:`Howto amoeba <Howto_amoeba>`
+* examples/amoeba
+* tools/amoeba
+* potentials/\*.amoeba
+* potentials/\*.hippo
+
+Description
+"""""""""""
+
+The *amoeba* style computes the AMOEBA polarizable field formulated
+by Jay Ponder's group at the U Washington at St Louis :ref:`(Ren)
+<amoeba-Ren>`, :ref:`(Shi) <amoeba-Shi>`.  The *hippo* style computes
+the HIPPO polarizable force field, an extension to AMOEBA, formulated
+by Josh Rackers and collaborators in the Ponder group :ref:`(Rackers)
+<amoeba-Rackers>`.
+
+These force fields can be used when polarization effects are desired
+in simulations of water, organic molecules, and biomolecules including
+proteins, provided that parameterizations (Tinker PRM force field
+files) are available for the systems you are interested in.  Files in
+the LAMMPS potentials directory with a "amoeba" or "hippo" suffix can
+be used.  The Tinker distribution and website have additional force
+field files as well.
+
+As discussed on the :doc:`Howto amoeba <Howto_amoeba>` doc page, the
+intermolecular (non-bonded) portion of the AMOEBA force field contains
+these terms:
+
+.. math::
+
+   U_{amoeba} = U_{multipole} + U_{polar} + U_{hal}
+
+while the HIPPO force field contains these terms:
+
+.. math::
+
+   U_{hippo} = U_{multipole} + U_{polar} + U_{qxfer} + U_{repulsion} + U_{dispersion}
+
+Conceptually, these terms compute the following interactions:
+
+* :math:`U_{hal}` = buffered 14-7 van der Waals with offsets applied to hydrogen atoms
+* :math:`U_{repulsion}` = Pauli repulsion due to rearrangement of electron density
+* :math:`U_{dispersion}` = dispersion between correlated, instantaneous induced dipole moments
+* :math:`U_{multipole}` = electrostatics between permanent point charges, dipoles, and quadrupoles
+* :math:`U_{polar}` = electronic polarization between induced point dipoles
+* :math:`U_{qxfer}` = charge transfer effects
+
+Note that the AMOEBA versus HIPPO force fields typically compute the
+same term differently using their own formulas.  The references on
+this doc page give full details for both force fields.
+
+The formulas for the AMOEBA energy terms are:
+
+.. math::
+
+   U_{hal} = \epsilon_{ij} \left( \frac{1.07}{\rho_{ij} + 0.07} \right)^7 \left( \frac{1.12}{\rho_{ij}^7 + 0.12} - 2 \right)
+   U_{multipole} = \vec{M_i}\bold{T_{ij}}\vec{M_j}
+      \vec{M} = \left( q, \vec{\mu_{perm}}, \bold{\Theta} \right)
+   U_{polar} = \frac{1}{2}\vec{\mu_i}^{ind} \vec{E_i}^{perm}
+
+The formulas for the HIPPO energy terms are:
+
+.. math::
+
+   U_{multipole} = Z_i \frac{1}{r_{ij}} Z_j + Z_i T_{ij}^{damp} \vec{M_j} + Z_j T_{ji}^{damp} \vec{M_i} + \vec{M_i} T_{ij}^{damp} \vec{M_j}
+      \vec{M} = \left( Q, \vec{\mu_{perm}}, \bold{\Theta} \right)
+   U_{polar} = \frac{1}{2}\vec{\mu_i}^{ind} \vec{E_i}^{perm}
+   U_{qxfer} = \epsilon_i e^{-\eta_j r_{ij}} + \epsilon_j e^{-\eta_i r_{ij}}
+   U_{repulsion} = \frac{K_i K_j}{r_{ij}} S^2
+      S^2 = \left( \int{\phi_i \phi_j} dv \right)^2 = \vec{M_i}\bold{T_{ij}^{repulsion}}\vec{M_j}
+   U_{dispersion} = -\frac{C_6^iC_6^j}{r_{ij}^6} \left( f_{damp}^{dispersion} \right)_{ij}^2
+
+.. note::
+
+  The AMOEBA and HIPPO force fields compute long-range charge, dipole,
+  and quadrupole interactions as well as long-range dispersion
+  effects.  However, unlike other models with long-range interactions
+  in LAMMPS, this does not require use of a KSpace style via the
+  :doc:`kspace_style <kspace_style>` command.  That is because for
+  AMOEBA and HIPPO the long-range computations are intertwined with
+  the pairwise computations.  So these pair style include both short-
+  and long-range computations.  This means the energy and virial
+  computed by the pair style as well as the "Pair" timing reported by
+  LAMMPS will include the long-range calculations.
+
+The implementation of the AMOEBA and HIPPO force fields in LAMMPS was
+done using F90 code provided by the Ponder group from their `Tinker MD
+code <https://dasher.wustl.edu/tinker/>`_.
+
+The current implementation (July 2022) of AMOEBA in LAMMPS matches the
+version discussed in :ref:`(Ponder) <amoeba-Ponder>`, :ref:`(Ren)
+<amoeba-Ren>`, and :ref:`(Shi) <amoeba-Shi>`.  Likewise the current
+implementation of HIPPO in LAMMPS matches the version discussed in
+:ref:`(Rackers) <amoeba-Rackers>`.
+
+----------
+
+Only a single pair_coeff command is used with either the *amoeba* and
+*hippo* styles which specifies two Tinker files, a PRM and KEY file.
+
+.. code-block:: LAMMPS
+
+   pair_coeff * * ../potentials/protein.prm.amoeba ../potentials/protein.key.amoeba
+   pair_coeff * * ../potentials/water.prm.hippo ../potentials/water.key.hippo
+
+Examples of the PRM files are in the potentials directory with an
+\*.amoeba or \*.hippo suffix.  The examples/amoeba directory has
+examples of both PRM and KEY files.
+
+A Tinker PRM file is composed of sections, each of which has multiple
+lines.  A Tinker KEY file is composed of lines, each of which has a
+keyword followed by zero or more parameters.
+
+The list of PRM sections and KEY keywords which LAMMPS recognizes are
+listed on the :doc:`Howto amoeba <Howto_amoeba>` doc page.  If not
+recognized, the section or keyword is skipped.
+
+Note that if the KEY file is specified as NULL, then no file is
+required; default values for various AMOEBA/HIPPO settings are used.
+The :doc:`Howto amoeba <Howto_amoeba>` doc page also gives the default
+settings.
+
+----------
+
+Mixing, shift, table, tail correction, restart, rRESPA info
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+These pair styles do not support the :doc:`pair_modify <pair_modify>`
+mix, shift, table, and tail options.
+
+These pair styles do not write their information to :doc:`binary
+restart files <restart>`, since it is stored in potential files.
+Thus, you need to re-specify the pair_style and pair_coeff commands in
+an input script that reads a restart file.
+
+These pair styles 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
+""""""""""""
+
+These pair styles are part of the AMOEBA package.  They are only
+enabled if LAMMPS was built with that package.  See the :doc:`Build
+package <Build_package>` doc page for more info.
+
+The AMOEBA and HIPPO potential (PRM) and KEY files provided with
+LAMMPS in the potentials and examples/amoeba directories are Tinker
+files parameterized for Tinker units.  Their numeric parameters are
+converted by LAMMPS to its real units :doc:`units <units>`.  Thus you
+can only use these pair styles with real units.
+
+These potentials do not yet calculate per-atom energy or virial
+contributions.
+
+As explained on the :doc:`AMOEBA and HIPPO howto <Howto_amoeba>` page,
+use of these pair styles to run a simulation with the AMOEBA or HIPPO
+force fields requires several things.
+
+The first is a data file generated by the tools/tinker/tinker2lmp.py
+conversion script which uses Tinker file force field file input to
+create a data file compatible with LAMMPS.
+
+The second is use of these commands:
+
+* :doc:`atom_style amoeba <atom_style>`
+* :doc:`fix property/atom <fix_property_atom>`
+* :doc:`special_bonds one/five <special_bonds>`
+
+And third, depending on the model being simulated, these
+commands for intramolecular interactions may also be required:
+
+* :doc:`bond_style class2 <bond_class2>`
+* :doc:`angle_style amoeba <angle_amoeba>`
+* :doc:`dihedral_style fourier <dihedral_fourier>`
+* :doc:`improper_style amoeba <improper_amoeba>`
+* :doc:`fix amoeba/pitorsion <fix_amoeba_pitorsion>`
+* :doc:`fix amoeba/bitorsion <fix_amoeba_bitorsion>`
+
+----------
+
+Related commands
+""""""""""""""""
+
+:doc:`atom_style amoeba <atom_style>`,
+:doc:`bond_style class2 <bond_class2>`,
+:doc:`angle_style amoeba <angle_amoeba>`,
+:doc:`dihedral_style fourier <dihedral_fourier>`,
+:doc:`improper_style amoeba <improper_amoeba>`,
+:doc:`fix amoeba/pitorsion <fix_amoeba_pitorsion>`,
+:doc:`fix amoeba/bitorsion <fix_amoeba_bitorsion>`,
+:doc:`special_bonds one/five <special_bonds>`,
+:doc:`fix property/atom <fix_property_atom>`
+
+Default
+"""""""
+
+none
+
+----------
+
+.. _amoeba-Ponder:
+
+**(Ponder)** Ponder, Wu, Ren, Pande, Chodera, Schnieders, Haque, Mobley, Lambrecht, DiStasio Jr, M. Head-Gordon, Clark, Johnson, T. Head-Gordon, J Phys Chem B, 114, 2549-2564 (2010).
+
+.. _amoeba-Rackers:
+
+**(Rackers)** Rackers, Silva, Wang, Ponder, J Chem Theory Comput, 17, 7056-7084 (2021).
+
+.. _amoeba-Ren:
+
+**(Ren)** Ren and Ponder, J Phys Chem B, 107, 5933 (2003).
+
+.. _amoeba-Shi:
+
+**(Shi)** Shi, Xia, Zhang, Best, Wu, Ponder, Ren, J Chem Theory Comp, 9, 4046, 2013.
+

doc/src/pair_dipole.rst

@@ -78,12 +78,12 @@ Examples
 Description
 """""""""""
 
-Style *lj/cut/dipole/cut* computes interactions between pairs of particles
-that each have a charge and/or a point dipole moment.  In addition to
-the usual Lennard-Jones interaction between the particles (Elj) the
-charge-charge (Eqq), charge-dipole (Eqp), and dipole-dipole (Epp)
-interactions are computed by these formulas for the energy (E), force
-(F), and torque (T) between particles I and J.
+Style *lj/cut/dipole/cut* computes interactions between pairs of
+particles that each have a charge and/or a point dipole moment.  In
+addition to the usual Lennard-Jones interaction between the particles
+(Elj) the charge-charge (Eqq), charge-dipole (Eqp), and dipole-dipole
+(Epp) interactions are computed by these formulas for the energy (E),
+force (F), and torque (T) between particles I and J.
 
 .. math::
 
@@ -112,18 +112,18 @@ interactions are computed by these formulas for the energy (E), force
                       \frac{3}{r^5} (\vec{p_i} \bullet \vec{r})
                       (\vec{p_j} \times \vec{r})
 
-where :math:`q_i` and :math:`q_j` are the charges on the two particles,
-:math:`\vec{p_i}` and :math:`\vec{p_j}` are the dipole moment vectors of
-the two particles, r is their separation distance, and the vector r =
-Ri - Rj is the separation vector between the two particles.  Note that
-Eqq and Fqq are simply Coulombic energy and force, Fij = -Fji as
-symmetric forces, and Tij != -Tji since the torques do not act
-symmetrically.  These formulas are discussed in :ref:`(Allen) <Allen2>`
-and in :ref:`(Toukmaji) <Toukmaji2>`.
+where :math:`q_i` and :math:`q_j` are the charges on the two
+particles, :math:`\vec{p_i}` and :math:`\vec{p_j}` are the dipole
+moment vectors of the two particles, r is their separation distance,
+and the vector r = Ri - Rj is the separation vector between the two
+particles.  Note that Eqq and Fqq are simply Coulombic energy and
+force, Fij = -Fji as symmetric forces, and Tij != -Tji since the
+torques do not act symmetrically.  These formulas are discussed in
+:ref:`(Allen) <Allen2>` and in :ref:`(Toukmaji) <Toukmaji2>`.
 
 Also note, that in the code, all of these terms (except Elj) have a
-:math:`C/\epsilon` prefactor, the same as the Coulombic term in the LJ +
-Coulombic pair styles discussed :doc:`here <pair_lj>`.  C is an
+:math:`C/\epsilon` prefactor, the same as the Coulombic term in the
+LJ + Coulombic pair styles discussed :doc:`here <pair_lj>`.  C is an
 energy-conversion constant and epsilon is the dielectric constant
 which can be set by the :doc:`dielectric <dielectric>` command.  The
 same is true of the equations that follow for other dipole pair
@@ -135,11 +135,11 @@ moment. In general, a shifted-force potential is a (slightly) modified
 potential containing extra terms that make both the energy and its
 derivative go to zero at the cutoff distance; this removes
 (cutoff-related) problems in energy conservation and any numerical
-instability in the equations of motion :ref:`(Allen) <Allen2>`. Shifted-force
-interactions for the Lennard-Jones (E_LJ), charge-charge (Eqq),
-charge-dipole (Eqp), dipole-charge (Epq) and dipole-dipole (Epp)
-potentials are computed by these formulas for the energy (E), force
-(F), and torque (T) between particles I and J:
+instability in the equations of motion :ref:`(Allen)
+<Allen2>`. Shifted-force interactions for the Lennard-Jones (E_LJ),
+charge-charge (Eqq), charge-dipole (Eqp), dipole-charge (Epq) and
+dipole-dipole (Epp) potentials are computed by these formulas for the
+energy (E), force (F), and torque (T) between particles I and J:
 
 .. math::
 
@@ -207,65 +207,52 @@ potentials are computed by these formulas for the energy (E), force
 where :math:`\epsilon` and :math:`\sigma` are the standard LJ
 parameters, :math:`r_c` is the cutoff, :math:`q_i` and :math:`q_j` are
 the charges on the two particles, :math:`\vec{p_i}` and
-:math:`\vec{p_j}` are the dipole moment vectors of the two particles, r
-is their separation distance, and the vector r = Ri - Rj is the
-separation vector between the two particles.  Note that Eqq and Fqq are
-simply Coulombic energy and force, Fij = -Fji as symmetric forces, and
-Tij != -Tji since the torques do not act symmetrically.  The
+:math:`\vec{p_j}` are the dipole moment vectors of the two particles,
+r is their separation distance, and the vector r = Ri - Rj is the
+separation vector between the two particles.  Note that Eqq and Fqq
+are simply Coulombic energy and force, Fij = -Fji as symmetric forces,
+and Tij != -Tji since the torques do not act symmetrically.  The
 shifted-force formula for the Lennard-Jones potential is reported in
 :ref:`(Stoddard) <Stoddard>`.  The original (non-shifted) formulas for
 the electrostatic potentials, forces and torques can be found in
-:ref:`(Price) <Price2>`. The shifted-force electrostatic potentials have
-been obtained by applying equation 5.13 of :ref:`(Allen) <Allen2>`. The
-formulas for the corresponding forces and torques have been obtained by
-applying the 'chain rule' as in appendix C.3 of :ref:`(Allen) <Allen2>`.
+:ref:`(Price) <Price2>`. The shifted-force electrostatic potentials
+have been obtained by applying equation 5.13 of :ref:`(Allen)
+<Allen2>`. The formulas for the corresponding forces and torques have
+been obtained by applying the 'chain rule' as in appendix C.3 of
+:ref:`(Allen) <Allen2>`.
 
 If one cutoff is specified in the pair_style command, it is used for
 both the LJ and Coulombic (q,p) terms.  If two cutoffs are specified,
 they are used as cutoffs for the LJ and Coulombic (q,p) terms
-respectively. This pair style also supports an optional *scale* keyword
-as part of a pair_coeff statement, where the interactions can be
-scaled according to this factor. This scale factor is also made available
-for use with fix adapt.
+respectively. This pair style also supports an optional *scale*
+keyword as part of a pair_coeff statement, where the interactions can
+be scaled according to this factor. This scale factor is also made
+available for use with fix adapt.
 
-Style *lj/cut/dipole/long* computes long-range point-dipole
-interactions as discussed in :ref:`(Toukmaji) <Toukmaji2>`. Dipole-dipole,
-dipole-charge, and charge-charge interactions are all supported, along
-with the standard 12/6 Lennard-Jones interactions, which are computed
-with a cutoff.  A :doc:`kspace_style <kspace_style>` must be defined to
-use this pair style.  Currently, only :doc:`kspace_style ewald/disp <kspace_style>` support long-range point-dipole
-interactions.
+Style *lj/cut/dipole/long* computes the short-range portion of
+point-dipole interactions as discussed in :ref:`(Toukmaji)
+<Toukmaji2>`. Dipole-dipole, dipole-charge, and charge-charge
+interactions are all supported, along with the standard 12/6
+Lennard-Jones interactions, which are computed with a cutoff.  A
+:doc:`kspace_style <kspace_style>` must be defined to use this pair
+style.  If only dipoles (not point charges) are included in the model,
+the kspace style can be one of these 3 options, all of which compute
+the long-range portion of dipole-dipole interactions.  If the model
+includes point charges (in addition to dipoles), then only the first
+of these kspace styles can be used:
 
-Style *lj/long/dipole/long* also computes point-dipole interactions as
-discussed in :ref:`(Toukmaji) <Toukmaji2>`. Long-range dipole-dipole,
-dipole-charge, and charge-charge interactions are all supported, along
-with the standard 12/6 Lennard-Jones interactions.  LJ interactions
-can be cutoff or long-ranged.
+* :doc:`kspace_style ewald/disp <kspace_style>`
+* :doc:`kspace_style ewald/dipole <kspace_style>`
+* :doc:`kspace_style pppm/dipole <kspace_style>`
 
-For style *lj/long/dipole/long*, if *flag_lj* is set to *long*, no
-cutoff is used on the LJ 1/r\^6 dispersion term.  The long-range
-portion is calculated by using the :doc:`kspace_style ewald_disp <kspace_style>` command.  The specified LJ cutoff then
-determines which portion of the LJ interactions are computed directly
-by the pair potential versus which part is computed in reciprocal
-space via the Kspace style.  If *flag_lj* is set to *cut*, the LJ
-interactions are simply cutoff, as with :doc:`pair_style lj/cut <pair_lj>`.  If *flag_lj* is set to *off*, LJ interactions
-are not computed at all.
+Style *lj/long/dipole/long* has the same functionality as style
+*lj/cut/dipole/long*, except it also has an option to compute 12/6
+Lennard-Jones interactions for use with a long-range dispersion kspace
+style.  This is done by setting its *flag_lj* argument to *long*.  For
+long-range LJ interactions, the doc:`kspace_style ewald/disp
+<kspace_style>` command must be used.
 
-If *flag_coul* is set to *long*, no cutoff is used on the Coulombic or
-dipole interactions.  The long-range portion is calculated by using
-*ewald_disp* of the :doc:`kspace_style <kspace_style>` command. If
-*flag_coul* is set to *off*, Coulombic and dipole interactions are not
-computed at all.
-
-Atoms with dipole moments should be integrated using the :doc:`fix nve/sphere update dipole <fix_nve_sphere>` or the :doc:`fix nvt/sphere update dipole <fix_nvt_sphere>` command to rotate the
-dipole moments.  The *omega* option on the :doc:`fix langevin <fix_langevin>` command can be used to thermostat the
-rotational motion.  The :doc:`compute temp/sphere <compute_temp_sphere>`
-command can be used to monitor the temperature, since it includes
-rotational degrees of freedom.  The :doc:`atom_style hybrid dipole sphere <atom_style>` command should be used since
-it defines the point dipoles and their rotational state.
-The magnitude and orientation of the dipole moment for each particle
-can be defined by the :doc:`set <set>` command or in the "Atoms" section
-of the data file read in by the :doc:`read_data <read_data>` command.
+----------
 
 The following coefficients must be defined for each pair of atoms
 types via the :doc:`pair_coeff <pair_coeff>` command as in the examples
@@ -287,6 +274,40 @@ type pair.
 
 ----------
 
+Note that for systems using these pair styles, typically particles
+should be able to exert torque on each other via their dipole moments
+so that the particle and its dipole moment can rotate.  This requires
+they not be point particles, but finite-size spheres.  Thus you should
+use a command like :doc:`atom_style hybrid sphere dipole <atom_style>`
+to use particles with both attributes.
+
+The magnitude and orientation of the dipole moment for each particle
+can be defined by the :doc:`set <set>` command or in the "Atoms"
+section of the data file read in by the :doc:`read_data <read_data>`
+command.
+
+Rotating finite-size particles have 6 degrees of freedom (DOFs),
+translation and rotational.  You can use the :doc:`compute temp/sphere
+<compute_temp_sphere>` command to monitor a temperature which includes
+all these DOFs.
+
+Finite-size particles with dipole moments should be integrated using
+one of these options:
+
+* :doc:`fix nve/sphere update dipole <fix_nve_sphere>`
+* :doc:`fix nve/sphere update dipole <fix_nve_sphere>` plus :doc:`fix langevin omega yes <fix_langevin>`
+* :doc:`fix nvt/sphere update dipole <fix_nvt_sphere>`
+* :doc:`fix npt/sphere update dipole <fix_npt_sphere>`
+
+In all cases the "update dipole" setting insures the dipole moments
+are also rotated when the finite-size spheres rotate.  The 2nd and 3rd
+bullets perform thermostatting; in the case of a Langevin thermostat
+the "omega yes" option also thermostats the rotational degrees of
+freedom (if desired).  The 4th bullet performs thermostatting and
+barostatting.
+
+----------
+
 .. include:: accel_styles.rst
 
 ----------

doc/src/pair_meam.rst

@@ -1,8 +1,11 @@
 .. index:: pair_style meam
+.. index:: pair_style meam/kk
 
 pair_style meam command
 =========================
 
+Accelerator Variants: *meam/kk*
+
 Syntax
 """"""
 
@@ -347,6 +350,12 @@ Most published MEAM parameter sets use the default values *attrac* = *repulse* =
 Setting *repuls* = *attrac* = *delta* corresponds to the form used in several
 recent published MEAM parameter sets, such as :ref:`(Valone) <Valone>`
 
+----------
+
+.. include:: accel_styles.rst
+
+----------
+
 .. note::
 
    The default form of the *erose* expression in LAMMPS was corrected

doc/src/pair_modify.rst

@@ -13,7 +13,7 @@ Syntax
 * one or more keyword/value pairs may be listed
 * keyword = *pair* or *shift* or *mix* or *table* or *table/disp* or *tabinner*
   or *tabinner/disp* or *tail* or *compute* or *nofdotr* or *special* or
-  *compute/tally*
+  *compute/tally* or *neigh/trim*
 
   .. parsed-literal::
 
@@ -37,6 +37,7 @@ Syntax
           which = *lj/coul* or *lj* or *coul*
           w1,w2,w3 = 1-2, 1-3, 1-4 weights from 0.0 to 1.0 inclusive
        *compute/tally* value = *yes* or *no*
+       *neigh/trim* value = *yes* or *no*
 
 Examples
 """"""""
@@ -283,6 +284,31 @@ the *pair* keyword.  Use *no* to disable, or *yes* to enable.
    The "pair_modify pair compute/tally" command must be issued
    **before** the corresponding compute style is defined.
 
+The *neigh/trim* keyword controls whether an explicit cutoff is set for
+each neighbor list request issued by individual pair sub-styles when
+using :doc:`pair hybrid/overlay <pair_hybrid>`.  When this keyword is
+set to *no*, then the cutoff of each pair sub-style neighbor list will
+be set equal to the largest cutoff, even if a shorter cutoff is
+specified for a particular sub-style.  If possible the neighbor list
+will be copied directly from another list.  When this keyword is set to
+*yes* then the cutoff of the neighbor list will be explicitly set to the
+value requested by the pair sub-style, and if possible the list will be
+created by trimming neighbors from another list with a longer cutoff,
+otherwise a new neighbor list will be created with the specified cutoff.
+The *yes* option can be faster when there are multiple pair styles with
+different cutoffs since the number of pair-wise distance checks between
+neighbors is reduced (but the time required to build the neighbor lists
+is increased). The *no* option could be faster when two or more neighbor
+lists have similar (but not exactly the same) cutoffs.
+
+.. note::
+
+   The "pair_modify neigh/trim" command *only* applies when there are
+   multiple pair sub-styles for the same atoms with different cutoffs,
+   i.e. when using pair style hybrid/overlay.  If you have different
+   cutoffs for different pairs for atoms type, the :doc:`neighbor style
+   multi <neighbor>` should be used to create optimized neighbor lists.
+
 ----------
 
 Restrictions
@@ -298,13 +324,13 @@ Related commands
 
 :doc:`pair_style <pair_style>`, :doc:`pair_style hybrid <pair_hybrid>`,
 :doc:`pair_coeff <pair_coeff>`, :doc:`thermo_style <thermo_style>`,
-:doc:`compute \*/tally <compute_tally>`
+:doc:`compute \*/tally <compute_tally>`, :doc:`neighbor multi <neighbor>`
 
 Default
 """""""
 
 The option defaults are mix = geometric, shift = no, table = 12,
-tabinner = sqrt(2.0), tail = no, and compute = yes.
+tabinner = sqrt(2.0), tail = no, compute = yes, and neigh/trim yes.
 
 Note that some pair styles perform mixing, but only a certain style of
 mixing.  See the doc pages for individual pair styles for details.

doc/src/pair_smatb.rst

@@ -14,13 +14,9 @@ Syntax
 
    pair_style style args
 
-* style = *smatb*
+* style = *smatb* or *smatb/single*
 * args = none
 
-.. parsed-literal::
-
-     *smatb*
-
 Examples
 """"""""
 
@@ -29,13 +25,18 @@ Examples
    pair_style smatb
    pair_coeff 1 1 2.88 10.35 4.178 0.210 1.818 4.07293506 4.9883063257983666
 
+   pair_style smatb/single
+   pair_coeff 1 1 2.88 10.35 4.178 0.210 1.818 4.07293506 4.9883063257983666
+
 
 Description
 """""""""""
 
-The *smatb* styles compute the Second Moment Approximation to the Tight Binding
-:ref:`(Cyrot) <Cyrot>`, :ref:`(Gupta) <Gupta>`, :ref:`(Rosato) <Rosato>`,
-given by
+.. versionadded:: 4May2022
+
+The *smatb* and *smatb/single* styles compute the Second Moment
+Approximation to the Tight Binding :ref:`(Cyrot) <Cyrot>`,
+:ref:`(Gupta) <Gupta>`, :ref:`(Rosato) <Rosato>`, given by
 
 .. math::
       E_{i}  = \sum_{j,R_{ij}\leq R_{c}} \alpha(R_{ij}) - \sqrt{\sum_{j,R_{ij}\leq R_{c}}\Xi^2(R_{ij})}
@@ -66,6 +67,8 @@ exponential terms and their first and second derivatives are smoothly
 reduced to zero, from the inner cutoff :math:`R_{sc}` to the outer
 cutoff :math:`R_{c}`.
 
+The *smatb/single* style is an optimization when using only a single atom type.
+
 Coefficients
 """"""""""""
 
@@ -100,10 +103,10 @@ For atom type pairs I,J and I != J the coefficients are not automatically mixed.
 Restrictions
 """"""""""""
 
-This pair style is part of the SMTBQ package and is only enabled
+These pair styles are part of the SMTBQ package and are only enabled
 if LAMMPS is built with that package.  See the :doc:`Build package <Build_package>` page for more info.
 
-These pair potentials require the :doc:`newton <newton>` setting to be "on" for pair interactions.
+These pair styles require the :doc:`newton <newton>` setting to be "on" for pair interactions.
 
 Related commands
 """"""""""""""""

doc/src/pair_spica.rst

@@ -1,27 +1,27 @@
-.. index:: pair_style lj/sdk
-.. index:: pair_style lj/sdk/gpu
-.. index:: pair_style lj/sdk/kk
-.. index:: pair_style lj/sdk/omp
-.. index:: pair_style lj/sdk/coul/long
-.. index:: pair_style lj/sdk/coul/long/gpu
-.. index:: pair_style lj/sdk/coul/long/omp
-.. index:: pair_style lj/sdk/coul/msm
-.. index:: pair_style lj/sdk/coul/msm/omp
+.. index:: pair_style lj/spica
+.. index:: pair_style lj/spica/gpu
+.. index:: pair_style lj/spica/kk
+.. index:: pair_style lj/spica/omp
+.. index:: pair_style lj/spica/coul/long
+.. index:: pair_style lj/spica/coul/long/gpu
+.. index:: pair_style lj/spica/coul/long/omp
+.. index:: pair_style lj/spica/coul/msm
+.. index:: pair_style lj/spica/coul/msm/omp
 
-pair_style lj/sdk command
-=========================
+pair_style lj/spica command
+===========================
 
-Accelerator Variants: *lj/sdk/gpu*, *lj/sdk/kk*, *lj/sdk/omp*
+Accelerator Variants: *lj/spica/gpu*, *lj/spica/kk*, *lj/spica/omp*
 
-pair_style lj/sdk/coul/long command
-===================================
+pair_style lj/spica/coul/long command
+=====================================
 
-Accelerator Variants: *lj/sdk/coul/long/gpu*, *lj/sdk/coul/long/omp*
+Accelerator Variants: *lj/spica/coul/long/gpu*, *lj/spica/coul/long/omp*
 
-pair_style lj/sdk/coul/msm command
-==================================
+pair_style lj/spica/coul/msm command
+====================================
 
-Accelerator Variants: *lj/sdk/coul/msm/omp*
+Accelerator Variants: *lj/spica/coul/msm/omp*
 
 Syntax
 """"""
@@ -30,14 +30,14 @@ Syntax
 
    pair_style style args
 
-* style = *lj/sdk* or *lj/sdk/coul/long*
+* style = *lj/spica* or *lj/spica/coul/long*
 * args = list of arguments for a particular style
 
 .. parsed-literal::
 
-     *lj/sdk* args = cutoff
+     *lj/spica* args = cutoff
        cutoff = global cutoff for Lennard Jones interactions (distance units)
-     *lj/sdk/coul/long* args = cutoff (cutoff2)
+     *lj/spica/coul/long* args = cutoff (cutoff2)
        cutoff = global cutoff for LJ (and Coulombic if only 1 arg) (distance units)
        cutoff2 = global cutoff for Coulombic (optional) (distance units)
 
@@ -46,21 +46,21 @@ Examples
 
 .. code-block:: LAMMPS
 
-   pair_style lj/sdk 2.5
+   pair_style lj/spica 2.5
    pair_coeff 1 1 lj12_6 1 1.1 2.8
 
-   pair_style lj/sdk/coul/long 10.0
-   pair_style lj/sdk/coul/long 10.0 12.0
+   pair_style lj/spica/coul/long 10.0
+   pair_style lj/spica/coul/long 10.0 12.0
    pair_coeff 1 1 lj9_6 100.0 3.5 12.0
 
-   pair_style lj/sdk/coul/msm 10.0
-   pair_style lj/sdk/coul/msm 10.0 12.0
+   pair_style lj/spica/coul/msm 10.0
+   pair_style lj/spica/coul/msm 10.0 12.0
    pair_coeff 1 1 lj9_6 100.0 3.5 12.0
 
 Description
 """""""""""
 
-The *lj/sdk* styles compute a 9/6, 12/4, or 12/6 Lennard-Jones potential,
+The *lj/spica* styles compute a 9/6, 12/4, 12/5, or 12/6 Lennard-Jones potential,
 given by
 
 .. math::
@@ -71,14 +71,20 @@ given by
    E = & \frac{3\sqrt{3}}{2} \epsilon \left[ \left(\frac{\sigma}{r}\right)^{12} -
                          \left(\frac{\sigma}{r}\right)^4 \right]
                          \qquad r < r_c \\
+   E = & \frac{12}{7}\left(\frac{12}{5}\right)^{\left(\frac{5}{7}\right)} \epsilon
+                         \left[ \left(\frac{\sigma}{r}\right)^{12} -
+                         \left(\frac{\sigma}{r}\right)^5 \right]
+                         \qquad r < r_c \\
    E = &  4 \epsilon  \left[ \left(\frac{\sigma}{r}\right)^{12} -
                          \left(\frac{\sigma}{r}\right)^6 \right]
                          \qquad r < r_c
 
-as required for the SDK Coarse-grained MD parameterization discussed in
-:ref:`(Shinoda) <Shinoda3>` and :ref:`(DeVane) <DeVane>`.  Rc is the cutoff.
+as required for the SPICA (formerly called SDK) and the pSPICA Coarse-grained MD parameterization discussed in
+:ref:`(Shinoda) <Shinoda3>`, :ref:`(DeVane) <DeVane>`, :ref:`(Seo) <Seo>`, and :ref:`(Miyazaki) <Miyazaki>`.
+Rc is the cutoff.
+Summary information on these force fields can be found at https://www.spica-ff.org
 
-Style *lj/sdk/coul/long* computes the adds Coulombic interactions
+Style *lj/spica/coul/long* computes the adds Coulombic interactions
 with an additional damping factor applied so it can be used in
 conjunction with the :doc:`kspace_style <kspace_style>` command and
 its *ewald* or *pppm* or *pppm/cg* option.  The Coulombic cutoff
@@ -92,7 +98,7 @@ above, or in the data file or restart files read by the
 :doc:`read_data <read_data>` or :doc:`read_restart <read_restart>`
 commands, or by mixing as described below:
 
-* cg_type (lj9_6, lj12_4, or lj12_6)
+* cg_type (lj9_6, lj12_4, lj12_5, or lj12_6)
 * epsilon (energy units)
 * sigma (distance units)
 * cutoff1 (distance units)
@@ -108,11 +114,15 @@ and Coulombic interactions for this type pair.  If both coefficients
 are specified, they are used as the LJ and Coulombic cutoffs for this
 type pair.
 
-For *lj/sdk/coul/long* and *lj/sdk/coul/msm* only the LJ cutoff can be
+For *lj/spica/coul/long* and *lj/spica/coul/msm* only the LJ cutoff can be
 specified since a Coulombic cutoff cannot be specified for an
 individual I,J type pair.  All type pairs use the same global
 Coulombic cutoff specified in the pair_style command.
 
+The original implementation of the above styles are
+style *lj/sdk*, *lj/sdk/coul/long*, and *lj/sdk/coul/msm*,
+and available for backward compatibility.
+
 ----------
 
 .. include:: accel_styles.rst
@@ -123,24 +133,24 @@ Mixing, shift, table, tail correction, restart, rRESPA info
 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
 
 For atom type pairs I,J and I != J, the epsilon and sigma coefficients
-and cutoff distance for all of the lj/sdk pair styles *cannot* be mixed,
+and cutoff distance for all of the lj/spica pair styles *cannot* be mixed,
 since different pairs may have different exponents. So all parameters
 for all pairs have to be specified explicitly through the "pair_coeff"
 command. Defining then in a data file is also not supported, due to
 limitations of that file format.
 
-All of the lj/sdk pair styles support the
+All of the lj/spica pair styles support the
 :doc:`pair_modify <pair_modify>` shift option for the energy of the
 Lennard-Jones portion of the pair interaction.
 
-The *lj/sdk/coul/long* pair styles support the
+The *lj/spica/coul/long* pair styles support the
 :doc:`pair_modify <pair_modify>` table option since they can tabulate
 the short-range portion of the long-range Coulombic interaction.
 
-All of the lj/sdk pair styles write their information to :doc:`binary restart files <restart>`, so pair_style and pair_coeff commands do
+All of the lj/spica pair styles write their 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.
 
-The lj/sdk and lj/cut/coul/long pair styles do not support
+The lj/spica and lj/cut/coul/long pair styles do not support
 the use of the *inner*, *middle*, and *outer* keywords of the :doc:`run_style respa <run_style>` command.
 
 ----------
@@ -148,8 +158,8 @@ the use of the *inner*, *middle*, and *outer* keywords of the :doc:`run_style re
 Restrictions
 """"""""""""
 
-All of the lj/sdk pair styles are part of the CG-SDK package.  The
-*lj/sdk/coul/long* style also requires the KSPACE package to be built
+All of the lj/spica pair styles are part of the CG-SPICA package.  The
+*lj/spica/coul/long* style also requires the KSPACE package to be built
 (which is enabled by default).  They are only enabled if LAMMPS was
 built with that package.  See the :doc:`Build package <Build_package>`
 doc page for more info.
@@ -157,7 +167,7 @@ doc page for more info.
 Related commands
 """"""""""""""""
 
-:doc:`pair_coeff <pair_coeff>`, :doc:`angle_style sdk <angle_sdk>`
+:doc:`pair_coeff <pair_coeff>`, :doc:`angle_style spica <angle_spica>`
 
 Default
 """""""
@@ -168,8 +178,16 @@ none
 
 .. _Shinoda3:
 
-**(Shinoda)** Shinoda, DeVane, Klein, Mol Sim, 33, 27 (2007).
+**(Shinoda)** Shinoda, DeVane, Klein, Mol Sim, 33, 27-36 (2007).
 
 .. _DeVane:
 
 **(DeVane)**  Shinoda, DeVane, Klein, Soft Matter, 4, 2453-2462 (2008).
+
+.. _Seo:
+
+**(Seo)**  Seo, Shinoda, J Chem Theory Comput, 15, 762-774 (2019).
+
+.. _Miyazaki:
+
+**(Miyazaki)**  Miyazaki, Okazaki, Shinoda, J Chem Theory Comput, 16, 782-793 (2020).

doc/src/pair_srp.rst

@@ -1,18 +1,23 @@
 .. index:: pair_style srp
+.. index:: pair_style srp/react
 
 pair_style srp command
 ======================
 
+pair_style srp/react command
+============================
 Syntax
 """"""
 
 .. code-block:: LAMMPS
 
    pair_style srp cutoff btype dist keyword value ...
+   pair_style srp/react cutoff btype dist react-id keyword value ...
 
 * cutoff = global cutoff for SRP interactions (distance units)
 * btype = bond type to apply SRP interactions to (can be wildcard, see below)
 * distance = *min* or *mid*
+* react-id = id of either fix bond/break or fix bond/create
 * zero or more keyword/value pairs may be appended
 * keyword = *exclude*
 
@@ -36,13 +41,19 @@ Examples
    pair_coeff 1 2 none
    pair_coeff 2 2 srp 40.0
 
+   fix        create all bond/create   100  1  2  1.0 1 prob  0.2  19852
+   pair_style hybrid dpd 1.0 1.0 12345 srp/react 0.8 * min create exclude yes
+   pair_coeff 1 1 dpd 60.0 50 1.0
+   pair_coeff 1 2 none
+   pair_coeff 2 2 srp/react 40.0
+
    pair_style hybrid srp 0.8 2 mid
    pair_coeff 1 1 none
    pair_coeff 1 2 none
    pair_coeff 2 2 srp 100.0 0.8
 
 Description
-"""""""""""
+
 
 Style *srp* computes a soft segmental repulsive potential (SRP) that
 acts between pairs of bonds. This potential is useful for preventing
@@ -121,6 +132,18 @@ at the cutoff distance :math:`r_c`.
 
 ----------
 
+Pair style *srp/react* interfaces the pair style *srp* with the
+bond breaking and formation mechanisms provided by fix *bond/break*
+and fix *bond/create*, respectively. When using this pair style, whenever a
+bond breaking (or formation) reaction occurs, the corresponding fictitious
+particle is deleted (or inserted) during the same simulation time step as
+the reaction. This is useful in the simulation of reactive systems involving
+large polymeric molecules :ref:`(Palkar) <Palkar>`  where the segmental repulsive
+potential is necessary to minimize topological violations, and also needs to be
+turned on and off according to the progress of the reaction.
+
+----------
+
 Mixing, shift, table, tail correction, restart, rRESPA info
 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
 
@@ -178,3 +201,8 @@ The default keyword value is exclude = yes.
 
 **(Sirk)** Sirk TW, Sliozberg YR, Brennan JK, Lisal M, Andzelm JW, J
 Chem Phys, 136 (13) 134903, 2012.
+
+.. _Palkar:
+
+**(Palkar)** Palkar V, Kuksenok O, J. Phys. Chem. B, 126 (1), 336-346, 2022
+

doc/src/pair_style.rst

@@ -116,6 +116,7 @@ accelerated styles exist.
 * :doc:`agni <pair_agni>` - AGNI machine-learning potential
 * :doc:`airebo <pair_airebo>` - AIREBO potential of Stuart
 * :doc:`airebo/morse <pair_airebo>` - AIREBO with Morse instead of LJ
+* :doc:`amoeba <pair_amoeba>` -
 * :doc:`atm <pair_atm>` - Axilrod-Teller-Muto potential
 * :doc:`awpmd/cut <pair_awpmd>` - Antisymmetrized Wave Packet MD potential for atoms and electrons
 * :doc:`beck <pair_beck>` - Beck potential
@@ -202,6 +203,7 @@ accelerated styles exist.
 * :doc:`hbond/dreiding/lj <pair_hbond_dreiding>` - DREIDING hydrogen bonding LJ potential
 * :doc:`hbond/dreiding/morse <pair_hbond_dreiding>` - DREIDING hydrogen bonding Morse potential
 * :doc:`hdnnp <pair_hdnnp>` - High-dimensional neural network potential
+* :doc:`hippo <pair_amoeba>` -
 * :doc:`ilp/graphene/hbn <pair_ilp_graphene_hbn>` - registry-dependent interlayer potential (ILP)
 * :doc:`ilp/tmd <pair_ilp_tmd>` - interlayer potential (ILP) potential for transition metal dichalcogenides (TMD)
 * :doc:`kim <pair_kim>` - interface to potentials provided by KIM project
@@ -258,9 +260,9 @@ accelerated styles exist.
 * :doc:`lj/long/tip4p/long <pair_lj_long>` - long-range LJ and long-range Coulomb for TIP4P water
 * :doc:`lj/mdf <pair_mdf>` - LJ potential with a taper function
 * :doc:`lj/relres <pair_lj_relres>` - LJ using multiscale Relative Resolution (RelRes) methodology :ref:`(Chaimovich) <Chaimovich2>`.
-* :doc:`lj/sdk <pair_sdk>` - LJ for SDK coarse-graining
-* :doc:`lj/sdk/coul/long <pair_sdk>` - LJ for SDK coarse-graining with long-range Coulomb
-* :doc:`lj/sdk/coul/msm <pair_sdk>` - LJ for SDK coarse-graining with long-range Coulomb via MSM
+* :doc:`lj/spica <pair_spica>` - LJ for SPICA coarse-graining
+* :doc:`lj/spica/coul/long <pair_spica>` - LJ for SPICA coarse-graining with long-range Coulomb
+* :doc:`lj/spica/coul/msm <pair_spica>` - LJ for SPICA coarse-graining with long-range Coulomb via MSM
 * :doc:`lj/sf/dipole/sf <pair_dipole>` - LJ with dipole interaction with shifted forces
 * :doc:`lj/smooth <pair_lj_smooth>` - smoothed Lennard-Jones potential
 * :doc:`lj/smooth/linear <pair_lj_smooth_linear>` - linear smoothed LJ potential
@@ -347,6 +349,7 @@ accelerated styles exist.
 * :doc:`spin/magelec <pair_spin_magelec>` -
 * :doc:`spin/neel <pair_spin_neel>` -
 * :doc:`srp <pair_srp>` -
+* :doc:`srp/react <pair_srp>` -
 * :doc:`sw <pair_sw>` - Stillinger-Weber 3-body potential
 * :doc:`sw/angle/table <pair_sw_angle_table>` - Stillinger-Weber potential with tabulated angular term
 * :doc:`sw/mod <pair_sw>` - modified Stillinger-Weber 3-body potential

doc/src/pair_sw.rst

@@ -24,13 +24,16 @@ Syntax
    pair_style style keyword values
 
 * style = *sw* or *sw/mod*
-* keyword = *maxdelcs*
+* keyword = *maxdelcs* or *threebody*
 
   .. parsed-literal::
 
-       *maxdelcs* value = delta1 delta2 (optional)
+       *maxdelcs* value = delta1 delta2 (optional, sw/mod only)
          delta1 = The minimum thershold for the variation of cosine of three-body angle
          delta2 = The maximum threshold for the variation of cosine of three-body angle
+       *threebody* value = *on* or *off* (optional, sw only)
+         on (default) = Compute both the three-body and two-body terms of the potential
+         off = Compute only the two-body term of the potential
 
 Examples
 """"""""
@@ -44,6 +47,11 @@ Examples
    pair_style sw/mod maxdelcs 0.25 0.35
    pair_coeff * * tmd.sw.mod Mo S S
 
+   pair_style hybrid sw threebody on sw threebody off
+   pair_coeff * * sw 1 mW_xL.sw mW NULL
+   pair_coeff 1 2 sw 2 mW_xL.sw mW xL
+   pair_coeff 2 2 sw 2 mW_xL.sw mW xL
+
 Description
 """""""""""
 
@@ -67,23 +75,28 @@ where :math:`\phi_2` is a two-body term and :math:`\phi_3` is a
 three-body term.  The summations in the formula are over all neighbors J
 and K of atom I within a cutoff distance :math:`a `\sigma`.
 
+.. versionadded:: 14Dec2021
+
 The *sw/mod* style is designed for simulations of materials when
-distinguishing three-body angles are necessary, such as borophene
-and transition metal dichalcogenides, which cannot be described
-by the original code for the Stillinger-Weber potential.
-For instance, there are several types of angles around each Mo atom in `MoS_2`,
-and some unnecessary angle types should be excluded in the three-body interaction.
-Such exclusion may be realized by selecting proper angle types directly.
-The exclusion of unnecessary angles is achieved here by the cut-off function (`f_C(\delta)`),
-which induces only minimum modifications for LAMMPS.
+distinguishing three-body angles are necessary, such as borophene and
+transition metal dichalcogenides, which cannot be described by the
+original code for the Stillinger-Weber potential.  For instance, there
+are several types of angles around each Mo atom in `MoS_2`, and some
+unnecessary angle types should be excluded in the three-body
+interaction.  Such exclusion may be realized by selecting proper angle
+types directly.  The exclusion of unnecessary angles is achieved here by
+the cut-off function (`f_C(\delta)`), which induces only minimum
+modifications for LAMMPS.
 
 Validation, benchmark tests, and applications of the *sw/mod* style
 can be found in :ref:`(Jiang2) <Jiang2>` and :ref:`(Jiang3) <Jiang3>`.
 
-The *sw/mod* style computes the energy E of a system of atoms, whose potential
-function is mostly the same as the Stillinger-Weber potential. The only modification
-is in the three-body term, where the value of :math:`\delta = \cos \theta_{ijk} - \cos \theta_{0ijk}`
-used in the original energy and force expression is scaled by a switching factor :math:`f_C(\delta)`:
+The *sw/mod* style computes the energy E of a system of atoms, whose
+potential function is mostly the same as the Stillinger-Weber
+potential. The only modification is in the three-body term, where the
+value of :math:`\delta = \cos \theta_{ijk} - \cos \theta_{0ijk}` used in
+the original energy and force expression is scaled by a switching factor
+:math:`f_C(\delta)`:
 
 .. math::
 
@@ -94,28 +107,46 @@ used in the original energy and force expression is scaled by a switching factor
     0 & \left| \delta \right| > \delta_2
     \end{array} \right. \\
 
-This cut-off function decreases smoothly from 1 to 0 over the range :math:`[\delta_1, \delta_2]`.
-This smoothly turns off the energy and force contributions for :math:`\left| \delta \right| > \delta_2`.
-It is suggested that :math:`\delta 1` and :math:`\delta_2` to be the value around
-:math:`0.5 \left| \cos \theta_1 - \cos \theta_2 \right|`, with
-:math:`\theta_1` and :math:`\theta_2` as the different types of angles around an atom.
-For borophene and transition metal dichalcogenides, :math:`\delta_1 = 0.25` and :math:`\delta_2 = 0.35`.
-This value enables the cut-off function to exclude unnecessary angles in the three-body SW terms.
+This cut-off function decreases smoothly from 1 to 0 over the range
+:math:`[\delta_1, \delta_2]`.  This smoothly turns off the energy and
+force contributions for :math:`\left| \delta \right| > \delta_2`.  It is
+suggested that :math:`\delta 1` and :math:`\delta_2` to be the value
+around :math:`0.5 \left| \cos \theta_1 - \cos \theta_2 \right|`, with
+:math:`\theta_1` and :math:`\theta_2` as the different types of angles
+around an atom.  For borophene and transition metal dichalcogenides,
+:math:`\delta_1 = 0.25` and :math:`\delta_2 = 0.35`.  This value enables
+the cut-off function to exclude unnecessary angles in the three-body SW
+terms.
 
 .. note::
 
-   The cut-off function is just to be used as a technique to exclude some unnecessary angles,
-   and it has no physical meaning. It should be noted that the force and potential are inconsistent
-   with each other in the decaying range of the cut-off function, as the angle dependence for the
-   cut-off function is not implemented in the force (first derivation of potential).
-   However, the angle variation is much smaller than the given threshold value for actual simulations,
-   so the inconsistency between potential and force can be neglected in actual simulations.
+   The cut-off function is just to be used as a technique to exclude
+   some unnecessary angles, and it has no physical meaning. It should be
+   noted that the force and potential are inconsistent with each other
+   in the decaying range of the cut-off function, as the angle
+   dependence for the cut-off function is not implemented in the force
+   (first derivation of potential).  However, the angle variation is
+   much smaller than the given threshold value for actual simulations,
+   so the inconsistency between potential and force can be neglected in
+   actual simulations.
 
-Only a single pair_coeff command is used with the *sw* and *sw/mod* styles
-which specifies a Stillinger-Weber potential file with parameters for all
-needed elements.  These are mapped to LAMMPS atom types by specifying
-N additional arguments after the filename in the pair_coeff command,
-where N is the number of LAMMPS atom types:
+The *threebody* keyword is optional and determines whether or not the
+three-body term of the potential is calculated.  The default value is
+"on" and it is only available for the plain *sw* pair style variants,
+but not available for the *sw/mod* and :doc:`sw/angle/table
+<pair_sw_angle_table>` pair style variants.  To turn off the threebody
+contributions all :math:`\lambda_{ijk}` parameters from the potential
+file are forcibly set to 0.  In addition the pair style implementation
+may employ code optimizations for the *threebody off* setting that can
+result in significant speedups versus the default.  These code optimizations
+are currently only available for the MANYBODY and OPENMP packages.
+
+Only a single pair_coeff command is used with the *sw* and *sw/mod*
+styles which specifies a Stillinger-Weber potential file with parameters
+for all needed elements, except for when the *threebody off* setting is
+used (see note below).  These are mapped to LAMMPS atom types by
+specifying N additional arguments after the filename in the pair_coeff
+command, where N is the number of LAMMPS atom types:
 
 * filename
 * N element names = mapping of SW elements to atom types
@@ -130,16 +161,22 @@ pair_coeff command:
 
 .. code-block:: LAMMPS
 
+   pair_style sw
    pair_coeff * * SiC.sw Si Si Si C
 
 The first 2 arguments must be \* \* so as to span all LAMMPS atom types.
-The first three Si arguments map LAMMPS atom types 1,2,3 to the Si
-element in the SW file.  The final C argument maps LAMMPS atom type 4
-to the C element in the SW file.  If a mapping value is specified as
-NULL, the mapping is not performed.  This can be used when a *sw*
+The first three Si arguments map LAMMPS atom types 1, 2, and 3 to the Si
+element in the SW file.  The final C argument maps LAMMPS atom type 4 to
+the C element in the SW file.  If an argument value is specified as
+NULL, the mapping is not performed.  This can be used when an *sw*
 potential is used as part of the *hybrid* pair style.  The NULL values
-are placeholders for atom types that will be used with other
-potentials.
+are placeholders for atom types that will be used with other potentials.
+
+.. note::
+
+   When the *threebody off* keyword is used, multiple pair_coeff commands may
+   be used to specific the pairs of atoms which don't require three-body term.
+   In these cases, the first 2 arguments are not required to be \* \*.
 
 Stillinger-Weber files in the *potentials* directory of the LAMMPS
 distribution have a ".sw" suffix.  Lines that are not blank or
@@ -243,30 +280,39 @@ 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 files.  Thus, you
-need to re-specify the pair_style and pair_coeff commands in an input
-script that reads a restart file.
+This pair style does not write its information to :doc:`binary restart
+files <restart>`, since it is stored in potential files.  Thus, you need
+to re-specify the pair_style and pair_coeff commands in an input script
+that reads a restart file.
 
 This pair style can only be used via the *pair* keyword of the
 :doc:`run_style respa <run_style>` command.  It does not support the
 *inner*, *middle*, *outer* keywords.
 
+The single() function of the *sw* pair style is only enabled and
+supported for the case of the *threebody off* setting.
+
 ----------
 
 Restrictions
 """"""""""""
 
-This pair style is part of the MANYBODY 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 is part of the MANYBODY 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"
 for pair interactions.
 
 The Stillinger-Weber potential files provided with LAMMPS (see the
 potentials directory) are parameterized for metal :doc:`units <units>`.
-You can use the SW potential with any LAMMPS units, but you would need
-to create your own SW potential file with coefficients listed in the
-appropriate units if your simulation does not use "metal" units.
+You can use the sw or sw/mod pair styles with any LAMMPS units, but you
+would need to create your own SW potential file with coefficients listed
+in the appropriate units if your simulation does not use "metal" units.
+If the potential file contains a 'UNITS:' metadata tag in the first line
+of the potential file, then LAMMPS can convert it transparently between
+"metal" and "real" units.
+
 
 Related commands
 """"""""""""""""
@@ -276,8 +322,9 @@ Related commands
 Default
 """""""
 
-The default values for the *maxdelcs* setting of the *sw/mod* pair
-style are *delta1* = 0.25 and *delta2* = 0.35`.
+The default value for the *threebody* setting of the "sw" pair style is
+"on", the default values for the "*maxdelcs* setting of the *sw/mod*
+pair style are *delta1* = 0.25 and *delta2* = 0.35`.
 
 ----------
 

doc/src/pair_sw_angle_table.rst

@@ -31,6 +31,8 @@ Used in example input script:
 Description
 """""""""""
 
+.. versionadded:: 2Jun2022
+
 The *sw/angle/table* style is a modification of the original
 :doc:`pair_style sw <pair_sw>`. It has been developed for coarse-grained
 simulations (of water) (:ref:`Scherer1 <Scherer1>`), but can be employed
@@ -296,7 +298,8 @@ for pair interactions.
 Related commands
 """"""""""""""""
 
-:doc:`pair_coeff <pair_coeff>`, :doc:`pair_style sw <pair_sw>`, :doc:`pair_style threebody/table <pair_threebody_table>`
+:doc:`pair_coeff <pair_coeff>`, :doc:`pair_style sw <pair_sw>`,
+:doc:`pair_style threebody/table <pair_threebody_table>`
 
 
 ----------

doc/src/pair_threebody_table.rst

@@ -35,6 +35,8 @@ Used in example input scripts:
 Description
 """""""""""
 
+.. versionadded:: 2Jun2022
+
 The *threebody/table* style is a pair style for generic tabulated
 three-body interactions.  It has been developed for (coarse-grained)
 simulations (of water) with Kernel-based machine learning (ML)

doc/src/read_dump.rst

@@ -24,12 +24,13 @@ Syntax
        *fx*,\ *fy*,\ *fz* = force components
 
 * zero or more keyword/value pairs may be appended
-* keyword = *nfile* or *box* or *replace* or *purge* or *trim* or *add* or *label* or *scaled* or *wrapped* or *format*
+* keyword = *nfile* or *box* or *timestep* or *replace* or *purge* or *trim* or *add* or *label* or *scaled* or *wrapped* or *format*
 
   .. parsed-literal::
 
        *nfile* value = Nfiles = how many parallel dump files exist
        *box* value = *yes* or *no* = replace simulation box with dump box
+       *timestep* value = *yes* or *no* = reset simulation timestep with dump timestep
        *replace* value = *yes* or *no* = overwrite atoms with dump atoms
        *purge* value = *yes* or *no* = delete all atoms before adding dump atoms
        *trim* value = *yes* or *no* = trim atoms not in dump snapshot
@@ -60,6 +61,7 @@ Examples
    read_dump dump.dcd 0 x y z box yes format molfile dcd
    read_dump dump.file 1000 x y z vx vy vz box yes format molfile lammpstrj /usr/local/lib/vmd/plugins/LINUXAMD64/plugins/molfile
    read_dump dump.file 5000 x y vx vy trim yes
+   read_dump dump.file 5000 x y vx vy add yes box no timestep no
    read_dump ../run7/dump.file.gz 10000 x y z box yes
    read_dump dump.xyz 10 x y z box no format molfile xyz ../plugins
    read_dump dump.dcd 0 x y z format molfile dcd
@@ -71,9 +73,9 @@ Description
 """""""""""
 
 Read atom information from a dump file to overwrite the current atom
-coordinates, and optionally the atom velocities and image flags and
-the simulation box dimensions.  This is useful for restarting a run
-from a particular snapshot in a dump file.  See the
+coordinates, and optionally the atom velocities and image flags, the
+simulation timestep, and the simulation box dimensions.  This is useful
+for restarting a run from a particular snapshot in a dump file.  See the
 :doc:`read_restart <read_restart>` and :doc:`read_data <read_data>`
 commands for alternative methods to do this.  Also see the
 :doc:`rerun <rerun>` command for a means of reading multiple snapshots
@@ -89,9 +91,9 @@ Also note that reading per-atom information from a dump snapshot is
 limited to the atom coordinates, velocities and image flags, as
 explained below.  Other atom properties, which may be necessary to run
 a valid simulation, such as atom charge, or bond topology information
-for a molecular system, are not read from (or even contained in) dump
-files.  Thus this auxiliary information should be defined in the usual
-way, e.g. in a data file read in by a :doc:`read_data <read_data>`
+for a molecular system, are not read from (or may not even be contained
+in) dump files.  Thus this auxiliary information should be defined in
+the usual way, e.g. in a data file read in by a :doc:`read_data <read_data>`
 command, before using the read_dump command, or by the :doc:`set <set>`
 command, after the dump snapshot is read.
 
@@ -165,11 +167,10 @@ variable *ntimestep*:
     uint64_t  ntimestep  5*scalar
       (0)    0 50 100 150 200
 
-Note that the *xyz*
-and *molfile* formats do not store the timestep.  For these formats,
-timesteps are numbered logically, in a sequential manner, starting
-from 0.  Thus to access the 10th snapshot in an *xyz* or *mofile*
-formatted dump file, use *Nstep* = 9.
+Note that the *xyz* and *molfile* formats do not store the timestep.
+For these formats, timesteps are numbered logically, in a sequential
+manner, starting from 0.  Thus to access the 10th snapshot in an *xyz*
+or *mofile* formatted dump file, use *Nstep* = 9.
 
 The dimensions of the simulation box for the selected snapshot are
 also read; see the *box* keyword discussion below.  For the *native*
@@ -266,8 +267,10 @@ for how this is done, determined by the specified fields and optional
 keywords.
 
 The timestep of the snapshot becomes the current timestep for the
-simulation.  See the :doc:`reset_timestep <reset_timestep>` command if
-you wish to change this after the dump snapshot is read.
+simulation unless the *timestep* keyword is specified with a *no* value
+(default setting is *yes*).  See the :doc:`reset_timestep <reset_timestep>`
+command if you wish to change this to a different value after the dump
+snapshot is read.
 
 If the *box* keyword is specified with a *yes* value, then the current
 simulation box dimensions are replaced by the dump snapshot box
@@ -391,7 +394,7 @@ Related commands
 Default
 """""""
 
-The option defaults are box = yes, replace = yes, purge = no, trim =
-no, add = no, scaled = no, wrapped = yes, and format = native.
+The option defaults are box = yes, timestep = yes, replace = yes, purge = no,
+trim = no, add = no, scaled = no, wrapped = yes, and format = native.
 
 .. _vmd: http://www.ks.uiuc.edu/Research/vmd

doc/src/special_bonds.rst

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

doc/src/thermo_modify.rst

@@ -153,6 +153,8 @@ containing the timestep and CPU time ("multi"), or in a YAML format
 block ("yaml").  This modify option overrides the *one*, *multi*, or
 *yaml* thermo_style settings.
 
+.. versionadded:: 4May2022
+
 The *colname* keyword can be used to change the default header keyword
 for a column or field of thermodynamic output.  The setting for *ID
 string* replaces the default text with the provided string.  *ID* can be

doc/src/variable.rst

@@ -685,7 +685,9 @@ of a run, according to this formula:
 The run begins on startstep and ends on stopstep.  Startstep and
 stopstep can span multiple runs, using the *start* and *stop* keywords
 of the :doc:`run <run>` command.  See the :doc:`run <run>` command for
-details of how to do this.
+details of how to do this.  If called in between runs or during a
+:doc:`run 0 <run>` command, the ramp(x,y) function will return the
+value of x.
 
 The stagger(x,y) function uses the current timestep to generate a new
 timestep.  X,y > 0 and x > y are required.  The generated timesteps
@@ -781,10 +783,14 @@ according to this formula:
 where dt = the timestep size.
 
 The run begins on startstep.  Startstep can span multiple runs, using
-the *start* keyword of the :doc:`run <run>` command.  See the
-:doc:`run <run>` command for details of how to do this.  Note that the
-:doc:`thermo_style <thermo_style>` keyword elaplong =
-timestep-startstep.
+the *start* keyword of the :doc:`run <run>` command.  See the :doc:`run
+<run>` command for details of how to do this.  Note that the
+:doc:`thermo_style <thermo_style>` keyword elaplong = timestep-startstep.
+If used between runs this function will return
+the value according to the end of the last run or the value of x if
+used before *any* runs.  This function assumes the length of the time
+step does not change and thus may not be used in combination with
+:doc:`fix dt/reset <fix_dt_reset>`.
 
 The swiggle(x,y,z) and cwiggle(x,y,z) functions each take 3 arguments:
 x = value0, y = amplitude, z = period.  They use the elapsed time to
@@ -799,10 +805,14 @@ run, according to one of these formulas, where omega = 2 PI / period:
 where dt = the timestep size.
 
 The run begins on startstep.  Startstep can span multiple runs, using
-the *start* keyword of the :doc:`run <run>` command.  See the
-:doc:`run <run>` command for details of how to do this.  Note that the
-:doc:`thermo_style <thermo_style>` keyword elaplong =
-timestep-startstep.
+the *start* keyword of the :doc:`run <run>` command.  See the :doc:`run
+<run>` command for details of how to do this.  Note that the
+:doc:`thermo_style <thermo_style>` keyword elaplong = timestep-startstep.
+If used between runs these functions will return
+the value according to the end of the last run or the value of x if
+used before *any* runs.  These functions assume the length of the time
+step does not change and thus may not be used in combination with
+:doc:`fix dt/reset <fix_dt_reset>`.
 
 ----------
 

doc/utils/check-styles.py

@@ -60,7 +60,7 @@ reader = {}
 region = {}
 total = 0
 
-index_pattern = re.compile(r"^.. index:: (compute|fix|pair_style|angle_style|bond_style|dihedral_style|improper_style|kspace_style)\s+([a-zA-Z0-9/_]+)$")
+index_pattern = re.compile(r"^.. index:: (compute|fix|pair_style|angle_style|bond_style|dihedral_style|improper_style|kspace_style|dump)\s+([a-zA-Z0-9/_]+)$")
 style_pattern = re.compile(r"(.+)Style\((.+),(.+)\)")
 upper = re.compile("[A-Z]+")
 gpu = re.compile("(.+)/gpu$")
@@ -84,7 +84,8 @@ def load_index_entries_in_file(path):
 
 def load_index_entries():
     index = {'compute': set(), 'fix': set(), 'pair_style': set(), 'angle_style': set(),
-             'bond_style': set(), 'dihedral_style': set(), 'improper_style': set(), 'kspace_style': set()}
+             'bond_style': set(), 'dihedral_style': set(), 'improper_style': set(),
+             'kspace_style': set(), 'dump': set()}
     rst_files = glob(os.path.join(doc_dir, '*.rst'))
     for f in rst_files:
         for command_type, style in load_index_entries_in_file(f):
@@ -254,8 +255,9 @@ for command_type, entries in index.items():
 
 print("Total number of style index entries:", total_index)
 
+skip_angle = ('sdk')
 skip_fix = ('python', 'NEIGH_HISTORY/omp','acks2/reax','qeq/reax','reax/c/bonds','reax/c/species')
-skip_pair = ('meam/c','lj/sf','reax/c')
+skip_pair = ('meam/c','lj/sf','reax/c','lj/sdk','lj/sdk/coul/long','lj/sdk/coul/msm')
 skip_compute = ('pressure/cylinder')
 
 counter = 0
@@ -269,13 +271,14 @@ counter += check_style('Commands_pair.rst', doc_dir, ":doc:`(.+) <pair.+>`",pair
 counter += check_style('pair_style.rst', doc_dir, ":doc:`(.+) <pair.+>` -",pair,'Pair',skip=skip_pair,suffix=False)
 counter += check_style('Commands_bond.rst', doc_dir, ":doc:`(.+) <bond.+>`",bond,'Bond',suffix=True)
 counter += check_style('bond_style.rst', doc_dir, ":doc:`(.+) <bond.+>` -",bond,'Bond',suffix=False)
-counter += check_style('Commands_bond.rst', doc_dir, ":doc:`(.+) <angle.+>`",angle,'Angle',suffix=True)
-counter += check_style('angle_style.rst', doc_dir, ":doc:`(.+) <angle.+>` -",angle,'Angle',suffix=False)
+counter += check_style('Commands_bond.rst', doc_dir, ":doc:`(.+) <angle.+>`",angle,'Angle',skip=skip_angle,suffix=True)
+counter += check_style('angle_style.rst', doc_dir, ":doc:`(.+) <angle.+>` -",angle,'Angle',skip=skip_angle,suffix=False)
 counter += check_style('Commands_bond.rst', doc_dir, ":doc:`(.+) <dihedral.+>`",dihedral,'Dihedral',suffix=True)
 counter += check_style('dihedral_style.rst', doc_dir, ":doc:`(.+) <dihedral.+>` -",dihedral,'Dihedral',suffix=False)
 counter += check_style('Commands_bond.rst', doc_dir, ":doc:`(.+) <improper.+>`",improper,'Improper',suffix=True)
 counter += check_style('improper_style.rst', doc_dir, ":doc:`(.+) <improper.+>` -",improper,'Improper',suffix=False)
 counter += check_style('Commands_kspace.rst', doc_dir, ":doc:`(.+) <kspace_style>`",kspace,'KSpace',suffix=True)
+counter += check_style('Commands_dump.rst', doc_dir, ":doc:`(.+) <dump.*>`",dump,'Dump',suffix=True)
 
 if counter:
     print(f"Found {counter} issue(s) with style lists")
@@ -284,12 +287,13 @@ counter = 0
 
 counter += check_style_index("compute", compute, index["compute"], skip=['pressure/cylinder'])
 counter += check_style_index("fix", fix, index["fix"], skip=['python','acks2/reax','qeq/reax','reax/c/bonds','reax/c/species'])
-counter += check_style_index("angle_style", angle, index["angle_style"])
+counter += check_style_index("angle_style", angle, index["angle_style"], skip=['sdk'])
 counter += check_style_index("bond_style", bond, index["bond_style"])
 counter += check_style_index("dihedral_style", dihedral, index["dihedral_style"])
 counter += check_style_index("improper_style", improper, index["improper_style"])
 counter += check_style_index("kspace_style", kspace, index["kspace_style"])
-counter += check_style_index("pair_style", pair, index["pair_style"], skip=['meam/c', 'lj/sf','reax/c'])
+counter += check_style_index("dump", dump, index["dump"])
+counter += check_style_index("pair_style", pair, index["pair_style"], skip=['meam/c','lj/sf','reax/c','lj/sdk','lj/sdk/coul/long','lj/sdk/coul/msm'])
 
 if counter:
     print(f"Found {counter} issue(s) with style index")

doc/utils/sphinx-config/_themes/lammps_theme/.gitignore

@@ -0,0 +1 @@
+!*.html