Merge branch 'fix-mdi-aimd-enhance' of github.com:lammps/lammps into fix-mdi-aimd-enhance

.github/dependabot.yml

@@ -0,0 +1,6 @@
+version: 2
+updates:
+  - package-ecosystem: "github-actions"
+    directory: "/"
+    schedule:
+      interval: "weekly"

.github/workflows/codeql-analysis.yml

@@ -13,6 +13,11 @@ jobs:
     if: ${{ github.repository == 'lammps/lammps' }}
     runs-on: ubuntu-latest
 
+    permissions:
+      security-events: write
+      actions: read
+      contents: read
+
     strategy:
       fail-fast: false
       matrix:
@@ -20,17 +25,17 @@ jobs:
 
     steps:
     - name: Checkout repository
-      uses: actions/checkout@v2
+      uses: actions/checkout@v3
       with:
         fetch-depth: 2
 
     - name: Setup Python
-      uses: actions/setup-python@v2
+      uses: actions/setup-python@v4
       with:
         python-version: '3.x'
 
     - name: Initialize CodeQL
-      uses: github/codeql-action/init@v1
+      uses: github/codeql-action/init@v2
       with:
         languages: ${{ matrix.language }}
         config-file: ./.github/codeql/${{ matrix.language }}.yml
@@ -48,4 +53,4 @@ jobs:
         cmake --build . --parallel 2
 
     - name: Perform CodeQL Analysis
-      uses: github/codeql-action/analyze@v1
+      uses: github/codeql-action/analyze@v2

.github/workflows/compile-msvc.yml

@@ -15,12 +15,12 @@ jobs:
 
     steps:
     - name: Checkout repository
-      uses: actions/checkout@v2
+      uses: actions/checkout@v3
       with:
         fetch-depth: 2
 
     - name: Select Python version
-      uses: actions/setup-python@v2
+      uses: actions/setup-python@v4
       with:
         python-version: '3.10'
 

.github/workflows/unittest-macos.yml

@@ -17,7 +17,7 @@ jobs:
 
     steps:
     - name: Checkout repository
-      uses: actions/checkout@v2
+      uses: actions/checkout@v3
       with:
         fetch-depth: 2
 
@@ -28,7 +28,7 @@ jobs:
       run: mkdir build
 
     - name: Set up ccache
-      uses: actions/cache@v2
+      uses: actions/cache@v3
       with:
         path: ${{ env.CCACHE_DIR }}
         key: macos-ccache-${{ github.sha }}

cmake/CMakeLists.txt

@@ -135,11 +135,13 @@ set(CMAKE_CXX_EXTENSIONS OFF CACHE BOOL "Use compiler extensions")
 # ugly hacks for MSVC which by default always reports an old C++ standard in the __cplusplus macro
 # and prints lots of pointless warnings about "unsafe" functions
 if(MSVC)
-  add_compile_options(/Zc:__cplusplus)
+  if(CMAKE_CXX_COMPILER_ID STREQUAL "MSVC")
-  add_compile_options(/wd4244)
+    add_compile_options(/Zc:__cplusplus)
-  add_compile_options(/wd4267)
+    add_compile_options(/wd4244)
-  if(LAMMPS_EXCEPTIONS)
+    add_compile_options(/wd4267)
-    add_compile_options(/EHsc)
+    if(LAMMPS_EXCEPTIONS)
+      add_compile_options(/EHsc)
+    endif()
   endif()
   add_compile_definitions(_CRT_SECURE_NO_WARNINGS)
 endif()
@@ -784,14 +786,16 @@ if(BUILD_SHARED_LIBS)
     find_package(Python COMPONENTS Interpreter)
   endif()
   if(BUILD_IS_MULTI_CONFIG)
-    set(LIBLAMMPS_SHARED_BINARY ${CMAKE_BINARY_DIR}/$<CONFIG>/liblammps${LAMMPS_MACHINE}${CMAKE_SHARED_LIBRARY_SUFFIX})
+    set(MY_BUILD_DIR ${CMAKE_BINARY_DIR}/$<CONFIG>)
   else()
-    set(LIBLAMMPS_SHARED_BINARY ${CMAKE_BINARY_DIR}/liblammps${LAMMPS_MACHINE}${CMAKE_SHARED_LIBRARY_SUFFIX})
+    set(MY_BUILD_DIR ${CMAKE_BINARY_DIR})
   endif()
+  set(LIBLAMMPS_SHARED_BINARY ${MY_BUILD_DIR}/liblammps${LAMMPS_MACHINE}${CMAKE_SHARED_LIBRARY_SUFFIX})
   if(Python_EXECUTABLE)
     add_custom_target(
       install-python ${CMAKE_COMMAND} -E remove_directory build
-      COMMAND ${Python_EXECUTABLE} ${LAMMPS_PYTHON_DIR}/install.py -p ${LAMMPS_PYTHON_DIR}/lammps -l ${LIBLAMMPS_SHARED_BINARY}
+      COMMAND ${Python_EXECUTABLE} ${LAMMPS_PYTHON_DIR}/install.py -p ${LAMMPS_PYTHON_DIR}/lammps
+              -l ${LIBLAMMPS_SHARED_BINARY} -w ${MY_BUILD_DIR}
       COMMENT "Installing LAMMPS Python module")
   else()
     add_custom_target(

cmake/Modules/LAMMPSInterfacePlugin.cmake

@@ -0,0 +1,197 @@
+# CMake script code to define LAMMPS settings required for building LAMMPS plugins
+
+# enforce out-of-source build
+if(${CMAKE_SOURCE_DIR} STREQUAL ${CMAKE_BINARY_DIR})
+  message(FATAL_ERROR "In-source builds are not allowed. You must create and use a build directory. "
+    "Please remove CMakeCache.txt and CMakeFiles first.")
+endif()
+
+# global LAMMPS/plugin build settings
+set(LAMMPS_SOURCE_DIR ""  CACHE PATH "Location of LAMMPS sources folder")
+if(NOT LAMMPS_SOURCE_DIR)
+  message(FATAL_ERROR "Must set LAMMPS_SOURCE_DIR")
+endif()
+
+# by default, install into $HOME/.local (not /usr/local),
+# so that no root access (and sudo) is needed
+if(CMAKE_INSTALL_PREFIX_INITIALIZED_TO_DEFAULT)
+  set(CMAKE_INSTALL_PREFIX "$ENV{HOME}/.local" CACHE PATH "Default install path" FORCE)
+endif()
+
+# ugly hacks for MSVC which by default always reports an old C++ standard in the __cplusplus macro
+# and prints lots of pointless warnings about "unsafe" functions
+if(MSVC)
+  if(CMAKE_CXX_COMPILER_ID STREQUAL "MSVC")
+    add_compile_options(/Zc:__cplusplus)
+    add_compile_options(/wd4244)
+    add_compile_options(/wd4267)
+    if(LAMMPS_EXCEPTIONS)
+      add_compile_options(/EHsc)
+    endif()
+  endif()
+  add_compile_definitions(_CRT_SECURE_NO_WARNINGS)
+endif()
+
+# C++11 is required
+set(CMAKE_CXX_STANDARD 11)
+set(CMAKE_CXX_STANDARD_REQUIRED ON)
+
+# Need -restrict with Intel compilers
+if(CMAKE_CXX_COMPILER_ID STREQUAL "Intel")
+  set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -restrict")
+endif()
+set(CMAKE_POSITION_INDEPENDENT_CODE TRUE)
+
+#######
+# helper functions from LAMMPSUtils.cmake
+function(validate_option name values)
+  string(TOLOWER ${${name}} needle_lower)
+  string(TOUPPER ${${name}} needle_upper)
+  list(FIND ${values} ${needle_lower} IDX_LOWER)
+  list(FIND ${values} ${needle_upper} IDX_UPPER)
+  if(${IDX_LOWER} LESS 0 AND ${IDX_UPPER} LESS 0)
+    list_to_bulletpoints(POSSIBLE_VALUE_LIST ${${values}})
+    message(FATAL_ERROR "\n########################################################################\n"
+      "Invalid value '${${name}}' for option ${name}\n"
+      "\n"
+      "Possible values are:\n"
+      "${POSSIBLE_VALUE_LIST}"
+      "########################################################################")
+  endif()
+endfunction(validate_option)
+
+# helper function for getting the most recently modified file or folder from a glob pattern
+function(get_newest_file path variable)
+  file(GLOB _dirs ${path})
+  set(_besttime 2000-01-01T00:00:00)
+  set(_bestfile "<unknown>")
+  foreach(_dir ${_dirs})
+    file(TIMESTAMP ${_dir} _newtime)
+    if(_newtime IS_NEWER_THAN _besttime)
+      set(_bestfile ${_dir})
+      set(_besttime ${_newtime})
+    endif()
+  endforeach()
+  if(_bestfile STREQUAL "<unknown>")
+    message(FATAL_ERROR "Could not find valid path at: ${path}")
+  endif()
+  set(${variable} ${_bestfile} PARENT_SCOPE)
+endfunction()
+
+#################################################################################
+# LAMMPS C++ interface. We only need the header related parts except on windows.
+add_library(lammps INTERFACE)
+target_include_directories(lammps INTERFACE ${LAMMPS_SOURCE_DIR})
+if((CMAKE_SYSTEM_NAME STREQUAL "Windows") AND CMAKE_CROSSCOMPILING)
+  target_link_libraries(lammps INTERFACE ${CMAKE_BINARY_DIR}/../liblammps.dll.a)
+endif()
+
+################################################################################
+# MPI configuration
+if(NOT CMAKE_CROSSCOMPILING)
+  find_package(MPI QUIET)
+  option(BUILD_MPI "Build MPI version" ${MPI_FOUND})
+else()
+  option(BUILD_MPI "Build MPI version" OFF)
+endif()
+
+if(BUILD_MPI)
+  # do not include the (obsolete) MPI C++ bindings which makes
+  # for leaner object files and avoids namespace conflicts
+  set(MPI_CXX_SKIP_MPICXX TRUE)
+  # We use a non-standard procedure to cross-compile with MPI on Windows
+  if((CMAKE_SYSTEM_NAME STREQUAL "Windows") AND CMAKE_CROSSCOMPILING)
+    # Download and configure custom MPICH files for Windows
+    message(STATUS "Downloading and configuring MPICH-1.4.1 for Windows")
+    set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win64-devel.tar.gz" CACHE STRING "URL for MPICH2 (win64) tarball")
+    set(MPICH2_WIN32_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win32-devel.tar.gz" CACHE STRING "URL for MPICH2 (win32) tarball")
+    set(MPICH2_WIN64_DEVEL_MD5 "4939fdb59d13182fd5dd65211e469f14" CACHE STRING "MD5 checksum of MPICH2 (win64) tarball")
+    set(MPICH2_WIN32_DEVEL_MD5 "a61d153500dce44e21b755ee7257e031" CACHE STRING "MD5 checksum of MPICH2 (win32) tarball")
+    mark_as_advanced(MPICH2_WIN64_DEVEL_URL)
+    mark_as_advanced(MPICH2_WIN32_DEVEL_URL)
+    mark_as_advanced(MPICH2_WIN64_DEVEL_MD5)
+    mark_as_advanced(MPICH2_WIN32_DEVEL_MD5)
+
+    include(ExternalProject)
+    if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64")
+      ExternalProject_Add(mpi4win_build
+        URL     ${MPICH2_WIN64_DEVEL_URL}
+        URL_MD5 ${MPICH2_WIN64_DEVEL_MD5}
+        CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND ""
+        BUILD_BYPRODUCTS <SOURCE_DIR>/lib/libmpi.a)
+    else()
+      ExternalProject_Add(mpi4win_build
+        URL     ${MPICH2_WIN32_DEVEL_URL}
+        URL_MD5 ${MPICH2_WIN32_DEVEL_MD5}
+        CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND ""
+        BUILD_BYPRODUCTS <SOURCE_DIR>/lib/libmpi.a)
+    endif()
+
+    ExternalProject_get_property(mpi4win_build SOURCE_DIR)
+    file(MAKE_DIRECTORY "${SOURCE_DIR}/include")
+    add_library(MPI::MPI_CXX UNKNOWN IMPORTED)
+    set_target_properties(MPI::MPI_CXX PROPERTIES
+      IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmpi.a"
+      INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include"
+      INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX")
+    add_dependencies(MPI::MPI_CXX mpi4win_build)
+
+    # set variables for status reporting at the end of CMake run
+    set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include")
+    set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX")
+    set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmpi.a")
+  else()
+    find_package(MPI REQUIRED)
+    option(LAMMPS_LONGLONG_TO_LONG "Workaround if your system or MPI version does not recognize 'long long' data types" OFF)
+    if(LAMMPS_LONGLONG_TO_LONG)
+      target_compile_definitions(lammps INTERFACE -DLAMMPS_LONGLONG_TO_LONG)
+    endif()
+  endif()
+  target_link_libraries(lammps INTERFACE MPI::MPI_CXX)
+else()
+  add_library(mpi_stubs INTERFACE)
+  target_include_directories(mpi_stubs INTERFACE $<BUILD_INTERFACE:${LAMMPS_SOURCE_DIR}/STUBS>)
+  target_link_libraries(lammps INTERFACE mpi_stubs)
+endif()
+
+################################################################################
+# detect if we may enable OpenMP support by default
+set(BUILD_OMP_DEFAULT OFF)
+find_package(OpenMP QUIET)
+if(OpenMP_FOUND)
+  check_include_file_cxx(omp.h HAVE_OMP_H_INCLUDE)
+  if(HAVE_OMP_H_INCLUDE)
+    set(BUILD_OMP_DEFAULT ON)
+  endif()
+endif()
+
+option(BUILD_OMP "Build with OpenMP support" ${BUILD_OMP_DEFAULT})
+
+if(BUILD_OMP)
+  find_package(OpenMP REQUIRED)
+  check_include_file_cxx(omp.h HAVE_OMP_H_INCLUDE)
+  if(NOT HAVE_OMP_H_INCLUDE)
+    message(FATAL_ERROR "Cannot find the 'omp.h' header file required for full OpenMP support")
+  endif()
+
+  if (((CMAKE_CXX_COMPILER_ID STREQUAL "GNU") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 9.0)) OR
+      (CMAKE_CXX_COMPILER_ID STREQUAL "PGI") OR
+      ((CMAKE_CXX_COMPILER_ID STREQUAL "Clang") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 10.0)) OR
+      ((CMAKE_CXX_COMPILER_ID STREQUAL "Intel") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 19.0)))
+    # GCC 9.x and later plus Clang 10.x and later implement strict OpenMP 4.0 semantics for consts.
+    # Intel 18.0 was tested to support both, so we switch to OpenMP 4+ from 19.x onward to be safe.
+    target_compile_definitions(lammps INTERFACE -DLAMMPS_OMP_COMPAT=4)
+  else()
+    target_compile_definitions(lammps INTERFACE -DLAMMPS_OMP_COMPAT=3)
+  endif()
+  target_link_libraries(lammps INTERFACE OpenMP::OpenMP_CXX)
+endif()
+
+################
+# integer size selection
+set(LAMMPS_SIZES "smallbig" CACHE STRING "LAMMPS integer sizes (smallsmall: all 32-bit, smallbig: 64-bit #atoms #timesteps, bigbig: also 64-bit imageint, 64-bit atom ids)")
+set(LAMMPS_SIZES_VALUES smallbig bigbig smallsmall)
+set_property(CACHE LAMMPS_SIZES PROPERTY STRINGS ${LAMMPS_SIZES_VALUES})
+validate_option(LAMMPS_SIZES LAMMPS_SIZES_VALUES)
+string(TOUPPER ${LAMMPS_SIZES} LAMMPS_SIZES)
+target_compile_definitions(lammps INTERFACE -DLAMMPS_${LAMMPS_SIZES})

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_URL "https://github.com/MolSSI-MDI/MDI_Library/archive/v1.4.1.tar.gz" CACHE STRING "URL for MDI tarball")
-  set(MDI_MD5 "836f5da400d8cff0f0e4435640f9454f" CACHE STRING "MD5 checksum 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)

cmake/Modules/Packages/ML-PACE.cmake

@@ -32,5 +32,6 @@ target_include_directories(pace PUBLIC ${PACE_EVALUATOR_INCLUDE_DIR} ${YAML_CPP_
 
 
 target_link_libraries(pace PRIVATE yaml-cpp-pace)
-
+if(CMAKE_PROJECT_NAME STREQUAL "lammps")
-target_link_libraries(lammps PRIVATE pace)
+  target_link_libraries(lammps PRIVATE pace)
+endif()

cmake/presets/mingw-cross.cmake

@@ -46,8 +46,8 @@ set(WIN_PACKAGES
   MISC
   ML-HDNNP
   ML-IAP
-  ML-SNAP
   ML-RANN
+  ML-SNAP
   MOFFF
   MOLECULE
   MOLFILE
@@ -56,6 +56,7 @@ set(WIN_PACKAGES
   ORIENT
   PERI
   PHONON
+  PLUGIN
   POEMS
   PTM
   QEQ

doc/Makefile

@@ -13,7 +13,7 @@ VENV           = $(BUILDDIR)/docenv
 ANCHORCHECK    = $(VENV)/bin/rst_anchor_check
 SPHINXCONFIG   = $(BUILDDIR)/utils/sphinx-config
 MATHJAX        = $(SPHINXCONFIG)/_static/mathjax
-MATHJAXTAG     = 3.2.1
+MATHJAXTAG     = 3.2.2
 
 PYTHON         = $(word 3,$(shell type python3))
 DOXYGEN        = $(word 3,$(shell type doxygen))

doc/lammps.1

@@ -1,7 +1,7 @@
-.TH LAMMPS "1" "2 June 2022" "2022-6-2"
+.TH LAMMPS "1" "23 June 2022" "2022-6-23"
 .SH NAME
 .B LAMMPS
-\- Molecular Dynamics Simulator.  Version 2 June 2022
+\- Molecular Dynamics Simulator.  Version 23 June 2022
 
 .SH SYNOPSIS
 .B lmp

doc/src/Commands_compute.rst

@@ -138,6 +138,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/Developer_plugins.rst

@@ -276,10 +276,27 @@ Compilation of the plugin can be managed via both, CMake or traditional
 GNU makefiles.  Some examples that can be used as a template are in the
 ``examples/plugins`` folder.  The CMake script code has some small
 adjustments to allow building the plugins for running unit tests with
-them. Another example that converts the KIM package into a plugin can be
+them.
-found in the ``examples/kim/plugin`` folder.  No changes to the sources
+
-of the KIM package themselves are needed; only the plugin interface and
+Another example that converts the KIM package into a plugin can be found
-loader code needs to be added.  This example only supports building with
+in the ``examples/kim/plugin`` folder.  No changes to the sources of the
-CMake, but is probably a more typical example. To compile you need to
+KIM package themselves are needed; only the plugin interface and loader
-run CMake with -DLAMMPS_SOURCE_DIR=<path/to/lammps/src/folder>.  Other
+code needs to be added.  This example only supports building with CMake,
+but is probably a more typical example. To compile you need to run CMake
+with -DLAMMPS_SOURCE_DIR=<path/to/lammps/src/folder>.  Other
 configuration setting are identical to those for compiling LAMMPS.
+
+A second example for a plugin from a package is in the
+``examples/PACKAGES/pace/plugin`` folder that will create a plugin from
+the ML-PACE package.  In this case the bulk of the code is in a static
+external library that is being downloaded and compiled first and then
+combined with the pair style wrapper and the plugin loader.  This
+example also contains a NSIS script that can be used to create an
+Installer package for Windows (the mutual licensing terms of the
+external library and LAMMPS conflict when distributing binaries, so the
+ML-PACE package cannot be linked statically, but the LAMMPS headers
+required to build the plugin are also available under a less restrictive
+license).  This will automatically set the required environment variable
+and launching a (compatible) LAMMPS binary will load and register the
+plugin and the ML-PACE package can then be used as it was linked into
+LAMMPS.

doc/src/Errors_warnings.rst

@@ -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.
 

doc/src/Howto_mdi.rst

@@ -72,7 +72,7 @@ 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 evalute the quantum energy and forces for a series of
+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.
 

doc/src/Howto_structured_data.rst

@@ -184,7 +184,7 @@ frame.
 
 .. code-block:: python
 
-   import re, yaml
+   import yaml
    import pandas as pd
 
    try:
@@ -193,7 +193,7 @@ frame.
        from yaml import SafeLoader as Loader
 
    with open("ave.yaml") as f:
-       ave = yaml.load(docs, Loader=Loader)
+       ave = yaml.load(f, Loader=Loader)
 
    keys = ave['keywords']
    df = {}

doc/src/Packages_details.rst

@@ -657,7 +657,7 @@ advection-diffusion-reaction systems. The equations of motion of these
 DPD extensions are integrated through a modified velocity-Verlet (MVV)
 algorithm.
 
-**Author:** Zhen Li (Division of Applied Mathematics, Brown University)
+**Author:** Zhen Li (Department of Mechanical Engineering, Clemson University)
 
 **Supporting info:**
 
@@ -1479,7 +1479,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 +1801,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/Speed_measure.rst

@@ -42,5 +42,4 @@ inaccurate relative timing data, because processors have to wait when
 communication occurs for other processors to catch up.  Thus the
 reported times for "Communication" or "Other" may be higher than they
 really are, due to load-imbalance.  If this is an issue, you can
-uncomment the MPI_Barrier() lines in src/timer.cpp, and re-compile
+use the :doc:`timer sync <timer>` command to obtain synchronized timings.
-LAMMPS, to obtain synchronized timings.

doc/src/Tools.rst

@@ -95,7 +95,7 @@ Miscellaneous tools
    * :ref:`LAMMPS shell <lammps_shell>`
    * :ref:`LAMMPS magic patterns for file(1) <magic>`
    * :ref:`Offline build tool <offline>`
-   * :ref:`singularity <singularity_tool>`
+   * :ref:`singularity/apptainer <singularity_tool>`
    * :ref:`SWIG interface <swig>`
    * :ref:`vim <vim>`
 
@@ -1007,14 +1007,15 @@ Ivanov, at University of Iceland (ali5 at hi.is).
 
 .. _singularity_tool:
 
-singularity tool
+singularity/apptainer tool
-----------------------------------------
+--------------------------
 
-The singularity sub-directory contains container definitions files
+The singularity sub-directory contains container definitions files that
-that can be used to build container images for building and testing
+can be used to build container images for building and testing LAMMPS on
-LAMMPS on specific OS variants using the `Singularity <https://sylabs.io>`_
+specific OS variants using the `Apptainer <https://apptainer.org>`_ or
-container software. Contributions for additional variants are welcome.
+`Singularity <https://sylabs.io>`_ container software. Contributions for
-For more details please see the README.md file in that folder.
+additional variants are welcome.  For more details please see the
+README.md file in that folder.
 
 ----------
 

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_ave_sphere_atom.rst

@@ -35,16 +35,24 @@ Examples
 Description
 """""""""""
 
-Define a computation that calculates the local density and temperature
+Define a computation that calculates the local mass density and
-for each atom and neighbors inside a spherical cutoff.
+temperature for each atom based on its neighbors inside a spherical
+cutoff.  If an atom has M neighbors, then its local mass density is
+calculated as the sum of its mass and its M neighbor masses, divided
+by the volume of the cutoff sphere (or circle in 2d).  The local
+temperature of the atom is calculated as the temperature of the
+collection of M+1 atoms, after subtracting the center-of-mass velocity
+of the M+1 atoms from each of the M+1 atom's velocities.  This is
+effectively the thermal velocity of the neighborhood of the central
+atom, similar to :doc:`compute temp/com <compute_temp_com>`.
 
-The optional keyword *cutoff* defines the distance cutoff
+The optional keyword *cutoff* defines the distance cutoff used when
-used when searching for neighbors. The default value is the cutoff
+searching for neighbors. The default value is the cutoff specified by
-specified by the pair style. If no pair style is defined, then a cutoff
+the pair style. If no pair style is defined, then a cutoff must be
-must be defined using this keyword. If the specified cutoff is larger than
+defined using this keyword. If the specified cutoff is larger than
-that of the pair_style plus neighbor skin (or no pair style is defined),
+that of the pair_style plus neighbor skin (or no pair style is
-the *comm_modify cutoff* option must also be set to match that of the
+defined), the *comm_modify cutoff* option must also be set to match
-*cutoff* keyword.
+that of the *cutoff* keyword.
 
 The neighbor list needed to compute this quantity is constructed each
 time the calculation is performed (i.e. each time a snapshot of atoms
@@ -55,16 +63,16 @@ too frequently.
 
    If you have a bonded system, then the settings of
    :doc:`special_bonds <special_bonds>` command can remove pairwise
-   interactions between atoms in the same bond, angle, or dihedral.  This
+   interactions between atoms in the same bond, angle, or dihedral.
-   is the default setting for the :doc:`special_bonds <special_bonds>`
+   This is the default setting for the :doc:`special_bonds
-   command, and means those pairwise interactions do not appear in the
+   <special_bonds>` command, and means those pairwise interactions do
-   neighbor list.  Because this fix uses the neighbor list, it also means
+   not appear in the neighbor list.  Because this compute uses the
-   those pairs will not be included in the order parameter.  This
+   neighbor list, it also means those pairs will not be included in
-   difficulty can be circumvented by writing a dump file, and using the
+   the order parameter.  This difficulty can be circumvented by
-   :doc:`rerun <rerun>` command to compute the order parameter for
+   writing a dump file, and using the :doc:`rerun <rerun>` command to
-   snapshots in the dump file.  The rerun script can use a
+   compute the order parameter for snapshots in the dump file.  The
-   :doc:`special_bonds <special_bonds>` command that includes all pairs in
+   rerun script can use a :doc:`special_bonds <special_bonds>` command
-   the neighbor list.
+   that includes all pairs in the neighbor list.
 
 ----------
 
@@ -77,17 +85,20 @@ too frequently.
 Output info
 """""""""""
 
-This compute calculates a per-atom array with two columns: density and temperature.
+This compute calculates a per-atom array with two columns: mass
+density in density :doc:`units <units>` and temperature in temperature
+:doc:`units <units>`.
 
 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
+from a compute as input.  See the :doc:`Howto output <Howto_output>`
-page for an overview of LAMMPS output options.
+doc page for an overview of LAMMPS output options.
 
 Restrictions
 """"""""""""
 
-This compute is part of the EXTRA-COMPUTE package.  It is only enabled if
+This compute is part of the EXTRA-COMPUTE package.  It is only enabled
-LAMMPS was built with that package.  See the :doc:`Build package <Build_package>` page for more info.
+if LAMMPS was built with that package.  See the :doc:`Build package
+<Build_package>` page for more info.
 
 Related commands
 """"""""""""""""
@@ -97,5 +108,5 @@ Related commands
 Default
 """""""
 
-The option defaults are *cutoff* = pair style cutoff
+The option defaults are *cutoff* = pair style cutoff.
 

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,6 +43,7 @@ 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*
 
@@ -78,6 +90,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 +225,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.
@@ -414,6 +467,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*,
@@ -464,8 +532,7 @@ The optional keyword defaults are *rmin0* = 0,
 
 .. _Thompson20141:
 
-**(Thompson)** Thompson, Swiler, Trott, Foiles, Tucker, under review, preprint
+**(Thompson)** Thompson, Swiler, Trott, Foiles, Tucker, J Comp Phys, 285, 316, (2015).
-available at `arXiv:1409.3880 <http://arxiv.org/abs/1409.3880>`_
 
 .. _Bartok20101:
 
@@ -486,4 +553,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/fix_mdi_qm.rst

@@ -57,9 +57,9 @@ explained below.
 
 These are example use cases for this fix, discussed further below:
 
-* perform an ab intitio MD (AIMD) simulation with quantum forces
+* perform an ab initio MD (AIMD) simulation with quantum forces
-* perform an energy minimziation with quantum forces
+* perform an energy minimization with quantum forces
-* perform a nudged elatic band (NEB) calculation 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
@@ -71,7 +71,7 @@ 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 dicussed below.  In each case, two instances of
+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
@@ -141,7 +141,7 @@ 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.  
+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
@@ -157,7 +157,7 @@ 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 includesf :doc:`fix box relax <fix_box_relaxq>` to
+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.
@@ -172,7 +172,7 @@ 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 evaulation of energy and forces for a series of *N*
+(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
@@ -180,7 +180,7 @@ 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 <deletea_atoms>`, and/or :doc:`displace_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

doc/src/kspace_style.rst

@@ -129,8 +129,8 @@ Examples
 
    kspace_style pppm 1.0e-4
    kspace_style pppm/cg 1.0e-5 1.0e-6
-   kspace style msm 1.0e-4
+   kspace_style msm 1.0e-4
-   kspace style scafacos fmm 1.0e-4
+   kspace_style scafacos fmm 1.0e-4
    kspace_style none
 
 Used in input scripts:

doc/src/mdi.rst

@@ -326,15 +326,15 @@ 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
-issuedd multiple times, e.g. in a loop that issues a :doc:`clear
+issued multiple times, e.g. in a loop that issues a :doc:`clear
-<clear>` command, then it cannot initiate/terminate the connection
+<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 intiate/terminate the connection.
+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 enery and force evaluation.  A *clear* command inside the
+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

doc/src/pair_e3b.rst

@@ -50,6 +50,12 @@ Examples
    pair_style hybrid/overlay e3b 1 lj/cut/tip4p/long 1 2 1 1 0.15 8.5
    pair_coeff * * e3b preset 2011
 
+Used in example input script:
+
+.. parsed-literal::
+
+   examples/PACKAGES/e3b/in.e3b-tip4p2005
+
 Description
 """""""""""
 
@@ -68,21 +74,27 @@ The *e3b* style computes an \"explicit three-body\" (E3B) potential for water :r
    0 & r>R_f\\
    \end{cases}
 
-This potential was developed as a water model that includes the three-body cooperativity of hydrogen bonding explicitly.
+This potential was developed as a water model that includes the
-To use it in this way, it must be applied in conjunction with a conventional two-body water model, through *pair_style hybrid/overlay*.
+three-body cooperativity of hydrogen bonding explicitly.  To use it in
-The three body interactions are split into three types: A, B, and C.
+this way, it must be applied in conjunction with a conventional two-body
-Type A corresponds to anti-cooperative double hydrogen bond donor interactions.
+water model, through pair style :doc:`hybrid/overlay <pair_hybrid>`.  The
-Type B corresponds to the cooperative interaction of molecules that both donate and accept a hydrogen bond.
+three body interactions are split into three types: A, B, and C.  Type A
-Type C corresponds to anti-cooperative double hydrogen bond acceptor interactions.
+corresponds to anti-cooperative double hydrogen bond donor interactions.
-The three-body interactions are smoothly cutoff by the switching function s(r) between Rs and Rc3.
+Type B corresponds to the cooperative interaction of molecules that both
-The two-body interactions are designed to correct for the effective many-body interactions implicitly included in the conventional two-body potential.
+donate and accept a hydrogen bond.  Type C corresponds to
-The two-body interactions are cut off sharply at Rc2, because K3 is typically significantly smaller than K2.
+anti-cooperative double hydrogen bond acceptor interactions.  The
-See :ref:`(Kumar 2008) <Kumar>` for more details.
+three-body interactions are smoothly cutoff by the switching function
+s(r) between Rs and Rc3.  The two-body interactions are designed to
+correct for the effective many-body interactions implicitly included in
+the conventional two-body potential.  The two-body interactions are cut
+off sharply at Rc2, because K3 is typically significantly smaller than
+K2.  See :ref:`(Kumar 2008) <Kumar>` for more details.
 
-Only a single *pair_coeff* command is used with the *e3b* style.
+Only a single :doc:`pair_coeff <pair_coeff>` command is used with the
-The first two arguments must be \* \*.
+*e3b* style and the first two arguments must be \* \*.  The oxygen atom
-The oxygen atom type for the pair style is passed as the only argument to the *pair_style* command, not in the *pair_coeff* command.
+type for the pair style is passed as the only argument to the
-The hydrogen atom type is inferred by the ordering of the atoms.
+*pair_style* command, not in the *pair_coeff* command.  The hydrogen
+atom type is inferred from the ordering of the atoms.
 
 .. note::
 
@@ -90,26 +102,41 @@ The hydrogen atom type is inferred by the ordering of the atoms.
    Each water molecule must have consecutive IDs with the oxygen first.
    This pair style does not test that this criteria is met.
 
-The *pair_coeff* command must have at least one keyword/value pair, as described above.
+The *pair_coeff* command must have at least one keyword/value pair, as
-The *preset* keyword sets the potential parameters to the values used in :ref:`(Tainter 2011) <Tainter2011>` or :ref:`(Tainter 2015) <Tainter2015>`.
+described above.  The *preset* keyword sets the potential parameters to
-To use the water models defined in those references, the *e3b* style should always be used in conjunction with an *lj/cut/tip4p/long* style through *pair_style hybrid/overlay*, as demonstrated in the second example above.
+the values used in :ref:`(Tainter 2011) <Tainter2011>` or
-The *preset 2011* option should be used with the :doc:`TIP4P water model <Howto_tip4p>`.
+:ref:`(Tainter 2015) <Tainter2015>`.  To use the water models defined in
-The *preset 2015* option should be used with the :doc:`TIP4P/2005 water model <Howto_tip4p>`.
+those references, the *e3b* style should always be used in conjunction
-If the *preset* keyword is used, no other keyword is needed.
+with an *lj/cut/tip4p/long* style through *pair_style hybrid/overlay*,
-Changes to the preset parameters can be made by specifying the *preset* keyword followed by the specific parameter to change, like *Ea*\ .
+as demonstrated in the second example above.  The *preset 2011* option
-Note that the other keywords must come after *preset* in the pair_style command.
+should be used with the :doc:`TIP4P water model <Howto_tip4p>`.  The
-The *e3b* style can also be used to implement any three-body potential of the same form by specifying all the keywords except *neigh*\ : *Ea*, *Eb*, *Ec*, *E2*, *K3*, *K2*, *Rc3*, *Rc2*, *Rs*, and *bondL*\ .
+*preset 2015* option should be used with the :doc:`TIP4P/2005 water
-The keyword *bondL* specifies the intramolecular OH bond length of the water model being used.
+model <Howto_tip4p>`.  If the *preset* keyword is used, no other keyword
-This is needed to include H atoms that are within the cutoff even when the attached oxygen atom is not.
+is needed.  Changes to the preset parameters can be made by specifying
+the *preset* keyword followed by the specific parameter to change, like
+*Ea*\ .  Note that the other keywords must come after *preset* in the
+pair_style command.  The *e3b* style can also be used to implement any
+three-body potential of the same form by specifying all the keywords
+except *neigh*\ : *Ea*, *Eb*, *Ec*, *E2*, *K3*, *K2*, *Rc3*, *Rc2*,
+*Rs*, and *bondL*\ .  The keyword *bondL* specifies the intramolecular
+OH bond length of the water model being used.  This is needed to include
+H atoms that are within the cutoff even when the attached oxygen atom is
+not.
 
-This pair style allocates arrays sized according to the number of pairwise interactions within Rc3.
+This pair style allocates arrays sized according to the number of
-To do this it needs an estimate for the number of water molecules within Rc3 of an oxygen atom.
+pairwise interactions within Rc3.  To do this it needs an estimate for
-This estimate defaults to 10 and can be changed using the *neigh* keyword, which takes an integer as an argument.
+the number of water molecules within Rc3 of an oxygen atom.  This
-If the neigh setting is too small, the simulation will fail with the error "neigh is too small".
+estimate defaults to 10 and can be changed using the *neigh* keyword,
-If the neigh setting is too large, the pair style will use more memory than necessary.
+which takes an integer as an argument.  If the neigh setting is too
+small, the simulation will fail with the error "neigh is too small".  If
+the neigh setting is too large, the pair style will use more memory than
+necessary.
 
-This pair style tallies a breakdown of the total E3B potential energy into sub-categories, which can be accessed via the :doc:`compute pair <compute_pair>` command as a vector of values of length 4.
+This pair style tallies a breakdown of the total E3B potential energy
-The 4 values correspond to the terms in the first equation above: the E2 term, the Ea term, the Eb term, and the Ec term.
+into sub-categories, which can be accessed via the :doc:`compute pair
+<compute_pair>` command as a vector of values of length 4.  The 4 values
+correspond to the terms in the first equation above: the E2 term, the Ea
+term, the Eb term, and the Ec term.
 
 See the examples/PACKAGES/e3b directory for a complete example script.
 

doc/src/pair_sw_angle_table.rst

@@ -23,9 +23,9 @@ Examples
 
 Used in example input script:
 
-  .. parsed-literal::
+.. parsed-literal::
 
-     examples/PACKAGES/manybody_table/in.spce_sw
+   examples/PACKAGES/manybody_table/in.spce_sw
 
 
 Description

doc/src/pair_threebody_table.rst

@@ -27,10 +27,10 @@ Examples
 
 Used in example input scripts:
 
-  .. parsed-literal::
+.. parsed-literal::
 
-     examples/PACKAGES/manybody_table/in.spce
+   examples/PACKAGES/manybody_table/in.spce
-     examples/PACKAGES/manybody_table/in.spce2
+   examples/PACKAGES/manybody_table/in.spce2
 
 Description
 """""""""""

doc/src/read_data.rst

@@ -56,7 +56,7 @@ Examples
    read_data ../run7/data.polymer.gz
    read_data data.protein fix mycmap crossterm CMAP
    read_data data.water add append offset 3 1 1 1 1 shift 0.0 0.0 50.0
-   read_data data.water add merge 1 group solvent
+   read_data data.water add merge group solvent
 
 Description
 """""""""""
@@ -622,6 +622,8 @@ of analysis.
      - atom-ID molecule-ID atom-type x y z
    * - charge
      - atom-ID atom-type q x y z
+   * - dielectric
+     - atom-ID atom-type q x y z normx normy normz area ed em epsilon curvature
    * - dipole
      - atom-ID atom-type q x y z mux muy muz
    * - dpd

doc/src/read_restart.rst

@@ -11,7 +11,6 @@ Syntax
    read_restart file flag
 
 * file = name of binary restart file to read in
-* flag = remap (optional)
 
 Examples
 """"""""
@@ -19,44 +18,40 @@ Examples
 .. code-block:: LAMMPS
 
    read_restart save.10000
-   read_restart save.10000 remap
    read_restart restart.*
    read_restart restart.*.mpiio
-   read_restart poly.*.% remap
 
 Description
 """""""""""
 
 Read in a previously saved system configuration from a restart file.
 This allows continuation of a previous run.  Details about what
-information is stored (and not stored) in a restart file is given
+information is stored (and not stored) in a restart file is given below.
-below.  Basically this operation will re-create the simulation box
+Basically this operation will re-create the simulation box with all its
-with all its atoms and their attributes as well as some related global
+atoms and their attributes as well as some related global settings, at
-settings, at the point in time it was written to the restart file by a
+the point in time it was written to the restart file by a previous
-previous simulation.  The simulation box will be partitioned into a
+simulation.  The simulation box will be partitioned into a regular 3d
-regular 3d grid of rectangular bricks, one per processor, based on the
+grid of rectangular bricks, one per processor, based on the number of
-number of processors in the current simulation and the settings of the
+processors in the current simulation and the settings of the
 :doc:`processors <processors>` command.  The partitioning can later be
-changed by the :doc:`balance <balance>` or :doc:`fix balance <fix_balance>` commands.
+changed by the :doc:`balance <balance>` or :doc:`fix balance
-
+<fix_balance>` commands.
-.. note::
-
-   Normally, restart files are written by the
-   :doc:`restart <restart>` or :doc:`write_restart <write_restart>` commands
-   so that all atoms in the restart file are inside the simulation box.
-   If this is not the case, the read_restart command will print an error
-   that atoms were "lost" when the file is read.  This error should be
-   reported to the LAMMPS developers so the invalid writing of the
-   restart file can be fixed.  If you still wish to use the restart file,
-   the optional *remap* flag can be appended to the read_restart command.
-   This should avoid the error, by explicitly remapping each atom back
-   into the simulation box, updating image flags for the atom
-   appropriately.
 
 Restart files are saved in binary format to enable exact restarts,
 meaning that the trajectories of a restarted run will precisely match
 those produced by the original run had it continued on.
 
+The binary restart file format was not designed with backward, forward,
+or cross-platform compatibility in mind, so the files are only expected
+to be read correctly by the same LAMMPS executable on the same platform.
+Changes to the architecture, compilation settings, or LAMMPS version can
+render a restart file unreadable or it may read the data incorrectly.
+If you want a more portable format, you can use the data file format as
+created by the :doc:`write_data <write_data>` command.  Binary restart
+files can also be converted into a data file from the command line by
+the LAMMPS executable that wrote them using the :ref:`-restart2data
+<restart2data>` command line flag.
+
 Several things can prevent exact restarts due to round-off effects, in
 which case the trajectories in the 2 runs will slowly diverge.  These
 include running on a different number of processors or changing
@@ -65,7 +60,8 @@ certain settings such as those set by the :doc:`newton <newton>` or
 these cases.
 
 Certain fixes will not restart exactly, though they should provide
-statistically similar results.  These include :doc:`fix shake <fix_shake>` and :doc:`fix langevin <fix_langevin>`.
+statistically similar results.  These include :doc:`fix shake
+<fix_shake>` and :doc:`fix langevin <fix_langevin>`.
 
 Certain pair styles will not restart exactly, though they should
 provide statistically similar results.  This is because the forces
@@ -81,18 +77,19 @@ produced the restart file, it could be a LAMMPS bug, so consider
 :doc:`reporting it <Errors_bugs>` if you think the behavior is a bug.
 
 Because restart files are binary, they may not be portable to other
-machines.  In this case, you can use the :doc:`-restart command-line switch <Run_options>` to convert a restart file to a data file.
+machines.  In this case, you can use the :doc:`-restart command-line
+switch <Run_options>` to convert a restart file to a data file.
 
-Similar to how restart files are written (see the
+Similar to how restart files are written (see the :doc:`write_restart
-:doc:`write_restart <write_restart>` and :doc:`restart <restart>`
+<write_restart>` and :doc:`restart <restart>` commands), the restart
-commands), the restart filename can contain two wild-card characters.
+filename can contain two wild-card characters.  If a "\*" appears in the
-If a "\*" appears in the filename, the directory is searched for all
+filename, the directory is searched for all filenames that match the
-filenames that match the pattern where "\*" is replaced with a timestep
+pattern where "\*" is replaced with a timestep value.  The file with the
-value.  The file with the largest timestep value is read in.  Thus,
+largest timestep value is read in.  Thus, this effectively means, read
-this effectively means, read the latest restart file.  It's useful if
+the latest restart file.  It's useful if you want your script to
-you want your script to continue a run from where it left off.  See
+continue a run from where it left off.  See the :doc:`run <run>` command
-the :doc:`run <run>` command and its "upto" option for how to specify
+and its "upto" option for how to specify the run command so it does not
-the run command so it does not need to be changed either.
+need to be changed either.
 
 If a "%" character appears in the restart filename, LAMMPS expects a
 set of multiple files to exist.  The :doc:`restart <restart>` and
@@ -222,17 +219,17 @@ its calculations in a consistent manner.
 
 .. note::
 
-   There are a handful of commands which can be used before or
+   There are a handful of commands which can be used before or between
-   between runs which may require a system initialization.  Examples
+   runs which may require a system initialization.  Examples include the
-   include the "balance", "displace_atoms", "delete_atoms", "set" (some
+   "balance", "displace_atoms", "delete_atoms", "set" (some options),
-   options), and "velocity" (some options) commands.  This is because
+   and "velocity" (some options) commands.  This is because they can
-   they can migrate atoms to new processors.  Thus they will also discard
+   migrate atoms to new processors.  Thus they will also discard unused
-   unused "state" information from fixes.  You will know the discard has
+   "state" information from fixes.  You will know the discard has
    occurred because a list of discarded fixes will be printed to the
    screen and log file, as explained above.  This means that if you wish
    to retain that info in a restarted run, you must re-specify the
-   relevant fixes and computes (which create fixes) before those commands
+   relevant fixes and computes (which create fixes) before those
-   are used.
+   commands are used.
 
 Some pair styles, like the :doc:`granular pair styles <pair_gran>`, also
 use a fix to store "state" information that persists from timestep to
@@ -245,18 +242,19 @@ LAMMPS allows bond interactions (angle, etc) to be turned off or
 deleted in various ways, which can affect how their info is stored in
 a restart file.
 
-If bonds (angles, etc) have been turned off by the :doc:`fix shake <fix_shake>` or :doc:`delete_bonds <delete_bonds>` command,
+If bonds (angles, etc) have been turned off by the :doc:`fix shake
-their info will be written to a restart file as if they are turned on.
+<fix_shake>` or :doc:`delete_bonds <delete_bonds>` command, their info
-This means they will need to be turned off again in a new run after
+will be written to a restart file as if they are turned on.  This means
-the restart file is read.
+they will need to be turned off again in a new run after the restart
+file is read.
 
 Bonds that are broken (e.g. by a bond-breaking potential) are written
 to the restart file as broken bonds with a type of 0.  Thus these
 bonds will still be broken when the restart file is read.
 
-Bonds that have been broken by the :doc:`fix bond/break <fix_bond_break>` command have disappeared from the
+Bonds that have been broken by the :doc:`fix bond/break
-system.  No information about these bonds is written to the restart
+<fix_bond_break>` command have disappeared from the system.  No
-file.
+information about these bonds is written to the restart file.
 
 ----------
 

doc/utils/sphinx-config/false_positives.txt

@@ -120,6 +120,7 @@ Antonelli
 api
 Apoorva
 Appl
+apptainer
 Apu
 arallel
 arccos
@@ -372,6 +373,7 @@ Caltech
 Camilloni
 Camiloni
 Campana
+Cangi
 Cao
 Capolungo
 Caro
@@ -2368,7 +2370,11 @@ nopreliminary
 Nord
 norder
 Nordlund
+noremap
 normals
+normx
+normy
+normz
 Noskov
 noslip
 noticable
@@ -2676,6 +2682,7 @@ polyhedra
 Polym
 polymorphism
 popen
+Popoola
 Popov
 popstore
 Poresag
@@ -2819,6 +2826,7 @@ radians
 Rafferty
 rahman
 Rahman
+Rajamanickam
 Ralf
 Raman
 ramped
@@ -3671,13 +3679,19 @@ vx
 Vx
 vxcm
 vxmu
+vxx
+vxy
+vxz
 vy
 Vy
 vycm
+vyy
+vyz
 vz
 Vz
 vzcm
 vzi
+vzz
 Waals
 Wadley
 wallstyle

examples/COUPLE/fortran2/LAMMPS-wrapper.cpp

@@ -12,8 +12,8 @@
 ------------------------------------------------------------------------- */
 
 /* ------------------------------------------------------------------------
-    Contributing author:  Karl D. Hammond <karlh@ugcs.caltech.edu>
+    Contributing author:  Karl D. Hammond <hammondkd@missouri.edu>
-                          University of Tennessee, Knoxville (USA), 2012
+                          University of Missouri (USA), 2012
 ------------------------------------------------------------------------- */
 
 /* This is set of "wrapper" functions to assist LAMMPS.F90, which itself
@@ -23,6 +23,7 @@
 
 #include <mpi.h>
 #include "LAMMPS-wrapper.h"
+#define LAMMPS_LIB_MPI 1
 #include <library.h>
 #include <lammps.h>
 #include <atom.h>
@@ -56,181 +57,40 @@ void lammps_error_all (void *ptr, const char *file, int line, const char *str)
 
 int lammps_extract_compute_vectorsize (void *ptr, char *id, int style)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *size;
-   int icompute = lmp->modify->find_compute(id);
+   size = (int *) lammps_extract_compute(ptr, id, style, LMP_SIZE_VECTOR);
-   if ( icompute < 0 ) return 0;
+   if (size) return *size;
-   class Compute *compute = lmp->modify->compute[icompute];
+   return 0;
-
-   if ( style == 0 )
-   {
-      if ( !compute->vector_flag )
-         return 0;
-      else
-         return compute->size_vector;
-   }
-   else if ( style == 1 )
-   {
-      return lammps_get_natoms (ptr);
-   }
-   else if ( style == 2 )
-   {
-      if ( !compute->local_flag )
-         return 0;
-      else
-         return compute->size_local_rows;
-   }
-   else
-      return 0;
 }
 
 void lammps_extract_compute_arraysize (void *ptr, char *id, int style,
       int *nrows, int *ncols)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *tmp;
-   int icompute = lmp->modify->find_compute(id);
+   tmp = (int *) lammps_extract_compute(ptr, id, style, LMP_SIZE_ROWS);
-   if ( icompute < 0 )
+   if (tmp) *nrows = *tmp;
-   {
+   tmp = (int *) lammps_extract_compute(ptr, id, style, LMP_SIZE_COLS);
-      *nrows = 0;
+   if (tmp) *ncols = *tmp;
-      *ncols = 0;
-   }
-   class Compute *compute = lmp->modify->compute[icompute];
-
-   if ( style == 0 )
-   {
-      if ( !compute->array_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = compute->size_array_rows;
-         *ncols = compute->size_array_cols;
-      }
-   }
-   else if ( style == 1 )
-   {
-      if ( !compute->peratom_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = lammps_get_natoms (ptr);
-         *ncols = compute->size_peratom_cols;
-      }
-   }
-   else if ( style == 2 )
-   {
-      if ( !compute->local_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = compute->size_local_rows;
-         *ncols = compute->size_local_cols;
-      }
-   }
-   else
-   {
-      *nrows = 0;
-      *ncols = 0;
-   }
-
    return;
 }
 
 int lammps_extract_fix_vectorsize (void *ptr, char *id, int style)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *size;
-   int ifix = lmp->modify->find_fix(id);
+   size = (int *) lammps_extract_fix(ptr, id, style, LMP_SIZE_VECTOR, 0, 0);
-   if ( ifix < 0 ) return 0;
+   if (size) return *size;
-   class Fix *fix = lmp->modify->fix[ifix];
+   return 0;
-
-   if ( style == 0 )
-   {
-      if ( !fix->vector_flag )
-         return 0;
-      else
-         return fix->size_vector;
-   }
-   else if ( style == 1 )
-   {
-      return lammps_get_natoms (ptr);
-   }
-   else if ( style == 2 )
-   {
-      if ( !fix->local_flag )
-         return 0;
-      else
-         return fix->size_local_rows;
-   }
-   else
-      return 0;
 }
 
 void lammps_extract_fix_arraysize (void *ptr, char *id, int style,
       int *nrows, int *ncols)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *tmp;
-   int ifix = lmp->modify->find_fix(id);
+   tmp = (int *) lammps_extract_fix(ptr, id, style, LMP_SIZE_ROWS, 0, 0);
-   if ( ifix < 0 )
+   if (tmp) *nrows = *tmp;
-   {
+   tmp = (int *) lammps_extract_fix(ptr, id, style, LMP_SIZE_COLS, 0, 0);
-      *nrows = 0;
+   if (tmp) *ncols = *tmp;
-      *ncols = 0;
-   }
-   class Fix *fix = lmp->modify->fix[ifix];
-
-   if ( style == 0 )
-   {
-      if ( !fix->array_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = fix->size_array_rows;
-         *ncols = fix->size_array_cols;
-      }
-   }
-   else if ( style == 1 )
-   {
-      if ( !fix->peratom_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = lammps_get_natoms (ptr);
-         *ncols = fix->size_peratom_cols;
-      }
-   }
-   else if ( style == 2 )
-   {
-      if ( !fix->local_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = fix->size_local_rows;
-         *ncols = fix->size_local_cols;
-      }
-   }
-   else
-   {
-      *nrows = 0;
-      *ncols = 0;
-   }
-
    return;
-
 }
 
 /* vim: set ts=3 sts=3 expandtab: */

examples/COUPLE/fortran2/LAMMPS-wrapper.h

@@ -12,8 +12,8 @@
 ------------------------------------------------------------------------- */
 
 /* ------------------------------------------------------------------------
-    Contributing author:  Karl D. Hammond <karlh@ugcs.caltech.edu>
+    Contributing author:  Karl D. Hammond <hammondkd@missouri.edu>
-                          University of Tennessee, Knoxville (USA), 2012
+                          University of Missouri (USA), 2012
 ------------------------------------------------------------------------- */
 
 /* This is set of "wrapper" functions to assist LAMMPS.F90, which itself

examples/COUPLE/fortran2/LAMMPS.F90

@@ -1292,7 +1292,7 @@ contains !! Wrapper functions local to this module {{{1
       Cname = string2Cstring (name)
       Ccount = size(data) / natoms
       if ( Ccount /= 1 .and. Ccount /= 3 ) &
-         call lammps_error_all (ptr, FLERR, 'lammps_gather_atoms requires&
+         call lammps_error_all (ptr, FLERR, 'lammps_scatter_atoms requires&
             & count to be either 1 or 3')
       Fdata = data
       Cdata = C_loc (Fdata(1))
@@ -1355,7 +1355,7 @@ contains !! Wrapper functions local to this module {{{1
       Cname = string2Cstring (name)
       Ccount = size(data) / ndata
       if ( Ccount /= 1 .and. Ccount /= 3 ) &
-         call lammps_error_all (ptr, FLERR, 'lammps_gather_atoms requires&
+         call lammps_error_all (ptr, FLERR, 'lammps_scatter_atoms requires&
             & count to be either 1 or 3')
       Fdata = data
       Cdata = C_loc (Fdata(1))

examples/COUPLE/fortran2/Makefile

@@ -14,7 +14,7 @@ CXXLIB = -lstdc++    # replace with your C++ runtime libs
 # Flags for Fortran compiler, C++ compiler, and C preprocessor, respectively
 FFLAGS = -O2 -fPIC
 CXXFLAGS = -O2 -fPIC
-CPPFLAGS = -DOMPI_SKIP_MPICXX=1 -DMPICH_SKIP_MPICXX -DLAMMPS_LIB_MPI
+CPPFLAGS = -DOMPI_SKIP_MPICXX=1 -DMPICH_SKIP_MPICXX
 
 all : liblammps_fortran.a liblammps_fortran.so
 

examples/COUPLE/fortran2/README

@@ -11,9 +11,8 @@ This interface was created by Karl Hammond who you can contact with
 questions:
 
 Karl D. Hammond
-University of Tennessee, Knoxville
+University of Missouri
-karlh at ugcs.caltech.edu
+hammondkd at missouri.edu
-karlh at utk.edu
 
 -------------------------------------
 

examples/COUPLE/fortran_dftb/LAMMPS-wrapper.cpp

@@ -12,8 +12,8 @@
 ------------------------------------------------------------------------- */
 
 /* ------------------------------------------------------------------------
-    Contributing author:  Karl D. Hammond <karlh@ugcs.caltech.edu>
+    Contributing author:  Karl D. Hammond <hammondkd@missouri.edu>
-                          University of Tennessee, Knoxville (USA), 2012
+                          University of Missouri (USA), 2012
 ------------------------------------------------------------------------- */
 
 /* This is set of "wrapper" functions to assist LAMMPS.F90, which itself
@@ -23,6 +23,7 @@
 
 #include <mpi.h>
 #include "LAMMPS-wrapper.h"
+#define LAMMPS_LIB_MPI 1
 #include <library.h>
 #include <lammps.h>
 #include <atom.h>
@@ -56,181 +57,40 @@ void lammps_error_all (void *ptr, const char *file, int line, const char *str)
 
 int lammps_extract_compute_vectorsize (void *ptr, char *id, int style)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *size;
-   int icompute = lmp->modify->find_compute(id);
+   size = (int *) lammps_extract_compute(ptr, id, style, LMP_SIZE_VECTOR);
-   if ( icompute < 0 ) return 0;
+   if (size) return *size;
-   class Compute *compute = lmp->modify->compute[icompute];
+   return 0;
-
-   if ( style == 0 )
-   {
-      if ( !compute->vector_flag )
-         return 0;
-      else
-         return compute->size_vector;
-   }
-   else if ( style == 1 )
-   {
-      return lammps_get_natoms (ptr);
-   }
-   else if ( style == 2 )
-   {
-      if ( !compute->local_flag )
-         return 0;
-      else
-         return compute->size_local_rows;
-   }
-   else
-      return 0;
 }
 
 void lammps_extract_compute_arraysize (void *ptr, char *id, int style,
       int *nrows, int *ncols)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *tmp;
-   int icompute = lmp->modify->find_compute(id);
+   tmp = (int *) lammps_extract_compute(ptr, id, style, LMP_SIZE_ROWS);
-   if ( icompute < 0 )
+   if (tmp) *nrows = *tmp;
-   {
+   tmp = (int *) lammps_extract_compute(ptr, id, style, LMP_SIZE_COLS);
-      *nrows = 0;
+   if (tmp) *ncols = *tmp;
-      *ncols = 0;
-   }
-   class Compute *compute = lmp->modify->compute[icompute];
-
-   if ( style == 0 )
-   {
-      if ( !compute->array_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = compute->size_array_rows;
-         *ncols = compute->size_array_cols;
-      }
-   }
-   else if ( style == 1 )
-   {
-      if ( !compute->peratom_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = lammps_get_natoms (ptr);
-         *ncols = compute->size_peratom_cols;
-      }
-   }
-   else if ( style == 2 )
-   {
-      if ( !compute->local_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = compute->size_local_rows;
-         *ncols = compute->size_local_cols;
-      }
-   }
-   else
-   {
-      *nrows = 0;
-      *ncols = 0;
-   }
-
    return;
 }
 
 int lammps_extract_fix_vectorsize (void *ptr, char *id, int style)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *size;
-   int ifix = lmp->modify->find_fix(id);
+   size = (int *) lammps_extract_fix(ptr, id, style, LMP_SIZE_VECTOR, 0, 0);
-   if ( ifix < 0 ) return 0;
+   if (size) return *size;
-   class Fix *fix = lmp->modify->fix[ifix];
+   return 0;
-
-   if ( style == 0 )
-   {
-      if ( !fix->vector_flag )
-         return 0;
-      else
-         return fix->size_vector;
-   }
-   else if ( style == 1 )
-   {
-      return lammps_get_natoms (ptr);
-   }
-   else if ( style == 2 )
-   {
-      if ( !fix->local_flag )
-         return 0;
-      else
-         return fix->size_local_rows;
-   }
-   else
-      return 0;
 }
 
 void lammps_extract_fix_arraysize (void *ptr, char *id, int style,
       int *nrows, int *ncols)
 {
-   class LAMMPS *lmp = (class LAMMPS *) ptr;
+   int *tmp;
-   int ifix = lmp->modify->find_fix(id);
+   tmp = (int *) lammps_extract_fix(ptr, id, style, LMP_SIZE_ROWS, 0, 0);
-   if ( ifix < 0 )
+   if (tmp) *nrows = *tmp;
-   {
+   tmp = (int *) lammps_extract_fix(ptr, id, style, LMP_SIZE_COLS, 0, 0);
-      *nrows = 0;
+   if (tmp) *ncols = *tmp;
-      *ncols = 0;
-   }
-   class Fix *fix = lmp->modify->fix[ifix];
-
-   if ( style == 0 )
-   {
-      if ( !fix->array_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = fix->size_array_rows;
-         *ncols = fix->size_array_cols;
-      }
-   }
-   else if ( style == 1 )
-   {
-      if ( !fix->peratom_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = lammps_get_natoms (ptr);
-         *ncols = fix->size_peratom_cols;
-      }
-   }
-   else if ( style == 2 )
-   {
-      if ( !fix->local_flag )
-      {
-         *nrows = 0;
-         *ncols = 0;
-      }
-      else
-      {
-         *nrows = fix->size_local_rows;
-         *ncols = fix->size_local_cols;
-      }
-   }
-   else
-   {
-      *nrows = 0;
-      *ncols = 0;
-   }
-
    return;
-
 }
 
 /* vim: set ts=3 sts=3 expandtab: */

examples/COUPLE/fortran_dftb/LAMMPS-wrapper.h

@@ -12,8 +12,8 @@
 ------------------------------------------------------------------------- */
 
 /* ------------------------------------------------------------------------
-    Contributing author:  Karl D. Hammond <karlh@ugcs.caltech.edu>
+    Contributing author:  Karl D. Hammond <hammondkd@missouri.edu>
-                          University of Tennessee, Knoxville (USA), 2012
+                          University of Missouri (USA), 2012
 ------------------------------------------------------------------------- */
 
 /* This is set of "wrapper" functions to assist LAMMPS.F90, which itself

examples/COUPLE/fortran_dftb/LAMMPS-wrapper2.cpp

@@ -12,8 +12,7 @@
 ------------------------------------------------------------------------- */
 
 /* ------------------------------------------------------------------------
-    Contributing author:  Karl D. Hammond <karlh@ugcs.caltech.edu>
+    Contributing author:  Nir Goldman, LLNL <ngoldman@llnl.gov>, 2016
-                          University of Tennessee, Knoxville (USA), 2012
 ------------------------------------------------------------------------- */
 
 /* This is set of "wrapper" functions to assist LAMMPS.F90, which itself

examples/COUPLE/fortran_dftb/LAMMPS-wrapper2.h

@@ -12,7 +12,7 @@
 ------------------------------------------------------------------------- */
 
 /* ------------------------------------------------------------------------
-    Contributing author:  Nir Goldman, ngoldman@llnl.gov, Oct. 19th, 2016
+    Contributing author:  Nir Goldman, LLNL <ngoldman@llnl.gov>, 2016
 ------------------------------------------------------------------------- */
 
 /* This is set of "wrapper" functions to assist LAMMPS.F90, which itself

examples/COUPLE/fortran_dftb/LAMMPS.F90

@@ -12,8 +12,8 @@
 !--------------------------------------------------------------------------
 
 !! ------------------------------------------------------------------------
-!   Contributing author:  Karl D. Hammond <karlh@ugcs.caltech.edu>
+!   Contributing author:  Karl D. Hammond <hammondkd@missouri.edu>
-!                         University of Tennessee, Knoxville (USA), 2012
+!                         University of Missouri (USA), 2012
 !--------------------------------------------------------------------------
 
 !! LAMMPS, a Fortran 2003 module containing an interface between Fortran

examples/ELASTIC_T/BORN_MATRIX/Argon/Analytical/compute_born.py

@@ -1,7 +1,6 @@
 #!/usr/bin/env python3
 # -*- coding: utf-8 -*-
 
-import sys
 import numpy as np
 
 def reduce_Born(Cf):

examples/ELASTIC_T/BORN_MATRIX/Argon/Numdiff/compute_born.py

@@ -1,7 +1,6 @@
 #!/usr/bin/env python3
 # -*- coding: utf-8 -*-
 
-import sys
 import numpy as np
 
 def reduce_Born(Cf):

examples/ELASTIC_T/BORN_MATRIX/Silicon/elastic_utils.py

@@ -1,5 +1,5 @@
 import numpy as np
-from lammps import lammps, LAMMPS_INT, LMP_STYLE_GLOBAL, LMP_VAR_EQUAL, LMP_VAR_ATOM
+from lammps import lammps, LMP_VAR_EQUAL
 
 # method for rotating elastic constants
 

examples/PACKAGES/pace/plugin/CMakeLists.txt

@@ -0,0 +1,36 @@
+##########################################
+# CMake build system for plugin examples.
+# The is meant to be used as a template for plugins that are
+# distributed independent from the LAMMPS package.
+##########################################
+
+cmake_minimum_required(VERSION 3.10)
+
+project(paceplugin VERSION 1.0 LANGUAGES CXX)
+
+set(CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR})
+include(CheckIncludeFileCXX)
+include(LAMMPSInterfacePlugin)
+include(ML-PACE)
+
+##########################
+# building the plugins
+
+add_library(paceplugin MODULE paceplugin.cpp ${LAMMPS_SOURCE_DIR}/ML-PACE/pair_pace.cpp)
+target_link_libraries(paceplugin PRIVATE pace)
+target_link_libraries(paceplugin PRIVATE lammps)
+target_include_directories(paceplugin PRIVATE ${LAMMPS_SOURCE_DIR}/ML-PACE)
+set_target_properties(paceplugin PROPERTIES PREFIX "" SUFFIX ".so")
+
+# MacOS seems to need this
+if(CMAKE_SYSTEM_NAME STREQUAL Darwin)
+  set_target_properties(paceplugin PROPERTIES LINK_FLAGS "-Wl,-undefined,dynamic_lookup")
+elseif(CMAKE_SYSTEM_NAME STREQUAL "Windows")
+# tell CMake to export all symbols to a .dll on Windows with special case for MinGW cross-compilers
+  set_target_properties(paceplugin PROPERTIES WINDOWS_EXPORT_ALL_SYMBOLS TRUE)
+  if(CMAKE_CROSSCOMPILING)
+    set_target_properties(paceplugin  PROPERTIES LINK_FLAGS "-Wl,--export-all-symbols")
+  endif()
+else()
+  set_target_properties(paceplugin PROPERTIES LINK_FLAGS "-rdynamic")
+endif()

examples/PACKAGES/pace/plugin/LAMMPSInterfacePlugin.cmake

@@ -0,0 +1 @@
+../../../../cmake/Modules/LAMMPSInterfacePlugin.cmake

examples/PACKAGES/pace/plugin/ML-PACE.cmake

@@ -0,0 +1 @@
+../../../../cmake/Modules/Packages/ML-PACE.cmake

examples/PACKAGES/pace/plugin/README.txt

@@ -0,0 +1,2 @@
+This folder contains a loader and support files to build the ML-PACE package as plugin.
+For more information please see: https://docs.lammps.org/Developer_plugins.html

examples/PACKAGES/pace/plugin/paceplugin.cpp

@@ -0,0 +1,28 @@
+
+#include "lammpsplugin.h"
+#include "version.h"
+
+#include "pair_pace.h"
+
+using namespace LAMMPS_NS;
+
+static Pair *pair_pace_creator(LAMMPS *lmp)
+{
+  return new PairPACE(lmp);
+}
+
+extern "C" void lammpsplugin_init(void *lmp, void *handle, void *regfunc)
+{
+  lammpsplugin_t plugin;
+  lammpsplugin_regfunc register_plugin = (lammpsplugin_regfunc) regfunc;
+
+  // register pace pair style
+  plugin.version = LAMMPS_VERSION;
+  plugin.style = "pair";
+  plugin.name = "pace";
+  plugin.info = "PACE plugin pair style v1.0";
+  plugin.author = "Axel Kohlmeyer (akohlmey@gmail.com)";
+  plugin.creator.v1 = (lammpsplugin_factory1 *) &pair_pace_creator;
+  plugin.handle = handle;
+  (*register_plugin)(&plugin, lmp);
+}

examples/PACKAGES/pace/plugin/paceplugin.nsis

@@ -0,0 +1,157 @@
+#!Nsis Installer Command Script
+#
+# The following external defines are recognized:
+# ${VERSION} = YYYYMMDD
+
+!include "MUI2.nsh"
+!include "FileFunc.nsh"
+
+!define MUI_ICON "lammps.ico"
+!define MUI_UNICON "lammps.ico"
+!define MUI_HEADERIMAGE
+!define MUI_HEADERIMAGE_BITMAP "lammps-text-logo-wide.bmp"
+!define MUI_HEADERIMAGE_RIGHT
+
+Unicode true
+XPStyle on
+
+!include "LogicLib.nsh"
+!addplugindir "envvar/Plugins/x86-unicode"
+!include "x64.nsh"
+
+RequestExecutionLevel user
+
+!macro VerifyUserIsAdmin
+UserInfo::GetAccountType
+pop $0
+${If} $0 != "admin"
+  messageBox mb_iconstop "Administrator rights required!"
+  setErrorLevel 740 ;ERROR_ELEVATION_REQUIRED
+  quit
+${EndIf}
+!macroend
+
+!define PACEPLUGIN "LAMMPS ML-PACE Plugin ${VERSION}"
+OutFile "LAMMPS-ML-PACE-plugin-${VERSION}.exe"
+
+Name "${PACEPLUGIN}"
+InstallDir "$LOCALAPPDATA\${PACEPLUGIN}"
+
+ShowInstDetails show
+ShowUninstDetails show
+SetCompressor lzma
+
+!define MUI_ABORTWARNING
+
+!insertmacro MUI_PAGE_DIRECTORY
+!insertmacro MUI_PAGE_INSTFILES
+
+!insertmacro MUI_UNPAGE_CONFIRM
+!insertmacro MUI_UNPAGE_INSTFILES
+
+!insertmacro MUI_LANGUAGE "English"
+
+function .onInit
+  # Determine if LAMMPS was already installed and check whether it was in 32-bit
+  # or 64-bit. Then look up path to uninstaller and offer to uninstall or quit
+  SetRegView 32
+  ReadRegDWORD $0 HKCU "Software\LAMMPS-ML-PACE" "Bits"
+  SetRegView LastUsed
+  ${If} $0 == "32"
+    SetRegView 32
+  ${ElseIf} $0 == "64"
+    SetRegView 64
+  ${Else}
+    SetRegView 64
+  ${EndIf}
+  ClearErrors
+  ReadRegStr $R0 HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" "UninstallString"
+  SetRegView LastUsed
+  ${If} ${Errors}
+    DetailPrint "LAMMPS ML-PACE plugin not (yet) installed"
+  ${Else}
+    MessageBox MB_YESNO "LAMMPS ML-PACE plugin ($0 bit) is already installed. Uninstall existing version?" /SD IDYES IDNO Quit
+      Pop $R1
+      StrCmp $R1 2 Quit +1
+      Exec $R0
+    Quit:
+      Quit
+  ${EndIf}
+  setShellVarContext all
+functionEnd
+
+Section "${PACEPLUGIN}" SecPaceplugin
+  SectionIn RO
+  # Write LAMMPS installation bitness marker. Always use 32-bit registry view
+  SetRegView 32
+  IntFmt $0 "0x%08X" 64
+  WriteRegDWORD HKCU "Software\LAMMPS-ML-PACE" "Bits" $0
+
+  # Switch to "native" registry view
+  SetRegView 64
+  SetShellVarContext current
+
+  SetOutPath "$INSTDIR"
+  File lammps.ico
+  File paceplugin.so
+
+  # Register Application and its uninstaller
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "DisplayName" "${PACEPLUGIN}"
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "Publisher" "The LAMMPS and PACE Developers"
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "URLInfoAbout" "lammps.org"
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "DisplayIcon" "$INSTDIR\lammps.ico"
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "DisplayVersion" "${VERSION}"
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "InstallLocation" "$INSTDIR"
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "UninstallString" "$\"$INSTDIR\uninstall.exe$\""
+  WriteRegStr HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "QuietUninstallString" "$\"$INSTDIR\uninstall.exe$\" /S"
+
+  ${GetSize} "$INSTDIR" "/S=0K" $0 $1 $2
+  IntFmt $0 "0x%08X" $0
+  WriteRegDWORD HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE" \
+                 "EstimatedSize" "$0"
+
+  # update path variables
+  EnVar::SetHKCU
+  # add to LAMMPS plugin search path
+  EnVar::AddValue "LAMMPS_PLUGIN_PATH" "$INSTDIR"
+
+  WriteUninstaller "$INSTDIR\Uninstall.exe"
+SectionEnd
+
+function un.onInit
+  SetShellVarContext current
+functionEnd
+
+Section "Uninstall"
+  # remove LAMMPS bitness/installation indicator always in 32-bit registry view
+  SetRegView 32
+  DeleteRegKey HKCU "Software\LAMMPS-ML-PACE"
+
+  # unregister extension, and uninstall info
+  SetRegView 64
+  SetShellVarContext current
+  # unregister installation
+  DeleteRegKey HKCU "Software\Microsoft\Windows\CurrentVersion\Uninstall\LAMMPS-ML-PACE"
+
+  # update path variables
+  EnVar::SetHKCU
+  # remove entry from LAMMPS plugin search path
+  EnVar::DeleteValue "LAMMPS_PLUGIN_PATH" "$INSTDIR"
+
+  Delete /REBOOTOK "$INSTDIR\paceplugin.so"
+  Delete /REBOOTOK "$INSTDIR\Uninstall.exe"
+  Delete /REBOOTOK "$INSTDIR\lammps.ico"
+  RMDir /REBOOTOK "$INSTDIR"
+SectionEnd
+
+# Local Variables:
+# mode: sh
+# End:

examples/kim/plugin/CMakeLists.txt

@@ -6,46 +6,11 @@
 
 cmake_minimum_required(VERSION 3.10)
 
-# enforce out-of-source build
-if(${CMAKE_SOURCE_DIR} STREQUAL ${CMAKE_BINARY_DIR})
-  message(FATAL_ERROR "In-source builds are not allowed. You must create and use a build directory. "
-    "Please remove CMakeCache.txt and CMakeFiles first.")
-endif()
-
 project(kimplugin VERSION 1.0 LANGUAGES CXX)
 
-set(LAMMPS_SOURCE_DIR ""  CACHE PATH "Location of LAMMPS sources folder")
-if(NOT LAMMPS_SOURCE_DIR)
-  message(FATAL_ERROR "Must set LAMMPS_SOURCE_DIR")
-endif()
-
-# by default, install into $HOME/.local (not /usr/local),
-# so that no root access (and sudo) is needed
-if(CMAKE_INSTALL_PREFIX_INITIALIZED_TO_DEFAULT)
-  set(CMAKE_INSTALL_PREFIX "$ENV{HOME}/.local" CACHE PATH "Default install path" FORCE)
-endif()
-
-# ugly hacks for MSVC which by default always reports an old C++ standard in the __cplusplus macro
-# and prints lots of pointless warnings about "unsafe" functions
-if(MSVC)
-  add_compile_options(/Zc:__cplusplus)
-  add_compile_options(/wd4244)
-  add_compile_options(/wd4267)
-  add_compile_definitions(_CRT_SECURE_NO_WARNINGS)
-endif()
-
-# C++11 is required
-set(CMAKE_CXX_STANDARD 11)
-set(CMAKE_CXX_STANDARD_REQUIRED ON)
-
-# Need -restrict with Intel compilers
-if(CMAKE_CXX_COMPILER_ID STREQUAL "Intel")
-  set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -restrict")
-endif()
-
 set(CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR})
 include(CheckIncludeFileCXX)
-include(LAMMPSInterfaceCXX)
+include(LAMMPSInterfacePlugin)
 
 ##########################
 # building the plugins
@@ -90,9 +55,9 @@ if(CMAKE_SYSTEM_NAME STREQUAL Darwin)
   set_target_properties(kimplugin PROPERTIES LINK_FLAGS "-Wl,-undefined,dynamic_lookup")
 elseif(CMAKE_SYSTEM_NAME STREQUAL "Windows")
 # tell CMake to export all symbols to a .dll on Windows with special case for MinGW cross-compilers
-  set_target_properties(kimplugin.so PROPERTIES WINDOWS_EXPORT_ALL_SYMBOLS TRUE)
+  set_target_properties(kimplugin PROPERTIES WINDOWS_EXPORT_ALL_SYMBOLS TRUE)
   if(CMAKE_CROSSCOMPILING)
-    set_target_properties(kimplugin  PROPERTIES LINK_FLAGS "-Wl,--export-all-symbols")
+    set_target_properties(kimplugin PROPERTIES LINK_FLAGS "-Wl,--export-all-symbols")
   endif()
 else()
   set_target_properties(kimplugin PROPERTIES LINK_FLAGS "-rdynamic")

examples/kim/plugin/LAMMPSInterfaceCXX.cmake

@@ -1,88 +0,0 @@
-# Cmake script code to define the LAMMPS C++ interface
-# settings required for building LAMMPS plugins
-
-################################################################################
-# helper function
-function(validate_option name values)
-  string(TOLOWER ${${name}} needle_lower)
-  string(TOUPPER ${${name}} needle_upper)
-  list(FIND ${values} ${needle_lower} IDX_LOWER)
-  list(FIND ${values} ${needle_upper} IDX_UPPER)
-  if(${IDX_LOWER} LESS 0 AND ${IDX_UPPER} LESS 0)
-    list_to_bulletpoints(POSSIBLE_VALUE_LIST ${${values}})
-    message(FATAL_ERROR "\n########################################################################\n"
-      "Invalid value '${${name}}' for option ${name}\n"
-      "\n"
-      "Possible values are:\n"
-      "${POSSIBLE_VALUE_LIST}"
-      "########################################################################")
-  endif()
-endfunction(validate_option)
-
-#################################################################################
-# LAMMPS C++ interface. We only need the header related parts.
-add_library(lammps INTERFACE)
-target_include_directories(lammps INTERFACE ${LAMMPS_SOURCE_DIR})
-if((CMAKE_SYSTEM_NAME STREQUAL "Windows") AND CMAKE_CROSSCOMPILING)
-  target_link_libraries(lammps INTERFACE ${CMAKE_BINARY_DIR}/../liblammps.dll.a)
-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()
-  option(BUILD_MPI "Build MPI version" OFF)
-endif()
-
-if(BUILD_MPI)
-  find_package(MPI REQUIRED)
-  option(LAMMPS_LONGLONG_TO_LONG "Workaround if your system or MPI version does not recognize 'long long' data types" OFF)
-  if(LAMMPS_LONGLONG_TO_LONG)
-    target_compile_definitions(lammps INTERFACE -DLAMMPS_LONGLONG_TO_LONG)
-  endif()
-  target_link_libraries(lammps INTERFACE MPI::MPI_CXX)
-else()
-  target_include_directories(lammps INTERFACE "${LAMMPS_SOURCE_DIR}/STUBS")
-endif()
-
-set(LAMMPS_SIZES "smallbig" CACHE STRING "LAMMPS integer sizes (smallsmall: all 32-bit, smallbig: 64-bit #atoms #timesteps, bigbig: also 64-bit imageint, 64-bit atom ids)")
-set(LAMMPS_SIZES_VALUES smallbig bigbig smallsmall)
-set_property(CACHE LAMMPS_SIZES PROPERTY STRINGS ${LAMMPS_SIZES_VALUES})
-validate_option(LAMMPS_SIZES LAMMPS_SIZES_VALUES)
-string(TOUPPER ${LAMMPS_SIZES} LAMMPS_SIZES)
-target_compile_definitions(lammps INTERFACE -DLAMMPS_${LAMMPS_SIZES})
-
-################################################################################
-# detect if we may enable OpenMP support by default
-set(BUILD_OMP_DEFAULT OFF)
-find_package(OpenMP QUIET)
-if(OpenMP_FOUND)
-  check_include_file_cxx(omp.h HAVE_OMP_H_INCLUDE)
-  if(HAVE_OMP_H_INCLUDE)
-    set(BUILD_OMP_DEFAULT ON)
-  endif()
-endif()
-
-option(BUILD_OMP "Build with OpenMP support" ${BUILD_OMP_DEFAULT})
-
-if(BUILD_OMP)
-  find_package(OpenMP REQUIRED)
-  check_include_file_cxx(omp.h HAVE_OMP_H_INCLUDE)
-  if(NOT HAVE_OMP_H_INCLUDE)
-    message(FATAL_ERROR "Cannot find the 'omp.h' header file required for full OpenMP support")
-  endif()
-
-  if (((CMAKE_CXX_COMPILER_ID STREQUAL "GNU") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 9.0)) OR
-      (CMAKE_CXX_COMPILER_ID STREQUAL "PGI") OR
-      ((CMAKE_CXX_COMPILER_ID STREQUAL "Clang") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 10.0)) OR
-      ((CMAKE_CXX_COMPILER_ID STREQUAL "Intel") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 19.0)))
-    # GCC 9.x and later plus Clang 10.x and later implement strict OpenMP 4.0 semantics for consts.
-    # Intel 18.0 was tested to support both, so we switch to OpenMP 4+ from 19.x onward to be safe.
-    target_compile_definitions(lammps INTERFACE -DLAMMPS_OMP_COMPAT=4)
-  else()
-    target_compile_definitions(lammps INTERFACE -DLAMMPS_OMP_COMPAT=3)
-  endif()
-  target_link_libraries(lammps INTERFACE OpenMP::OpenMP_CXX)
-endif()

examples/kim/plugin/LAMMPSInterfacePlugin.cmake

@@ -0,0 +1 @@
+../../../cmake/Modules/LAMMPSInterfacePlugin.cmake

examples/kim/plugin/README.txt

@@ -0,0 +1,2 @@
+This folder contains a loader and support files to build the KIM package as plugin.
+For more information please see: https://docs.lammps.org/Developer_plugins.html

examples/mdi/aimd_driver.py

@@ -26,7 +26,7 @@
 #   -nsteps 10
 #     number of timesteps, default = 10
 
-import sys,math,random
+import sys
 import mdi
 import numpy as np
 from mpi4py import MPI
@@ -42,10 +42,9 @@ def error(txt=None):
 def perform_aimd(world,mm_comm,qm_comm):
 
   me = world.Get_rank()
-  nprocs = world.Get_size()
 
   # receive number of atoms from the MM engine
-  
+
   mdi.MDI_Send_command("<NATOMS",mm_comm)
   natoms = mdi.MDI_Recv(1,mdi.MDI_INT,mm_comm)
   natoms = world.bcast(natoms,root=0)

examples/plugins/CMakeLists.txt

@@ -32,9 +32,14 @@ else()
   # ugly hacks for MSVC which by default always reports an old C++ standard in the __cplusplus macro
   # and prints lots of pointless warnings about "unsafe" functions
   if(MSVC)
-    add_compile_options(/Zc:__cplusplus)
+    if(CMAKE_CXX_COMPILER_ID STREQUAL "MSVC")
-    add_compile_options(/wd4244)
+      add_compile_options(/Zc:__cplusplus)
-    add_compile_options(/wd4267)
+      add_compile_options(/wd4244)
+      add_compile_options(/wd4267)
+      if(LAMMPS_EXCEPTIONS)
+        add_compile_options(/EHsc)
+      endif()
+    endif()
     add_compile_definitions(_CRT_SECURE_NO_WARNINGS)
   endif()
 endif()

examples/snap/in.grid.snap

@@ -0,0 +1,94 @@
+# Demonstrate calculation of SNAP bispectrum descriptors on a grid
+
+# CORRECTNESS: The two atom positions coincide with two of
+# the gridpoints, so c_b[2][1-5] should match c_mygrid[8][4-8].
+# The same is true for compute grid/local c_mygridlocal[8][4-11].
+# Local arrays can not be access directly in the script,
+# but they are printed out to file dump.blocal
+
+variable 	nrep index 1
+variable 	a index 3.316
+variable 	ngrid index 2
+
+units		metal
+atom_modify	map hash
+
+# generate the box and atom positions using a BCC lattice
+
+variable       	nx equal ${nrep}
+variable 	ny equal ${nrep}
+variable 	nz equal ${nrep}
+
+boundary	p p p
+
+lattice		custom $a &
+		a1 1 0 0 &
+		a2 0 1 0  &
+		a3 0 0 1 &
+		basis 0 0 0 &
+		basis 0.5 0.5 0.5 &
+
+region		box block 0 ${nx} 0 ${ny} 0 ${nz}
+create_box	1 box
+create_atoms	1 box
+
+mass 		1 180.88
+
+# define atom compute and grid compute
+
+group 		snapgroup type 1
+variable 	twojmax equal 2
+variable 	rcutfac equal 4.67637
+variable 	rfac0 equal 0.99363
+variable 	rmin0 equal 0
+variable 	wj equal 1
+variable 	radelem equal 0.5
+variable 	bzero equal 0
+variable 	quadratic equal 0
+variable 	switch equal 1
+
+variable 	snap_options string &
+		"${rcutfac} ${rfac0} ${twojmax} ${radelem} &
+		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} &
+		bzeroflag ${bzero} switchflag ${switch}"
+
+# build zero potential to satisfy compute sna/atom
+
+pair_style      zero ${rcutfac}
+pair_coeff      * *
+
+# define atom and grid computes
+
+compute       	b all sna/atom ${snap_options}
+compute 	mygrid all sna/grid grid ${ngrid} ${ngrid} ${ngrid} &
+	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid ${ngrid} ${ngrid} ${ngrid} &
+	 	${snap_options}
+
+# define output
+
+variable	B5atom equal c_b[2][5]
+variable	B5grid equal c_mygrid[8][8]
+
+variable	rmse_global equal "sqrt(   &
+	 (c_mygrid[8][1] - x[2])^2 +      &
+	 (c_mygrid[8][2] - y[2])^2 +      &
+	 (c_mygrid[8][3] - z[2])^2 +      &
+	 (c_mygrid[8][4] - c_b[2][1])^2 + &
+	 (c_mygrid[8][5] - c_b[2][2])^2 + &
+	 (c_mygrid[8][6] - c_b[2][3])^2 + &
+	 (c_mygrid[8][7] - c_b[2][4])^2 + &
+	 (c_mygrid[8][8] - c_b[2][5])^2   &
+	 )"
+
+thermo_style	custom step v_B5atom v_B5grid v_rmse_global
+
+# this is the only way to view the local grid
+
+dump 1 all local 1000 dump.blocal c_mygridlocal[*]
+dump 2 all custom 1000 dump.batom id x y z c_b[*]
+
+# run
+
+run		0
+

examples/snap/in.grid.tri

@@ -0,0 +1,114 @@
+# Demonstrate calculation of SNAP bispectrum
+# descriptors on a grid for triclinic cell
+
+# This triclinic cell has 6 times the volume of the single
+# unit cell used by in.grid
+# and contains 12 atoms.  It is a 3x2x1 supercell
+# with each unit cell containing 2 atoms and the
+# reduced lattice vectors are [1 0 0], [1 1 0], and [1 1 1].
+# The grid is listed in x-fastest order
+
+# CORRECTNESS: The atom positions coincide with certain
+# gridpoints, so c_b[1][1-5] should match c_mygrid[1][4-8]
+# and c_b[7][1-5] should match c_mygrid[13][4-8].
+# Local arrays can not be access directly in the script,
+# but they are printed out to file dump.blocal.tri
+
+# Initialize simulation
+
+variable     	nrep index 1
+variable     	a index 3.316
+variable     	ngrid index 2
+
+variable     	nrepx equal 3*${nrep}
+variable     	nrepy equal 2*${nrep}
+variable     	nrepz equal 1*${nrep}
+
+variable     	ngridx equal 3*${ngrid}
+variable     	ngridy equal 2*${ngrid}
+variable     	ngridz equal 1*${ngrid}
+
+units	     	metal
+atom_modify  	map hash sort 0 0
+
+# generate the box and atom positions using a BCC lattice
+
+variable     	nx equal ${nrepx}
+variable     	ny equal ${nrepy}
+variable     	nz equal ${nrepz}
+
+boundary     	p p p
+
+lattice		custom $a &
+		a1 1 0 0 &
+		a2 1 1 0  &
+		a3 1 1 1 &
+		basis 0 0 0 &
+		basis 0.0 0.0 0.5 &
+		spacing 1 1 1
+
+box 		tilt large
+region		box prism 0 ${nx} 0 ${ny} 0 ${nz} ${ny} ${nz} ${nz}
+create_box	1 box
+create_atoms	1 box
+
+mass 		1 180.88
+
+# define atom compute and grid compute
+
+group 		snapgroup type 1
+variable 	twojmax equal 2
+variable 	rcutfac equal 4.67637
+variable 	rfac0 equal 0.99363
+variable 	rmin0 equal 0
+variable 	wj equal 1
+variable 	radelem equal 0.5
+variable 	bzero equal 0
+variable 	quadratic equal 0
+variable 	switch equal 1
+
+variable 	snap_options string &
+		"${rcutfac} ${rfac0} ${twojmax} ${radelem} &
+		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} &
+		bzeroflag ${bzero} switchflag ${switch}"
+
+# build zero potential to satisfy compute sna/atom
+
+pair_style      zero ${rcutfac}
+pair_coeff      * *
+
+# define atom and grid computes
+
+compute       	b all sna/atom ${snap_options}
+compute 	mygrid all sna/grid grid ${ngridx} ${ngridy} ${ngridz} &
+	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid ${ngridx} ${ngridy} ${ngridz} &
+	 	${snap_options}
+
+# define output
+
+variable	B5atom equal c_b[7][5]
+variable	B5grid equal c_mygrid[13][8]
+
+# do not compare x,y,z because assignment of ids
+# to atoms is not unnique for different processor grids
+
+variable	rmse_global equal "sqrt(    &
+	 (c_mygrid[13][4] - c_b[7][1])^2 + &
+	 (c_mygrid[13][5] - c_b[7][2])^2 + &
+	 (c_mygrid[13][6] - c_b[7][3])^2 + &
+	 (c_mygrid[13][7] - c_b[7][4])^2 + &
+	 (c_mygrid[13][8] - c_b[7][5])^2   &
+	 )"
+
+thermo_style	custom step v_B5atom v_B5grid v_rmse_global
+
+# this is the only way to view the local grid
+
+dump 1 all local 1000 dump.blocal.tri c_mygridlocal[*]
+dump 2 all custom 1000 dump.batom.tri id x y z c_b[*]
+
+# run
+
+run		0
+

examples/snap/log.15Jun22.grid.snap.g++.1

@@ -0,0 +1,159 @@
+LAMMPS (2 Jun 2022)
+  using 1 OpenMP thread(s) per MPI task
+# Demonstrate calculation of SNAP bispectrum descriptors on a grid
+
+# CORRECTNESS: The two atom positions coincide with two of
+# the gridpoints, so c_b[2][1-5] should match c_mygrid[8][4-8].
+# The same is true for compute grid/local c_mygridlocal[8][4-11].
+# Local arrays can not be access directly in the script,
+# but they are printed out to file dump.blocal
+
+variable 	nrep index 1
+variable 	a index 3.316
+variable 	ngrid index 2
+
+units		metal
+atom_modify	map hash
+
+# generate the box and atom positions using a BCC lattice
+
+variable       	nx equal ${nrep}
+variable       	nx equal 1
+variable 	ny equal ${nrep}
+variable 	ny equal 1
+variable 	nz equal ${nrep}
+variable 	nz equal 1
+
+boundary	p p p
+
+lattice		custom $a 		a1 1 0 0 		a2 0 1 0  		a3 0 0 1 		basis 0 0 0 		basis 0.5 0.5 0.5
+lattice		custom 3.316 		a1 1 0 0 		a2 0 1 0  		a3 0 0 1 		basis 0 0 0 		basis 0.5 0.5 0.5
+Lattice spacing in x,y,z = 3.316 3.316 3.316
+region		box block 0 ${nx} 0 ${ny} 0 ${nz}
+region		box block 0 1 0 ${ny} 0 ${nz}
+region		box block 0 1 0 1 0 ${nz}
+region		box block 0 1 0 1 0 1
+create_box	1 box
+Created orthogonal box = (0 0 0) to (3.316 3.316 3.316)
+  1 by 1 by 1 MPI processor grid
+create_atoms	1 box
+Created 2 atoms
+  using lattice units in orthogonal box = (0 0 0) to (3.316 3.316 3.316)
+  create_atoms CPU = 0.000 seconds
+
+mass 		1 180.88
+
+# define atom compute and grid compute
+
+group 		snapgroup type 1
+2 atoms in group snapgroup
+variable 	twojmax equal 2
+variable 	rcutfac equal 4.67637
+variable 	rfac0 equal 0.99363
+variable 	rmin0 equal 0
+variable 	wj equal 1
+variable 	radelem equal 0.5
+variable 	bzero equal 0
+variable 	quadratic equal 0
+variable 	switch equal 1
+
+variable 	snap_options string 		"${rcutfac} ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}"
+4.67637 ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# build zero potential to satisfy compute sna/atom
+
+pair_style      zero ${rcutfac}
+pair_style      zero 4.67637
+pair_coeff      * *
+
+# define atom and grid computes
+
+compute       	b all sna/atom ${snap_options}
+compute       	b all sna/atom 4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygrid all sna/grid grid ${ngrid} ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 2 ${ngrid} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 2 2 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 2 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygridlocal all sna/grid/local grid ${ngrid} ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 2 ${ngrid} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 2 2 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 2 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# define output
+
+variable	B5atom equal c_b[2][5]
+variable	B5grid equal c_mygrid[8][8]
+
+variable	rmse_global equal "sqrt(   	 (c_mygrid[8][1] - x[2])^2 +      	 (c_mygrid[8][2] - y[2])^2 +      	 (c_mygrid[8][3] - z[2])^2 +      	 (c_mygrid[8][4] - c_b[2][1])^2 + 	 (c_mygrid[8][5] - c_b[2][2])^2 + 	 (c_mygrid[8][6] - c_b[2][3])^2 + 	 (c_mygrid[8][7] - c_b[2][4])^2 + 	 (c_mygrid[8][8] - c_b[2][5])^2   	 )"
+
+thermo_style	custom step v_B5atom v_B5grid v_rmse_global
+
+# this is the only way to view the local grid
+
+dump 1 all local 1000 dump.blocal c_mygridlocal[*]
+dump 2 all custom 1000 dump.batom id x y z c_b[*]
+
+# run
+
+run		0
+WARNING: No fixes with time integration, atoms won't move (src/verlet.cpp:60)
+Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
+Neighbor list info ...
+  update every 1 steps, delay 10 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 6.67637
+  ghost atom cutoff = 6.67637
+  binsize = 3.338185, bins = 1 1 1
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair zero, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d
+      bin: standard
+  (2) compute sna/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 7.127 | 7.127 | 7.127 Mbytes
+   Step        v_B5atom       v_B5grid    v_rmse_global 
+         0   1.0427295      1.0427295      9.1551336e-16
+Loop time of 1.43e-06 on 1 procs for 0 steps with 2 atoms
+
+139.9% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 1.43e-06   |            |       |100.00
+
+Nlocal:              2 ave           2 max           2 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:            339 ave         339 max         339 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:             64 ave          64 max          64 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:          128 ave         128 max         128 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 128
+Ave neighs/atom = 64
+Neighbor list builds = 0
+Dangerous builds = 0
+
+Total wall time: 0:00:00

examples/snap/log.15Jun22.grid.snap.g++.4

@@ -0,0 +1,160 @@
+LAMMPS (2 Jun 2022)
+  using 1 OpenMP thread(s) per MPI task
+# Demonstrate calculation of SNAP bispectrum descriptors on a grid
+
+# CORRECTNESS: The two atom positions coincide with two of
+# the gridpoints, so c_b[2][1-5] should match c_mygrid[8][4-8].
+# The same is true for compute grid/local c_mygridlocal[8][4-11].
+# Local arrays can not be access directly in the script,
+# but they are printed out to file dump.blocal
+
+variable 	nrep index 1
+variable 	a index 3.316
+variable 	ngrid index 2
+
+units		metal
+atom_modify	map hash
+
+# generate the box and atom positions using a BCC lattice
+
+variable       	nx equal ${nrep}
+variable       	nx equal 1
+variable 	ny equal ${nrep}
+variable 	ny equal 1
+variable 	nz equal ${nrep}
+variable 	nz equal 1
+
+boundary	p p p
+
+lattice		custom $a 		a1 1 0 0 		a2 0 1 0  		a3 0 0 1 		basis 0 0 0 		basis 0.5 0.5 0.5
+lattice		custom 3.316 		a1 1 0 0 		a2 0 1 0  		a3 0 0 1 		basis 0 0 0 		basis 0.5 0.5 0.5
+Lattice spacing in x,y,z = 3.316 3.316 3.316
+region		box block 0 ${nx} 0 ${ny} 0 ${nz}
+region		box block 0 1 0 ${ny} 0 ${nz}
+region		box block 0 1 0 1 0 ${nz}
+region		box block 0 1 0 1 0 1
+create_box	1 box
+Created orthogonal box = (0 0 0) to (3.316 3.316 3.316)
+  1 by 2 by 2 MPI processor grid
+create_atoms	1 box
+Created 2 atoms
+  using lattice units in orthogonal box = (0 0 0) to (3.316 3.316 3.316)
+  create_atoms CPU = 0.001 seconds
+
+mass 		1 180.88
+
+# define atom compute and grid compute
+
+group 		snapgroup type 1
+2 atoms in group snapgroup
+variable 	twojmax equal 2
+variable 	rcutfac equal 4.67637
+variable 	rfac0 equal 0.99363
+variable 	rmin0 equal 0
+variable 	wj equal 1
+variable 	radelem equal 0.5
+variable 	bzero equal 0
+variable 	quadratic equal 0
+variable 	switch equal 1
+
+variable 	snap_options string 		"${rcutfac} ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}"
+4.67637 ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# build zero potential to satisfy compute sna/atom
+
+pair_style      zero ${rcutfac}
+pair_style      zero 4.67637
+pair_coeff      * *
+
+# define atom and grid computes
+
+compute       	b all sna/atom ${snap_options}
+compute       	b all sna/atom 4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygrid all sna/grid grid ${ngrid} ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 2 ${ngrid} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 2 2 	 	${snap_options}
+compute 	mygrid all sna/grid grid 2 2 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygridlocal all sna/grid/local grid ${ngrid} ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 ${ngrid} ${ngrid} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 2 ${ngrid} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 2 2 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 2 2 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# define output
+
+variable	B5atom equal c_b[2][5]
+variable	B5grid equal c_mygrid[8][8]
+
+variable	rmse_global equal "sqrt(   	 (c_mygrid[8][1] - x[2])^2 +      	 (c_mygrid[8][2] - y[2])^2 +      	 (c_mygrid[8][3] - z[2])^2 +      	 (c_mygrid[8][4] - c_b[2][1])^2 + 	 (c_mygrid[8][5] - c_b[2][2])^2 + 	 (c_mygrid[8][6] - c_b[2][3])^2 + 	 (c_mygrid[8][7] - c_b[2][4])^2 + 	 (c_mygrid[8][8] - c_b[2][5])^2   	 )"
+
+thermo_style	custom step v_B5atom v_B5grid v_rmse_global
+
+# this is the only way to view the local grid
+
+dump 1 all local 1000 dump.blocal c_mygridlocal[*]
+dump 2 all custom 1000 dump.batom id x y z c_b[*]
+
+# run
+
+run		0
+WARNING: No fixes with time integration, atoms won't move (src/verlet.cpp:60)
+Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
+Neighbor list info ...
+  update every 1 steps, delay 10 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 6.67637
+  ghost atom cutoff = 6.67637
+  binsize = 3.338185, bins = 1 1 1
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair zero, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d
+      bin: standard
+  (2) compute sna/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+WARNING: Proc sub-domain size < neighbor skin, could lead to lost atoms (src/domain.cpp:970)
+Per MPI rank memory allocation (min/avg/max) = 6.123 | 6.631 | 7.139 Mbytes
+   Step        v_B5atom       v_B5grid    v_rmse_global 
+         0   1.0427295      1.0427295      1.6316879e-15
+Loop time of 2.57125e-06 on 4 procs for 0 steps with 2 atoms
+
+107.0% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 2.571e-06  |            |       |100.00
+
+Nlocal:            0.5 ave           1 max           0 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+Nghost:          274.5 ave         275 max         274 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+Neighs:             16 ave          40 max           0 min
+Histogram: 2 0 0 0 0 0 1 0 0 1
+FullNghs:           32 ave          64 max           0 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+
+Total # of neighbors = 128
+Ave neighs/atom = 64
+Neighbor list builds = 0
+Dangerous builds = 0
+
+Total wall time: 0:00:00

examples/snap/log.15Jun22.grid.tri.g++.1

@@ -0,0 +1,192 @@
+LAMMPS (2 Jun 2022)
+  using 1 OpenMP thread(s) per MPI task
+# Demonstrate calculation of SNAP bispectrum
+# descriptors on a grid for triclinic cell
+
+# This triclinic cell has 6 times the volume of the single
+# unit cell used by in.grid
+# and contains 12 atoms.  It is a 3x2x1 supercell
+# with each unit cell containing 2 atoms and the
+# reduced lattice vectors are [1 0 0], [1 1 0], and [1 1 1].
+# The grid is listed in x-fastest order
+
+# CORRECTNESS: The atom positions coincide with certain
+# gridpoints, so c_b[1][1-5] should match c_mygrid[1][4-8]
+# and c_b[7][1-5] should match c_mygrid[13][4-8].
+# Local arrays can not be access directly in the script,
+# but they are printed out to file dump.blocal.tri
+
+# Initialize simulation
+
+variable     	nrep index 1
+variable     	a index 3.316
+variable     	ngrid index 2
+
+variable     	nrepx equal 3*${nrep}
+variable     	nrepx equal 3*1
+variable     	nrepy equal 2*${nrep}
+variable     	nrepy equal 2*1
+variable     	nrepz equal 1*${nrep}
+variable     	nrepz equal 1*1
+
+variable     	ngridx equal 3*${ngrid}
+variable     	ngridx equal 3*2
+variable     	ngridy equal 2*${ngrid}
+variable     	ngridy equal 2*2
+variable     	ngridz equal 1*${ngrid}
+variable     	ngridz equal 1*2
+
+units	     	metal
+atom_modify  	map hash sort 0 0
+
+# generate the box and atom positions using a BCC lattice
+
+variable     	nx equal ${nrepx}
+variable     	nx equal 3
+variable     	ny equal ${nrepy}
+variable     	ny equal 2
+variable     	nz equal ${nrepz}
+variable     	nz equal 1
+
+boundary     	p p p
+
+lattice		custom $a 		a1 1 0 0 		a2 1 1 0  		a3 1 1 1 		basis 0 0 0 		basis 0.0 0.0 0.5 		spacing 1 1 1
+lattice		custom 3.316 		a1 1 0 0 		a2 1 1 0  		a3 1 1 1 		basis 0 0 0 		basis 0.0 0.0 0.5 		spacing 1 1 1
+Lattice spacing in x,y,z = 3.316 3.316 3.316
+
+box 		tilt large
+region		box prism 0 ${nx} 0 ${ny} 0 ${nz} ${ny} ${nz} ${nz}
+region		box prism 0 3 0 ${ny} 0 ${nz} ${ny} ${nz} ${nz}
+region		box prism 0 3 0 2 0 ${nz} ${ny} ${nz} ${nz}
+region		box prism 0 3 0 2 0 1 ${ny} ${nz} ${nz}
+region		box prism 0 3 0 2 0 1 2 ${nz} ${nz}
+region		box prism 0 3 0 2 0 1 2 1 ${nz}
+region		box prism 0 3 0 2 0 1 2 1 1
+create_box	1 box
+Created triclinic box = (0 0 0) to (9.948 6.632 3.316) with tilt (6.632 3.316 3.316)
+WARNING: Triclinic box skew is large (src/domain.cpp:224)
+  1 by 1 by 1 MPI processor grid
+create_atoms	1 box
+Created 12 atoms
+  using lattice units in triclinic box = (0 0 0) to (9.948 6.632 3.316) with tilt (6.632 3.316 3.316)
+  create_atoms CPU = 0.000 seconds
+
+mass 		1 180.88
+
+# define atom compute and grid compute
+
+group 		snapgroup type 1
+12 atoms in group snapgroup
+variable 	twojmax equal 2
+variable 	rcutfac equal 4.67637
+variable 	rfac0 equal 0.99363
+variable 	rmin0 equal 0
+variable 	wj equal 1
+variable 	radelem equal 0.5
+variable 	bzero equal 0
+variable 	quadratic equal 0
+variable 	switch equal 1
+
+variable 	snap_options string 		"${rcutfac} ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}"
+4.67637 ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# build zero potential to satisfy compute sna/atom
+
+pair_style      zero ${rcutfac}
+pair_style      zero 4.67637
+pair_coeff      * *
+
+# define atom and grid computes
+
+compute       	b all sna/atom ${snap_options}
+compute       	b all sna/atom 4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygrid all sna/grid grid ${ngridx} ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 4 ${ngridz} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 4 2 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 4 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygridlocal all sna/grid/local grid ${ngridx} ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 4 ${ngridz} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 4 2 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 4 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# define output
+
+variable	B5atom equal c_b[7][5]
+variable	B5grid equal c_mygrid[13][8]
+
+# do not compare x,y,z because assignment of ids
+# to atoms is not unnique for different processor grids
+
+variable	rmse_global equal "sqrt(    	 (c_mygrid[13][4] - c_b[7][1])^2 + 	 (c_mygrid[13][5] - c_b[7][2])^2 + 	 (c_mygrid[13][6] - c_b[7][3])^2 + 	 (c_mygrid[13][7] - c_b[7][4])^2 + 	 (c_mygrid[13][8] - c_b[7][5])^2   	 )"
+
+thermo_style	custom step v_B5atom v_B5grid v_rmse_global
+
+# this is the only way to view the local grid
+
+dump 1 all local 1000 dump.blocal.tri c_mygridlocal[*]
+dump 2 all custom 1000 dump.batom.tri id x y z c_b[*]
+
+# run
+
+run		0
+WARNING: No fixes with time integration, atoms won't move (src/verlet.cpp:60)
+Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
+Neighbor list info ...
+  update every 1 steps, delay 10 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 6.67637
+  ghost atom cutoff = 6.67637
+  binsize = 3.338185, bins = 6 3 1
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair zero, perpetual
+      attributes: half, newton on
+      pair build: half/bin/newton/tri
+      stencil: half/bin/3d/tri
+      bin: standard
+  (2) compute sna/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 7.183 | 7.183 | 7.183 Mbytes
+   Step        v_B5atom       v_B5grid    v_rmse_global 
+         0   1.0427295      1.0427295      7.2262471e-14
+Loop time of 1.414e-06 on 1 procs for 0 steps with 12 atoms
+
+70.7% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 1.414e-06  |            |       |100.00
+
+Nlocal:             12 ave          12 max          12 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:            604 ave         604 max         604 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:            384 ave         384 max         384 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:          768 ave         768 max         768 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 768
+Ave neighs/atom = 64
+Neighbor list builds = 0
+Dangerous builds = 0
+
+Total wall time: 0:00:00

examples/snap/log.15Jun22.grid.tri.g++.4

@@ -0,0 +1,192 @@
+LAMMPS (2 Jun 2022)
+  using 1 OpenMP thread(s) per MPI task
+# Demonstrate calculation of SNAP bispectrum
+# descriptors on a grid for triclinic cell
+
+# This triclinic cell has 6 times the volume of the single
+# unit cell used by in.grid
+# and contains 12 atoms.  It is a 3x2x1 supercell
+# with each unit cell containing 2 atoms and the
+# reduced lattice vectors are [1 0 0], [1 1 0], and [1 1 1].
+# The grid is listed in x-fastest order
+
+# CORRECTNESS: The atom positions coincide with certain
+# gridpoints, so c_b[1][1-5] should match c_mygrid[1][4-8]
+# and c_b[7][1-5] should match c_mygrid[13][4-8].
+# Local arrays can not be access directly in the script,
+# but they are printed out to file dump.blocal.tri
+
+# Initialize simulation
+
+variable     	nrep index 1
+variable     	a index 3.316
+variable     	ngrid index 2
+
+variable     	nrepx equal 3*${nrep}
+variable     	nrepx equal 3*1
+variable     	nrepy equal 2*${nrep}
+variable     	nrepy equal 2*1
+variable     	nrepz equal 1*${nrep}
+variable     	nrepz equal 1*1
+
+variable     	ngridx equal 3*${ngrid}
+variable     	ngridx equal 3*2
+variable     	ngridy equal 2*${ngrid}
+variable     	ngridy equal 2*2
+variable     	ngridz equal 1*${ngrid}
+variable     	ngridz equal 1*2
+
+units	     	metal
+atom_modify  	map hash sort 0 0
+
+# generate the box and atom positions using a BCC lattice
+
+variable     	nx equal ${nrepx}
+variable     	nx equal 3
+variable     	ny equal ${nrepy}
+variable     	ny equal 2
+variable     	nz equal ${nrepz}
+variable     	nz equal 1
+
+boundary     	p p p
+
+lattice		custom $a 		a1 1 0 0 		a2 1 1 0  		a3 1 1 1 		basis 0 0 0 		basis 0.0 0.0 0.5 		spacing 1 1 1
+lattice		custom 3.316 		a1 1 0 0 		a2 1 1 0  		a3 1 1 1 		basis 0 0 0 		basis 0.0 0.0 0.5 		spacing 1 1 1
+Lattice spacing in x,y,z = 3.316 3.316 3.316
+
+box 		tilt large
+region		box prism 0 ${nx} 0 ${ny} 0 ${nz} ${ny} ${nz} ${nz}
+region		box prism 0 3 0 ${ny} 0 ${nz} ${ny} ${nz} ${nz}
+region		box prism 0 3 0 2 0 ${nz} ${ny} ${nz} ${nz}
+region		box prism 0 3 0 2 0 1 ${ny} ${nz} ${nz}
+region		box prism 0 3 0 2 0 1 2 ${nz} ${nz}
+region		box prism 0 3 0 2 0 1 2 1 ${nz}
+region		box prism 0 3 0 2 0 1 2 1 1
+create_box	1 box
+Created triclinic box = (0 0 0) to (9.948 6.632 3.316) with tilt (6.632 3.316 3.316)
+WARNING: Triclinic box skew is large (src/domain.cpp:224)
+  2 by 2 by 1 MPI processor grid
+create_atoms	1 box
+Created 12 atoms
+  using lattice units in triclinic box = (0 0 0) to (9.948 6.632 3.316) with tilt (6.632 3.316 3.316)
+  create_atoms CPU = 0.000 seconds
+
+mass 		1 180.88
+
+# define atom compute and grid compute
+
+group 		snapgroup type 1
+12 atoms in group snapgroup
+variable 	twojmax equal 2
+variable 	rcutfac equal 4.67637
+variable 	rfac0 equal 0.99363
+variable 	rmin0 equal 0
+variable 	wj equal 1
+variable 	radelem equal 0.5
+variable 	bzero equal 0
+variable 	quadratic equal 0
+variable 	switch equal 1
+
+variable 	snap_options string 		"${rcutfac} ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}"
+4.67637 ${rfac0} ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 ${twojmax} ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 ${radelem} 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		${wj} rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 ${rmin0} quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag ${quadratic} 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag ${bzero} switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag ${switch}
+4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# build zero potential to satisfy compute sna/atom
+
+pair_style      zero ${rcutfac}
+pair_style      zero 4.67637
+pair_coeff      * *
+
+# define atom and grid computes
+
+compute       	b all sna/atom ${snap_options}
+compute       	b all sna/atom 4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygrid all sna/grid grid ${ngridx} ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 4 ${ngridz} 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 4 2 	 	${snap_options}
+compute 	mygrid all sna/grid grid 6 4 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+compute 	mygridlocal all sna/grid/local grid ${ngridx} ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 ${ngridy} ${ngridz} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 4 ${ngridz} 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 4 2 	 	${snap_options}
+compute 	mygridlocal all sna/grid/local grid 6 4 2 	 	4.67637 0.99363 2 0.5 		1 rmin0 0 quadraticflag 0 		bzeroflag 0 switchflag 1
+
+# define output
+
+variable	B5atom equal c_b[7][5]
+variable	B5grid equal c_mygrid[13][8]
+
+# do not compare x,y,z because assignment of ids
+# to atoms is not unnique for different processor grids
+
+variable	rmse_global equal "sqrt(    	 (c_mygrid[13][4] - c_b[7][1])^2 + 	 (c_mygrid[13][5] - c_b[7][2])^2 + 	 (c_mygrid[13][6] - c_b[7][3])^2 + 	 (c_mygrid[13][7] - c_b[7][4])^2 + 	 (c_mygrid[13][8] - c_b[7][5])^2   	 )"
+
+thermo_style	custom step v_B5atom v_B5grid v_rmse_global
+
+# this is the only way to view the local grid
+
+dump 1 all local 1000 dump.blocal.tri c_mygridlocal[*]
+dump 2 all custom 1000 dump.batom.tri id x y z c_b[*]
+
+# run
+
+run		0
+WARNING: No fixes with time integration, atoms won't move (src/verlet.cpp:60)
+Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
+Neighbor list info ...
+  update every 1 steps, delay 10 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 6.67637
+  ghost atom cutoff = 6.67637
+  binsize = 3.338185, bins = 6 3 1
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair zero, perpetual
+      attributes: half, newton on
+      pair build: half/bin/newton/tri
+      stencil: half/bin/3d/tri
+      bin: standard
+  (2) compute sna/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 7.15 | 7.15 | 7.15 Mbytes
+   Step        v_B5atom       v_B5grid    v_rmse_global 
+         0   1.0427295      1.0427295      1.9367585e-14
+Loop time of 2.65825e-06 on 4 procs for 0 steps with 12 atoms
+
+84.6% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 2.658e-06  |            |       |100.00
+
+Nlocal:              3 ave           4 max           2 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+Nghost:            459 ave         460 max         458 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+Neighs:             96 ave         128 max          64 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+FullNghs:          192 ave         256 max         128 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+
+Total # of neighbors = 768
+Ave neighs/atom = 64
+Neighbor list builds = 0
+Dangerous builds = 0
+
+Total wall time: 0:00:00

lib/install_helpers.py

@@ -1,4 +1,4 @@
-import hashlib,os,subprocess,sys
+import hashlib,os,subprocess
 
 # try to auto-detect the maximum number of available CPUs
 def get_cpus():

lib/latte/Install.py

@@ -71,7 +71,8 @@ buildflag = args.build
 pathflag = args.path is not None
 version = args.version
 suffixflag = args.machine is not None
-suffix = args.machine
+if suffixflag:
+  suffix = args.machine
 
 if pathflag:
   lattedir = args.path
@@ -132,8 +133,6 @@ os.symlink(os.path.join(lattedir, 'src', 'latte_c_bind.o'), 'filelink.o')
 # copy Makefile.lammps.suffix to Makefile.lammps
 
 if suffixflag or not os.path.exists("Makefile.lammps"):
-  if suffix is None:
-    suffix = 'gfortran'
   print("Creating Makefile.lammps")
   if os.path.exists("Makefile.lammps.%s" % suffix):
     shutil.copyfile("Makefile.lammps.%s" % suffix, 'Makefile.lammps')

lib/mdi/Install.py

@@ -39,7 +39,7 @@ url = "https://github.com/MolSSI-MDI/MDI_Library/archive/v%s.tar.gz" % version
 
 # known checksums for different MDI versions. used to validate the download.
 checksums = { \
-              '1.3.2' : '836f5da400d8cff0f0e4435640f9454f', \
+              '1.4.1' : 'f9505fccd4c79301a619f6452dad4ad9', \
               }
 
 # print error message or help
@@ -177,7 +177,7 @@ try:
 except:
   n_cpus = 1
 
-print("Building lib%s.so ..." % lib)
+print("Building lib%s.a ..." % lib)
 cmd = "make -f Makefile.auto clean; make -f Makefile.auto -j%d" % n_cpus
 txt = subprocess.check_output(cmd,shell=True,stderr=subprocess.STDOUT)
 print(txt.decode('UTF-8'))
@@ -201,10 +201,10 @@ makefile_lammps.write(str(rpath_option) + "\n")
 makefile_lammps.close()
 
 
-shared_files = glob.glob( os.path.join( homepath, "liblink", "lib%s.so*" % lib) )
+shared_files = glob.glob( os.path.join( homepath, "liblink", "lib%s.a" % lib) )
 if len(shared_files) > 0:
   print("Build was successful")
 else:
-  error("Build of lib/%s/lib%s.so was NOT successful" % (lib,lib))
+  error("Build of lib/%s/lib%s.a was NOT successful" % (lib,lib))
 if has_extramake and not os.path.exists("Makefile.lammps"):
   print("lib/%s/Makefile.lammps was NOT created" % lib)

lib/mdi/Makefile.g++

@@ -8,7 +8,7 @@ EXTRAMAKE = Makefile.lammps.empty
 
 lib: 	$(OBJ)
 	mkdir -p build
-	cd build; cmake -Dlibtype=SHARED -Dlanguage=CXX -D CMAKE_C_COMPILER=gcc -D CMAKE_CXX_COMPILER=g++ ../MDI_Library; make
+	cd build; cmake -Dlibtype=STATIC -Dlanguage=C -D CMAKE_POSITION_INDEPENDENT_CODE=yes -D CMAKE_C_COMPILER=gcc -D CMAKE_CXX_COMPILER=g++ ../MDI_Library; make
 	@cp $(EXTRAMAKE) Makefile.lammps
 
 # ------ CLEAN ------

lib/mdi/Makefile.gcc

@@ -8,7 +8,7 @@ EXTRAMAKE = Makefile.lammps.empty
 
 lib: 	$(OBJ)
 	mkdir -p build
-	cd build; cmake -Dlibtype=SHARED -D CMAKE_C_COMPILER=gcc ../MDI_Library; make
+	cd build; cmake -Dlibtype=STATIC -Dlanguage=C -D CMAKE_POSITION_INDEPENDENT_CODE=yes -D CMAKE_C_COMPILER=gcc ../MDI_Library; make
 	@cp $(EXTRAMAKE) Makefile.lammps
 
 # ------ CLEAN ------

lib/mdi/Makefile.icc

@@ -8,7 +8,7 @@ EXTRAMAKE = Makefile.lammps.empty
 
 lib: 	$(OBJ)
 	mkdir -p build
-	cd build; cmake -Dlibtype=SHARED -D CMAKE_C_COMPILER=icc -D CMAKE_CXX_COMPILER=icc ../MDI_Library; make
+	cd build; cmake -Dlibtype=STATIC -Dlanguage=C -D CMAKE_POSITION_INDEPENDENT_CODE=yes -D CMAKE_C_COMPILER=icc -D CMAKE_CXX_COMPILER=icpc ../MDI_Library; make
 	@cp $(EXTRAMAKE) Makefile.lammps
 
 # ------ CLEAN ------

lib/mdi/Makefile.mpi

@@ -8,7 +8,7 @@ EXTRAMAKE = Makefile.lammps.empty
 
 lib: 	$(OBJ)
 	mkdir -p build
-	cd build; cmake -Dlibtype=SHARED -D CMAKE_C_COMPILER=mpicc -D CMAKE_CXX_COMPILER=mpicxx ../MDI_Library; make
+	cd build; cmake -Dlibtype=STATIC  -Dlanguage=C -D CMAKE_POSITION_INDEPENDENT_CODE=yes -D CMAKE_C_COMPILER=mpicc -D CMAKE_CXX_COMPILER=mpicxx ../MDI_Library; make
 	@cp $(EXTRAMAKE) Makefile.lammps
 
 # ------ CLEAN ------

python/install.py

@@ -11,7 +11,7 @@ independently and used to build the wheel without installing it.
 """
 
 from __future__ import print_function
-import sys,os,shutil,time,glob,subprocess
+import sys,os,shutil,glob,subprocess
 from argparse import ArgumentParser
 
 parser = ArgumentParser(prog='install.py',
@@ -23,6 +23,8 @@ parser.add_argument("-l", "--lib", required=True,
                     help="path to the compiled LAMMPS shared library")
 parser.add_argument("-n", "--noinstall", action="store_true", default=False,
                     help="only build a binary wheel. Don't attempt to install it")
+parser.add_argument("-w", "--wheeldir", required=False,
+                    help="path to a directory where the created wheel will be stored")
 
 args = parser.parse_args()
 
@@ -30,7 +32,7 @@ args = parser.parse_args()
 
 if args.package:
   if not os.path.exists(args.package):
-    print( "ERROR: LAMMPS package %s does not exist" % args.package)
+    print("ERROR: LAMMPS package %s does not exist" % args.package)
     parser.print_help()
     sys.exit(1)
   else:
@@ -38,12 +40,20 @@ if args.package:
 
 if args.lib:
   if not os.path.exists(args.lib):
-    print( "ERROR: LAMMPS shared library %s does not exist" % args.lib)
+    print("ERROR: LAMMPS shared library %s does not exist" % args.lib)
     parser.print_help()
     sys.exit(1)
   else:
     args.lib = os.path.abspath(args.lib)
 
+if args.wheeldir:
+  if not os.path.exists(args.wheeldir):
+    print("ERROR: directory %s to store the wheel does not exist" % args.wheeldir)
+    parser.print_help()
+    sys.exit(1)
+  else:
+    args.wheeldir = os.path.abspath(args.wheeldir)
+
 # we need to switch to the folder of the python package
 olddir = os.path.abspath('.')
 os.chdir(os.path.dirname(args.package))
@@ -80,10 +90,18 @@ os.remove(os.path.join('lammps',os.path.basename(args.lib)))
 
 # stop here if we were asked not to install the wheel we created
 if args.noinstall:
-    exit(0)
+  if args.wheeldir:
+    for wheel in glob.glob('lammps-*.whl'):
+      shutil.copy(wheel, args.wheeldir)
+      os.remove(wheel)
+  exit(0)
 
 # install the wheel with pip. first try to install in the default environment.
 # that will be a virtual environment, if active, or the system folder.
+# if in a virtual environment, we must not use the python executable
+# that is running this script (configured by cmake), but use "python"
+# from the regular system path. The user may have changed to the virtual
+# environment *after* running cmake.
 # recent versions of pip will automatically drop to use the user folder
 # in case the system folder is not writable.
 
@@ -93,21 +111,35 @@ if args.noinstall:
 # must be uninstalled manually. We must not ignore this and drop
 # back to install into a (forced) user folder.
 
-print("Installing wheel")
+if "VIRTUAL_ENV" in os.environ:
+  print("Installing wheel into virtual environment")
+  py_exe = 'python'
+else:
+  print("Installing wheel into system site-packages folder")
+  py_exe = sys.executable
+
 for wheel in glob.glob('lammps-*.whl'):
   try:
-    txt = subprocess.check_output([sys.executable, '-m', 'pip', 'install', '--force-reinstall', wheel], stderr=subprocess.STDOUT, shell=False)
+    txt = subprocess.check_output([py_exe, '-m', 'pip', 'install', '--force-reinstall', wheel], stderr=subprocess.STDOUT, shell=False)
     print(txt.decode('UTF-8'))
+    if args.wheeldir:
+      shutil.copy(wheel, args.wheeldir)
+    else:
+      shutil.copy(wheel, olddir)
+    os.remove(wheel)
     continue
   except subprocess.CalledProcessError as err:
     errmsg = err.output.decode('UTF-8')
     if errmsg.find("distutils installed"):
       sys.exit(errmsg + "You need to uninstall the LAMMPS python module manually first.\n")
   try:
-    print('Installing wheel into standard site-packages folder failed. Trying user folder now')
+    print('Installing wheel into system site-packages folder failed. Trying user folder now')
     txt = subprocess.check_output([sys.executable, '-m', 'pip', 'install', '--user', '--force-reinstall', wheel], stderr=subprocess.STDOUT, shell=False)
     print(txt.decode('UTF-8'))
+    if args.wheeldir:
+      shutil.copy(wheel, args.wheeldir)
+    else:
+      shutil.copy(wheel, olddir)
+    os.remove(wheel)
   except:
     sys.exit('Failed to install wheel ' + wheel)
-  shutil.copy(wheel, olddir)
-  os.remove(wheel)

python/lammps/core.py

@@ -188,20 +188,17 @@ class lammps(object):
       [c_void_p,POINTER(c_double),POINTER(c_double),c_double,c_double,c_double]
     self.lib.lammps_reset_box.restype = None
 
-    self.lib.lammps_gather_atoms.argtypes = \
+    self.lib.lammps_gather_atoms.argtypes = [c_void_p,c_char_p,c_int,c_int,c_void_p]
-      [c_void_p,c_char_p,c_int,c_int,c_void_p]
     self.lib.lammps_gather_atoms.restype = None
 
-    self.lib.lammps_gather_atoms_concat.argtypes = \
+    self.lib.lammps_gather_atoms_concat.argtypes = [c_void_p,c_char_p,c_int,c_int,c_void_p]
-      [c_void_p,c_char_p,c_int,c_int,c_void_p]
     self.lib.lammps_gather_atoms_concat.restype = None
 
     self.lib.lammps_gather_atoms_subset.argtypes = \
       [c_void_p,c_char_p,c_int,c_int,c_int,POINTER(c_int),c_void_p]
     self.lib.lammps_gather_atoms_subset.restype = None
 
-    self.lib.lammps_scatter_atoms.argtypes = \
+    self.lib.lammps_scatter_atoms.argtypes = [c_void_p,c_char_p,c_int,c_int,c_void_p]
-      [c_void_p,c_char_p,c_int,c_int,c_void_p]
     self.lib.lammps_scatter_atoms.restype = None
 
     self.lib.lammps_scatter_atoms_subset.argtypes = \
@@ -211,20 +208,17 @@ class lammps(object):
     self.lib.lammps_gather_bonds.argtypes = [c_void_p,c_void_p]
     self.lib.lammps_gather_bonds.restype = None
 
-    self.lib.lammps_gather.argtypes = \
+    self.lib.lammps_gather.argtypes = [c_void_p,c_char_p,c_int,c_int,c_void_p]
-      [c_void_p,c_char_p,c_int,c_int,c_void_p]
     self.lib.lammps_gather.restype = None
 
-    self.lib.lammps_gather_concat.argtypes = \
+    self.lib.lammps_gather_concat.argtypes = [c_void_p,c_char_p,c_int,c_int,c_void_p]
-      [c_void_p,c_char_p,c_int,c_int,c_void_p]
     self.lib.lammps_gather_concat.restype = None
 
     self.lib.lammps_gather_subset.argtypes = \
       [c_void_p,c_char_p,c_int,c_int,c_int,POINTER(c_int),c_void_p]
     self.lib.lammps_gather_subset.restype = None
 
-    self.lib.lammps_scatter.argtypes = \
+    self.lib.lammps_scatter.argtypes = [c_void_p,c_char_p,c_int,c_int,c_void_p]
-      [c_void_p,c_char_p,c_int,c_int,c_void_p]
     self.lib.lammps_scatter.restype = None
 
     self.lib.lammps_scatter_subset.argtypes = \
@@ -244,7 +238,8 @@ class lammps(object):
     self.lib.lammps_neighlist_num_elements.argtypes = [c_void_p, c_int]
     self.lib.lammps_neighlist_num_elements.restype  = c_int
 
-    self.lib.lammps_neighlist_element_neighbors.argtypes = [c_void_p, c_int, c_int, POINTER(c_int), POINTER(c_int), POINTER(POINTER(c_int))]
+    self.lib.lammps_neighlist_element_neighbors.argtypes = \
+      [c_void_p, c_int, c_int, POINTER(c_int), POINTER(c_int), POINTER(POINTER(c_int))]
     self.lib.lammps_neighlist_element_neighbors.restype  = None
 
     self.lib.lammps_is_running.argtypes = [c_void_p]
@@ -368,11 +363,9 @@ class lammps(object):
             if type(cmdargs[i]) is str:
               cmdargs[i] = cmdargs[i].encode()
           cargs = (c_char_p*narg)(*cmdargs)
-          self.lib.lammps_open.argtypes = [c_int, c_char_p*narg, \
+          self.lib.lammps_open.argtypes = [c_int, c_char_p*narg, MPI_Comm, c_void_p]
-                                           MPI_Comm, c_void_p]
         else:
-          self.lib.lammps_open.argtypes = [c_int, c_char_p, \
+          self.lib.lammps_open.argtypes = [c_int, c_char_p, MPI_Comm, c_void_p]
-                                           MPI_Comm, c_void_p]
 
         self.opened = 1
         comm_ptr = self.MPI._addressof(comm)
@@ -390,8 +383,7 @@ class lammps(object):
             if type(cmdargs[i]) is str:
               cmdargs[i] = cmdargs[i].encode()
           cargs = (c_char_p*narg)(*cmdargs)
-          self.lib.lammps_open_no_mpi.argtypes = [c_int, c_char_p*narg, \
+          self.lib.lammps_open_no_mpi.argtypes = [c_int, c_char_p*narg, c_void_p]
-                                                  c_void_p]
           self.lmp = c_void_p(self.lib.lammps_open_no_mpi(narg,cargs,None))
         else:
           self.lib.lammps_open_no_mpi.argtypes = [c_int, c_char_p, c_void_p]
@@ -963,17 +955,14 @@ class lammps(object):
       return ptr
 
     elif ctype == LMP_SIZE_COLS:
-      if cstyle == LMP_STYLE_GLOBAL  \
+      if cstyle == LMP_STYLE_GLOBAL or cstyle == LMP_STYLE_ATOM or cstyle == LMP_STYLE_LOCAL:
-         or cstyle == LMP_STYLE_ATOM \
-         or cstyle == LMP_STYLE_LOCAL:
         self.lib.lammps_extract_compute.restype = POINTER(c_int)
         with ExceptionCheck(self):
           ptr = self.lib.lammps_extract_compute(self.lmp,cid,cstyle,ctype)
         return ptr[0]
 
     elif ctype == LMP_SIZE_VECTOR or ctype == LMP_SIZE_ROWS:
-      if cstyle == LMP_STYLE_GLOBAL  \
+      if cstyle == LMP_STYLE_GLOBAL or cstyle == LMP_STYLE_LOCAL:
-         or cstyle == LMP_STYLE_LOCAL:
         self.lib.lammps_extract_compute.restype = POINTER(c_int)
         with ExceptionCheck(self):
           ptr = self.lib.lammps_extract_compute(self.lmp,cid,cstyle,ctype)

python/lammps/numpy_wrapper.py

@@ -165,7 +165,7 @@ class numpy_wrapper:
     """
     value = self.lmp.extract_compute(cid, cstyle, ctype)
 
-    if cstyle in (LMP_STYLE_GLOBAL, LMP_STYLE_LOCAL):
+    if cstyle == LMP_STYLE_GLOBAL:
       if ctype == LMP_TYPE_VECTOR:
         nrows = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_VECTOR)
         return self.darray(value, nrows)
@@ -173,6 +173,13 @@ class numpy_wrapper:
         nrows = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_ROWS)
         ncols = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_COLS)
         return self.darray(value, nrows, ncols)
+    elif cstyle == LMP_STYLE_LOCAL:
+      nrows = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_ROWS)
+      ncols = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_COLS)
+      if ncols == 0:
+        return self.darray(value, nrows)
+      else:
+        return self.darray(value, nrows, ncols)
     elif cstyle == LMP_STYLE_ATOM:
       if ctype == LMP_TYPE_VECTOR:
         nlocal = self.lmp.extract_global("nlocal")

python/makewheel.py

@@ -1,14 +1,26 @@
 #!/usr/bin/env python
 
-import sys,os,shutil
+import sys,os,site
 
-# find python script to activate the virtual environment and source it
+base = os.path.abspath('buildwheel')
 if sys.platform == 'win32':
-  virtenv=os.path.join('buildwheel','Scripts','activate_this.py')
+  bin_dir=os.path.join(base,'Scripts')
 else:
-  virtenv=os.path.join('buildwheel','bin','activate_this.py')
+  bin_dir=os.path.join(base,'bin')
 
-exec(open(virtenv).read(), {'__file__': virtenv})
+# prepend bin to PATH, set venv path
+os.environ["PATH"] = os.pathsep.join([bin_dir] + os.environ.get("PATH", "").split(os.pathsep))
+os.environ["VIRTUAL_ENV"] = base
+
+# add the virtual environments libraries to the host python import mechanism
+prev_length = len(sys.path)
+for lib in "__LIB_FOLDERS__".split(os.pathsep):
+    path = os.path.realpath(os.path.join(bin_dir, lib))
+    site.addsitedir(path)
+sys.path[:] = sys.path[prev_length:] + sys.path[0:prev_length]
+
+sys.real_prefix = sys.prefix
+sys.prefix = base
 
 # update pip and install all requirements to build the wheel
 os.system('python -m pip install --upgrade pip')

python/setup.py

@@ -3,7 +3,7 @@
 from setuptools import setup
 from setuptools.dist import Distribution
 from sys import version_info
-import os,time,shutil
+import os,time
 LAMMPS_PYTHON_DIR = os.path.dirname(os.path.realpath(__file__))
 LAMMPS_DIR = os.path.dirname(LAMMPS_PYTHON_DIR)
 LAMMPS_SOURCE_DIR = os.path.join(LAMMPS_DIR, 'src')
@@ -24,7 +24,7 @@ def get_lammps_version():
 
 class BinaryDistribution(Distribution):
     """Wrapper to enforce creating a binary package"""
-    def has_ext_modules(foo):
+    def has_ext_modules(self):
         return True
 
 if version_info.major >= 3:

src/.gitignore

@@ -173,12 +173,20 @@
 /pair_tdpd.cpp
 /pair_tdpd.h
 
+/compute_grid.cpp
+/compute_grid.h
+/compute_grid_local.cpp
+/compute_grid_local.h
 /compute_sna_atom.cpp
 /compute_sna_atom.h
 /compute_snad_atom.cpp
 /compute_snad_atom.h
 /compute_snav_atom.cpp
 /compute_snav_atom.h
+/compute_sna_grid.cpp
+/compute_sna_grid.h
+/compute_sna_grid_local.cpp
+/compute_sna_grid_local.h
 /compute_snap.cpp
 /compute_snap.h
 /openmp_snap.h
@@ -997,8 +1005,8 @@
 /neb.h
 /netcdf_units.cpp
 /netcdf_units.h
-/pair_3b_table.cpp
+/pair_threebody_table.cpp
-/pair_3b_table.h
+/pair_threebody_table.h
 /pair_adp.cpp
 /pair_adp.h
 /pair_agni.cpp
@@ -1293,8 +1301,8 @@
 /pair_sph_taitwater_morris.h
 /pair_sw.cpp
 /pair_sw.h
-/pair_sw_3b_table.cpp
+/pair_sw_angle_table.cpp
-/pair_sw_3b_table.h
+/pair_sw_angle_table.h
 /pair_sw_mod.cpp
 /pair_sw_mod.h
 /pair_tersoff.cpp

src/ATC/fix_atc.cpp

@@ -284,7 +284,7 @@ FixATC::FixATC(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg),
   int me = ATC::LammpsInterface::instance()->comm_rank();
 
   string groupName(arg[1]);
-  int igroup = group->find(groupName.c_str());
+  int igroup = group->find(groupName);
   int atomCount = group->count(igroup);
 
   try {

src/DIELECTRIC/atom_vec_dielectric.cpp

@@ -78,9 +78,9 @@ AtomVecDielectric::AtomVecDielectric(LAMMPS *_lmp) : AtomVec(_lmp)
     "mu", "area", "ed", "em", "epsilon", "curvature", "q_unscaled"};
   fields_create = {"q", "molecule", "num_bond", "num_angle", "num_dihedral", "num_improper",
     "nspecial", "mu", "area", "ed", "em", "epsilon", "curvature", "q_unscaled"};
-  fields_data_atom = { "id", "molecule", "type", "q", "x", "mu3", "area", "ed", "em", "epsilon",
+  fields_data_atom = {"id", "molecule", "type", "q", "x", "mu3", "area", "ed", "em", "epsilon",
     "curvature"};
-  fields_data_vel = {"id v"};
+  fields_data_vel = {"id", "v"};
   // clang-format on
 
   setup_fields();

src/DPD-MESO/fix_mvv_dpd.cpp

@@ -17,8 +17,8 @@
    modified velocity-Verlet (MVV) algorithm.
    Setting verlet = 0.5 recovers the standard velocity-Verlet algorithm.
 
-   Contributing author: Zhen Li (Brown University)
+   Contributing author: Zhen Li (Clemson University)
-   Email: zhen_li@brown.edu
+   Email: zli7@clemson.edu
 ------------------------------------------------------------------------- */
 
 #include "fix_mvv_dpd.h"

src/DPD-MESO/fix_mvv_edpd.cpp

@@ -17,8 +17,8 @@
    v and edpd_T) using the modified velocity-Verlet (MVV) algorithm.
    Setting verlet = 0.5 recovers the standard velocity-Verlet algorithm.
 
-   Contributing author: Zhen Li (Brown University)
+   Contributing author: Zhen Li (Clemson University)
-   Email: zhen_li@brown.edu
+   Email: zli7@clemson.edu
 
    Please cite the related publication:
    Z. Li, Y.-H. Tang, H. Lei, B. Caswell and G.E. Karniadakis. "Energy-

src/DPD-MESO/fix_mvv_tdpd.cpp

@@ -17,8 +17,8 @@
    v and cc) using the modified velocity-Verlet (MVV) algorithm.
    Setting verlet = 0.5 recovers the standard velocity-Verlet algorithm.
 
-   Contributing author: Zhen Li (Brown University)
+   Contributing author: Zhen Li (Clemson University)
-   Email: zhen_li@brown.edu
+   Email: zli7@clemson.edu
 
    Please cite the related publication:
    Z. Li, A. Yazdani, A. Tartakovsky and G.E. Karniadakis. "Transport

src/DPD-MESO/pair_edpd.cpp

@@ -13,8 +13,8 @@
 ------------------------------------------------------------------------- */
 
 /* ----------------------------------------------------------------------
-   Contributing author: Zhen Li (Brown University)
+   Contributing author: Zhen Li (Clemson University)
-   Email: zhen_li@brown.edu
+   Email: zli7@clemson.edu
 ------------------------------------------------------------------------- */
 
 #include "pair_edpd.h"

src/DPD-MESO/pair_mdpd.cpp

@@ -13,8 +13,8 @@
 ------------------------------------------------------------------------- */
 
 /* ----------------------------------------------------------------------
-   Contributing author: Zhen Li (Brown University)
+   Contributing author: Zhen Li (Clemson University)
-   Email: zhen_li@brown.edu
+   Email: zli7@clemson.edu
 ------------------------------------------------------------------------- */
 
 #include "pair_mdpd.h"

src/DPD-MESO/pair_mdpd_rhosum.cpp

@@ -17,7 +17,7 @@
    before the force calculation.
    The code uses 3D Lucy kernel, it can be modified for other kernels.
 
-   Contributing author: Zhen Li (Brown University)
+   Contributing author: Zhen Li (Clemson University)
 ------------------------------------------------------------------------- */
 
 #include "pair_mdpd_rhosum.h"

src/DPD-MESO/pair_tdpd.cpp

@@ -13,8 +13,8 @@
 ------------------------------------------------------------------------- */
 
 /* ----------------------------------------------------------------------
-   Contributing author: Zhen Li (Brown University)
+   Contributing author: Zhen Li (Clemson University)
-   Email: zhen_li@brown.edu
+   Email: zli7@clemson.edu
 ------------------------------------------------------------------------- */
 
 #include "pair_tdpd.h"

src/DPD-SMOOTH/pair_sdpd_taitwater_isothermal.cpp

@@ -44,7 +44,7 @@ static const double sqrt_2_inv = std::sqrt(0.5);
 /* ---------------------------------------------------------------------- */
 
 PairSDPDTaitwaterIsothermal::PairSDPDTaitwaterIsothermal (LAMMPS *lmp)
-: Pair (lmp) {
+: Pair (lmp), random(nullptr) {
   restartinfo = 0;
   single_enable =0;
 }
@@ -61,6 +61,7 @@ PairSDPDTaitwaterIsothermal::~PairSDPDTaitwaterIsothermal () {
     memory->destroy (soundspeed);
     memory->destroy (B);
   }
+  delete random;
 }
 
 /* ---------------------------------------------------------------------- */

src/EXTRA-COMPUTE/compute_ave_sphere_atom.cpp

@@ -11,10 +11,15 @@
    See the README file in the top-level LAMMPS directory.
 ------------------------------------------------------------------------- */
 
+/* ----------------------------------------------------------------------
+   Contributing author: Stan Moore (SNL)
+------------------------------------------------------------------------- */
+
 #include "compute_ave_sphere_atom.h"
 
 #include "atom.h"
 #include "comm.h"
+#include "domain.h"
 #include "error.h"
 #include "force.h"
 #include "math_const.h"
@@ -98,7 +103,10 @@ void ComputeAveSphereAtom::init()
   }
 
   cutsq = cutoff * cutoff;
-  sphere_vol = 4.0 / 3.0 * MY_PI * cutsq * cutoff;
+  if (domain->dimension == 3)
+    volume = 4.0 / 3.0 * MY_PI * cutsq * cutoff;
+  else
+    volume = MY_PI * cutsq;
 
   // need an occasional full neighbor list
 
@@ -121,7 +129,7 @@ void ComputeAveSphereAtom::compute_peratom()
   double xtmp, ytmp, ztmp, delx, dely, delz, rsq;
   int *ilist, *jlist, *numneigh, **firstneigh;
   int count;
-  double vsum[3], vavg[3], vnet[3];
+  double p[3], vcom[3], vnet[3];
 
   invoked_peratom = update->ntimestep;
 
@@ -152,12 +160,26 @@ void ComputeAveSphereAtom::compute_peratom()
 
   double **x = atom->x;
   double **v = atom->v;
+  double *mass = atom->mass;
+  double *rmass = atom->rmass;
+  int *type = atom->type;
   int *mask = atom->mask;
+  double massone_i, massone_j, totalmass;
+
+  double adof = domain->dimension;
+  double mvv2e = force->mvv2e;
+  double mv2d = force->mv2d;
+  double boltz = force->boltz;
 
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
 
     if (mask[i] & groupbit) {
+      if (rmass)
+        massone_i = rmass[i];
+      else
+        massone_i = mass[type[i]];
+
       xtmp = x[i][0];
       ytmp = x[i][1];
       ztmp = x[i][2];
@@ -167,13 +189,18 @@ void ComputeAveSphereAtom::compute_peratom()
       // i atom contribution
 
       count = 1;
-      vsum[0] = v[i][0];
+      totalmass = massone_i;
-      vsum[1] = v[i][1];
+      p[0] = v[i][0] * massone_i;
-      vsum[2] = v[i][2];
+      p[1] = v[i][1] * massone_i;
+      p[2] = v[i][2] * massone_i;
 
       for (jj = 0; jj < jnum; jj++) {
         j = jlist[jj];
         j &= NEIGHMASK;
+        if (rmass)
+          massone_j = rmass[j];
+        else
+          massone_j = mass[type[j]];
 
         delx = xtmp - x[j][0];
         dely = ytmp - x[j][1];
@@ -181,42 +208,45 @@ void ComputeAveSphereAtom::compute_peratom()
         rsq = delx * delx + dely * dely + delz * delz;
         if (rsq < cutsq) {
           count++;
-          vsum[0] += v[j][0];
+          totalmass += massone_j;
-          vsum[1] += v[j][1];
+          p[0] += v[j][0] * massone_j;
-          vsum[2] += v[j][2];
+          p[1] += v[j][1] * massone_j;
+          p[2] += v[j][2] * massone_j;
         }
       }
 
-      vavg[0] = vsum[0] / count;
+      vcom[0] = p[0] / totalmass;
-      vavg[1] = vsum[1] / count;
+      vcom[1] = p[1] / totalmass;
-      vavg[2] = vsum[2] / count;
+      vcom[2] = p[2] / totalmass;
 
       // i atom contribution
 
-      count = 1;
+      vnet[0] = v[i][0] - vcom[0];
-      vnet[0] = v[i][0] - vavg[0];
+      vnet[1] = v[i][1] - vcom[1];
-      vnet[1] = v[i][1] - vavg[1];
+      vnet[2] = v[i][2] - vcom[2];
-      vnet[2] = v[i][2] - vavg[2];
+      double ke_sum = massone_i * (vnet[0] * vnet[0] + vnet[1] * vnet[1] + vnet[2] * vnet[2]);
-      double ke_sum = vnet[0] * vnet[0] + vnet[1] * vnet[1] + vnet[2] * vnet[2];
 
       for (jj = 0; jj < jnum; jj++) {
         j = jlist[jj];
         j &= NEIGHMASK;
+        if (rmass)
+          massone_j = rmass[j];
+        else
+          massone_j = mass[type[j]];
 
         delx = xtmp - x[j][0];
         dely = ytmp - x[j][1];
         delz = ztmp - x[j][2];
         rsq = delx * delx + dely * dely + delz * delz;
         if (rsq < cutsq) {
-          count++;
+          vnet[0] = v[j][0] - vcom[0];
-          vnet[0] = v[j][0] - vavg[0];
+          vnet[1] = v[j][1] - vcom[1];
-          vnet[1] = v[j][1] - vavg[1];
+          vnet[2] = v[j][2] - vcom[2];
-          vnet[2] = v[j][2] - vavg[2];
+          ke_sum += massone_j * (vnet[0] * vnet[0] + vnet[1] * vnet[1] + vnet[2] * vnet[2]);
-          ke_sum += vnet[0] * vnet[0] + vnet[1] * vnet[1] + vnet[2] * vnet[2];
         }
       }
-      double density = count / sphere_vol;
+      double density = mv2d * totalmass / volume;
-      double temp = ke_sum / 3.0 / count;
+      double temp = mvv2e * ke_sum / (adof * count * boltz);
       result[i][0] = density;
       result[i][1] = temp;
     }

src/EXTRA-COMPUTE/compute_ave_sphere_atom.h

@@ -37,7 +37,7 @@ class ComputeAveSphereAtom : public Compute {
 
  protected:
   int nmax;
-  double cutoff, cutsq, sphere_vol;
+  double cutoff, cutsq, volume;
   class NeighList *list;
 
   double **result;

src/EXTRA-COMPUTE/compute_stress_cartesian.cpp

@@ -266,7 +266,7 @@ void ComputeStressCartesian::compute_array()
   Pair *pair = force->pair;
   double **cutsq = force->pair->cutsq;
 
-  double xi1, xi2, xj1, xj2;
+  double xi1, xi2;
 
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
@@ -305,9 +305,6 @@ void ComputeStressCartesian::compute_array()
           }
         }
       }
-      xj1 = x[j][dir1];
-      xj2 = x[j][dir2];
-
       delx = x[j][0] - xtmp;
       dely = x[j][1] - ytmp;
       delz = x[j][2] - ztmp;
@@ -318,7 +315,7 @@ void ComputeStressCartesian::compute_array()
       // Check if inside cut-off
       if (rsq >= cutsq[itype][jtype]) continue;
       pair->single(i, j, itype, jtype, rsq, factor_coul, factor_lj, fpair);
-      compute_pressure(fpair, xi1, xi2, xj1, xj2, delx, dely, delz);
+      compute_pressure(fpair, xi1, xi2, delx, dely, delz);
     }
   }
 
@@ -356,8 +353,8 @@ void ComputeStressCartesian::compute_array()
   }
 }
 
-void ComputeStressCartesian::compute_pressure(double fpair, double xi, double yi, double xj,
+void ComputeStressCartesian::compute_pressure(double fpair, double xi, double yi, double delx,
-                                              double yj, double delx, double dely, double delz)
+                                              double dely, double delz)
 {
   int bin1, bin2, next_bin1, next_bin2;
   double la = 0.0, lb = 0.0, l_sum = 0.0;

src/EXTRA-COMPUTE/compute_stress_cartesian.h

@@ -41,7 +41,7 @@ class ComputeStressCartesian : public Compute {
   double *dens, *pkxx, *pkyy, *pkzz, *pcxx, *pcyy, *pczz;
   double *tdens, *tpkxx, *tpkyy, *tpkzz, *tpcxx, *tpcyy, *tpczz;
   class NeighList *list;
-  void compute_pressure(double, double, double, double, double, double, double, double);
+  void compute_pressure(double, double, double, double, double, double);
 };
 
 }    // namespace LAMMPS_NS

src/EXTRA-DUMP/dump_xtc.cpp

@@ -433,10 +433,10 @@ int xdropen(XDR *xdrs, const char *filename, const char *type)
     return 0;
   }
   if (*type == 'w' || *type == 'W') {
-    type = (char *) "w+";
+    type = (char *) "wb+";
     lmode = XDR_ENCODE;
   } else {
-    type = (char *) "r";
+    type = (char *) "rb";
     lmode = XDR_DECODE;
   }
   xdrfiles[xdrid] = fopen(filename, type);

src/EXTRA-DUMP/dump_yaml.cpp

@@ -124,6 +124,12 @@ void DumpYAML::write_data(int n, double *mybuf)
     }
     fputs("]\n", fp);
   }
+}
+
+/* ---------------------------------------------------------------------- */
+
+void DumpYAML::write_footer()
+{
   fputs("...\n", fp);
 }
 

src/EXTRA-DUMP/dump_yaml.h

@@ -35,6 +35,7 @@ class DumpYAML : public DumpCustom {
   void write() override;
   void write_header(bigint) override;
   void write_data(int, double *) override;
+  void write_footer() override;
 
   int modify_param(int, char **) override;
 };

src/INTERLAYER/pair_ilp_tmd.cpp

@@ -483,7 +483,7 @@ void PairILPTMD::calc_normal()
       }
     }
     //############################ For the edge atoms of TMD ################################
-    else if (cont > 1 && cont < Nnei) {
+    else if (cont < Nnei) {
       if (strcmp(elements[itype], "Mo") == 0 || strcmp(elements[itype], "W") == 0 ||
           strcmp(elements[itype], "S") == 0 || strcmp(elements[itype], "Se") == 0) {
         // derivatives of Ni[l] respect to the cont neighbors

src/KOKKOS/compute_ave_sphere_atom_kokkos.cpp

@@ -11,11 +11,16 @@
    See the README file in the top-level LAMMPS directory.
 ------------------------------------------------------------------------- */
 
+/* ----------------------------------------------------------------------
+   Contributing author: Stan Moore (SNL)
+------------------------------------------------------------------------- */
+
 #include "compute_ave_sphere_atom_kokkos.h"
 
 #include "atom_kokkos.h"
 #include "atom_masks.h"
 #include "comm.h"
+#include "domain.h"
 #include "error.h"
 #include "force.h"
 #include "memory_kokkos.h"
@@ -105,11 +110,19 @@ void ComputeAveSphereAtomKokkos<DeviceType>::compute_peratom()
   // compute properties for each atom in group
   // use full neighbor list to count atoms less than cutoff
 
-  atomKK->sync(execution_space,X_MASK|V_MASK|TYPE_MASK|MASK_MASK);
+  atomKK->sync(execution_space,X_MASK|V_MASK|RMASS_MASK|TYPE_MASK|MASK_MASK);
   x = atomKK->k_x.view<DeviceType>();
   v = atomKK->k_v.view<DeviceType>();
+  rmass = atomKK->k_rmass.view<DeviceType>();
+  mass = atomKK->k_mass.view<DeviceType>();
+  type = atomKK->k_type.view<DeviceType>();
   mask = atomKK->k_mask.view<DeviceType>();
 
+  adof = domain->dimension;
+  mvv2e = force->mvv2e;
+  mv2d = force->mv2d;
+  boltz = force->boltz;
+
   Kokkos::deep_copy(d_result,0.0);
 
   copymode = 1;
@@ -125,8 +138,13 @@ template<class DeviceType>
 KOKKOS_INLINE_FUNCTION
 void ComputeAveSphereAtomKokkos<DeviceType>::operator()(TagComputeAveSphereAtom, const int &ii) const
 {
+  double massone_i,massone_j;
+
   const int i = d_ilist[ii];
   if (mask[i] & groupbit) {
+    if (rmass.data()) massone_i = rmass[i];
+    else massone_i = mass[type[i]];
+
     const X_FLOAT xtmp = x(i,0);
     const X_FLOAT ytmp = x(i,1);
     const X_FLOAT ztmp = x(i,2);
@@ -135,14 +153,17 @@ void ComputeAveSphereAtomKokkos<DeviceType>::operator()(TagComputeAveSphereAtom,
     // i atom contribution
 
     int count = 1;
-    double vsum[3];
+    double totalmass = massone_i;
-    vsum[0] = v(i,0);
+    double p[3];
-    vsum[1] = v(i,1);
+    p[0] = v(i,0)*massone_i;
-    vsum[2] = v(i,2);
+    p[1] = v(i,1)*massone_i;
+    p[2] = v(i,2)*massone_i;
 
     for (int jj = 0; jj < jnum; jj++) {
       int j = d_neighbors(i,jj);
       j &= NEIGHMASK;
+      if (rmass.data()) massone_j = rmass[j];
+      else massone_j = mass[type[j]];
 
       const F_FLOAT delx = x(j,0) - xtmp;
       const F_FLOAT dely = x(j,1) - ytmp;
@@ -150,44 +171,45 @@ void ComputeAveSphereAtomKokkos<DeviceType>::operator()(TagComputeAveSphereAtom,
       const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
       if (rsq < cutsq) {
         count++;
-        vsum[0] += v(j,0);
+        totalmass += massone_j;
-        vsum[1] += v(j,1);
+        p[0] += v(j,0)*massone_j;
-        vsum[2] += v(j,2);
+        p[1] += v(j,1)*massone_j;
+        p[2] += v(j,2)*massone_j;
       }
     }
 
-    double vavg[3];
+    double vcom[3];
-    vavg[0] = vsum[0]/count;
+    vcom[0] = p[0]/totalmass;
-    vavg[1] = vsum[1]/count;
+    vcom[1] = p[1]/totalmass;
-    vavg[2] = vsum[2]/count;
+    vcom[2] = p[2]/totalmass;
 
     // i atom contribution
 
-    count = 1;
     double vnet[3];
-    vnet[0] = v(i,0) - vavg[0];
+    vnet[0] = v(i,0) - vcom[0];
-    vnet[1] = v(i,1) - vavg[1];
+    vnet[1] = v(i,1) - vcom[1];
-    vnet[2] = v(i,2) - vavg[2];
+    vnet[2] = v(i,2) - vcom[2];
-    double ke_sum = vnet[0]*vnet[0] + vnet[1]*vnet[1] + vnet[2]*vnet[2];
+    double ke_sum = massone_i * (vnet[0]*vnet[0] + vnet[1]*vnet[1] + vnet[2]*vnet[2]);
 
     for (int jj = 0; jj < jnum; jj++) {
       int j = d_neighbors(i,jj);
       j &= NEIGHMASK;
+      if (rmass.data()) massone_j = rmass[j];
+      else massone_j = mass[type[j]];
 
       const F_FLOAT delx = x(j,0) - xtmp;
       const F_FLOAT dely = x(j,1) - ytmp;
       const F_FLOAT delz = x(j,2) - ztmp;
       const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
       if (rsq < cutsq) {
-        count++;
+        vnet[0] = v(j,0) - vcom[0];
-        vnet[0] = v(j,0) - vavg[0];
+        vnet[1] = v(j,1) - vcom[1];
-        vnet[1] = v(j,1) - vavg[1];
+        vnet[2] = v(j,2) - vcom[2];
-        vnet[2] = v(j,2) - vavg[2];
+        ke_sum += massone_j * (vnet[0]*vnet[0] + vnet[1]*vnet[1] + vnet[2]*vnet[2]);
-        ke_sum += vnet[0]*vnet[0] + vnet[1]*vnet[1] + vnet[2]*vnet[2];
       }
     }
-    double density = count/sphere_vol;
+    double density = mv2d*totalmass/volume;
-    double temp = ke_sum/3.0/count;
+    double temp = mvv2e*ke_sum/(adof*count*boltz);
     d_result(i,0) = density;
     d_result(i,1) = temp;
   }

src/KOKKOS/compute_ave_sphere_atom_kokkos.h

@@ -46,13 +46,18 @@ template <class DeviceType> class ComputeAveSphereAtomKokkos : public ComputeAve
   void operator()(TagComputeAveSphereAtom, const int &) const;
 
  private:
-  typename AT::t_x_array_randomread x;
+  double adof,mvv2e,mv2d,boltz;
-  typename AT::t_v_array_randomread v;
+
+  typename AT::t_x_array x;
+  typename AT::t_v_array v;
+  typename ArrayTypes<DeviceType>::t_float_1d rmass;
+  typename ArrayTypes<DeviceType>::t_float_1d mass;
+  typename ArrayTypes<DeviceType>::t_int_1d type;
   typename ArrayTypes<DeviceType>::t_int_1d mask;
 
   typename AT::t_neighbors_2d d_neighbors;
-  typename AT::t_int_1d_randomread d_ilist;
+  typename AT::t_int_1d d_ilist;
-  typename AT::t_int_1d_randomread d_numneigh;
+  typename AT::t_int_1d d_numneigh;
 
   DAT::tdual_float_2d k_result;
   typename AT::t_float_2d d_result;

src/KOKKOS/pair_pace_kokkos.cpp

@@ -534,8 +534,7 @@ void PairPACEKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
   }
 
   copymode = 1;
-  int newton_pair = force->newton_pair;
+  if (!force->newton_pair)
-  if (newton_pair == false)
     error->all(FLERR,"PairPACEKokkos requires 'newton on'");
 
   if (recursive)

src/KOKKOS/pair_reaxff_kokkos.cpp

@@ -2628,7 +2628,7 @@ int PairReaxFFKokkos<DeviceType>::preprocess_angular(int i, int itype, int j_sta
 template<class DeviceType>
 template<bool POPULATE>
 KOKKOS_INLINE_FUNCTION
-int PairReaxFFKokkos<DeviceType>::preprocess_torsion(int i, int /*itype*/, int itag,
+int PairReaxFFKokkos<DeviceType>::preprocess_torsion(int i, int /*itype*/, tagint itag,
   F_FLOAT xtmp, F_FLOAT ytmp, F_FLOAT ztmp, int j_start, int j_end, int location_torsion) const {
 
   // in reaxff_torsion_angles: j = i, k = j, i = k;

src/KOKKOS/pair_reaxff_kokkos.h

@@ -257,7 +257,7 @@ class PairReaxFFKokkos : public PairReaxFF {
   // Abstraction for counting and populating torsion intermediated
   template<bool POPULATE>
   KOKKOS_INLINE_FUNCTION
-  int preprocess_torsion(int, int, int, F_FLOAT, F_FLOAT, F_FLOAT, int, int, int) const;
+  int preprocess_torsion(int, int, tagint, F_FLOAT, F_FLOAT, F_FLOAT, int, int, int) const;
 
   template<int NEIGHFLAG, int EVFLAG>
   KOKKOS_INLINE_FUNCTION