Merge pull request #3460 from akohlmey/maintenance-2022-06-23

Second round of maintenance fixes and backports for the stable release

cmake/CMakeLists.txt

@@ -7,6 +7,10 @@ cmake_minimum_required(VERSION 3.10)
 if(POLICY CMP0074)
   cmake_policy(SET CMP0074 NEW)
 endif()
+# set policy to silence warnings about ignoring  ${CMAKE_REQUIRED_LIBRARIES} but use it
+if(POLICY CMP0075)
+  cmake_policy(SET CMP0075 NEW)
+endif()
 # set policy to silence warnings about missing executable permissions in
 # pythonx.y-config when cross-compiling. review occasionally if it may be set to NEW
 if(POLICY CMP0109)
@@ -366,9 +370,9 @@ pkg_depends(ELECTRODE KSPACE)
 
 # 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)
+find_package(OpenMP COMPONENTS CXX QUIET)
+if(OpenMP_CXX_FOUND)
+  check_omp_h_include()
   if(HAVE_OMP_H_INCLUDE)
     set(BUILD_OMP_DEFAULT ON)
   endif()
@@ -377,8 +381,8 @@ 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)
+  find_package(OpenMP COMPONENTS CXX REQUIRED)
+  check_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()

cmake/Modules/FindCythonize.cmake

@@ -22,7 +22,7 @@ endif()
 if(Python_EXECUTABLE)
   get_filename_component(_python_path ${Python_EXECUTABLE} PATH)
   find_program(Cythonize_EXECUTABLE
-    NAMES cythonize3 cythonize cythonize.bat
+    NAMES cythonize-${Python_VERSION_MAJOR}.${Python_VERSION_MINOR} cythonize3 cythonize cythonize.bat
     HINTS ${_python_path})
 endif()
 

cmake/Modules/LAMMPSInterfacePlugin.cmake

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

cmake/Modules/LAMMPSUtils.cmake

@@ -24,6 +24,21 @@ function(validate_option name values)
     endif()
 endfunction(validate_option)
 
+# helper function to check for usable omp.h header
+function(check_omp_h_include)
+  find_package(OpenMP COMPONENTS CXX QUIET)
+  if(OpenMP_CXX_FOUND)
+    set(CMAKE_REQUIRED_FLAGS ${OpenMP_CXX_FLAGS})
+    set(CMAKE_REQUIRED_INCLUDES ${OpenMP_CXX_INCLUDE_DIRS})
+    set(CMAKE_REQUIRED_LINK_OPTIONS ${OpenMP_CXX_FLAGS})
+    set(CMAKE_REQUIRED_LIBRARIES ${OpenMP_CXX_LIBRARIES})
+    check_include_file_cxx(omp.h _have_omp_h)
+  else()
+    set(_have_omp_h FALSE)
+  endif()
+  set(HAVE_OMP_H_INCLUDE ${_have_omp_h} PARENT_SCOPE)
+endfunction()
+
 # 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})
@@ -110,14 +125,16 @@ function(FetchPotentials pkgfolder potfolder)
       math(EXPR plusone "${blank}+1")
       string(SUBSTRING ${line} 0 ${blank} pot)
       string(SUBSTRING ${line} ${plusone} -1 sum)
-      if(EXISTS ${LAMMPS_POTENTIALS_DIR}/${pot})
+      if(EXISTS "${LAMMPS_POTENTIALS_DIR}/${pot}")
         file(MD5 "${LAMMPS_POTENTIALS_DIR}/${pot}" oldsum)
       endif()
       if(NOT sum STREQUAL oldsum)
-        message(STATUS "Checking external potential ${pot} from ${LAMMPS_POTENTIALS_URL}")
-        file(DOWNLOAD "${LAMMPS_POTENTIALS_URL}/${pot}.${sum}" "${CMAKE_BINARY_DIR}/${pot}"
+        message(STATUS "Downloading external potential ${pot} from ${LAMMPS_POTENTIALS_URL}")
+        string(MD5 TMP_EXT "${CMAKE_BINARY_DIR}")
+        file(DOWNLOAD "${LAMMPS_POTENTIALS_URL}/${pot}.${sum}" "${CMAKE_BINARY_DIR}/${pot}.${TMP_EXT}"
           EXPECTED_HASH MD5=${sum} SHOW_PROGRESS)
-        file(COPY "${CMAKE_BINARY_DIR}/${pot}" DESTINATION ${LAMMPS_POTENTIALS_DIR})
+        file(COPY "${CMAKE_BINARY_DIR}/${pot}.${TMP_EXT}" DESTINATION "${LAMMPS_POTENTIALS_DIR}")
+        file(RENAME "${LAMMPS_POTENTIALS_DIR}/${pot}.${TMP_EXT}" "${LAMMPS_POTENTIALS_DIR}/${pot}")
       endif()
     endforeach()
   endif()

cmake/Modules/Packages/KOKKOS.cmake

@@ -15,8 +15,9 @@ if(Kokkos_ENABLE_OPENMP)
   if(NOT BUILD_OMP)
     message(FATAL_ERROR "Must enable BUILD_OMP with Kokkos_ENABLE_OPENMP")
   else()
-    if(LAMMPS_OMP_COMPAT_LEVEL LESS 4)
-      message(FATAL_ERROR "Compiler must support OpenMP 4.0 or later with Kokkos_ENABLE_OPENMP")
+    # NVHPC does not seem to provide a detectable OpenMP version, but is far beyond version 3.1
+    if((OpenMP_CXX_VERSION VERSION_LESS 3.1) AND NOT (CMAKE_CXX_COMPILER_ID STREQUAL "NVHPC"))
+      message(FATAL_ERROR "Compiler must support OpenMP 3.1 or later with Kokkos_ENABLE_OPENMP")
     endif()
   endif()
 endif()
@@ -47,8 +48,8 @@ if(DOWNLOAD_KOKKOS)
   list(APPEND KOKKOS_LIB_BUILD_ARGS "-DCMAKE_CXX_EXTENSIONS=${CMAKE_CXX_EXTENSIONS}")
   list(APPEND KOKKOS_LIB_BUILD_ARGS "-DCMAKE_TOOLCHAIN_FILE=${CMAKE_TOOLCHAIN_FILE}")
   include(ExternalProject)
-  set(KOKKOS_URL "https://github.com/kokkos/kokkos/archive/3.6.00.tar.gz" CACHE STRING "URL for KOKKOS tarball")
-  set(KOKKOS_MD5 "b5c44ea961031795f434002cd7b31c20" CACHE STRING "MD5 checksum of KOKKOS tarball")
+  set(KOKKOS_URL "https://github.com/kokkos/kokkos/archive/3.6.01.tar.gz" CACHE STRING "URL for KOKKOS tarball")
+  set(KOKKOS_MD5 "0ec97fc0c356dd65bd2487defe81a7bf" CACHE STRING "MD5 checksum of KOKKOS tarball")
   mark_as_advanced(KOKKOS_URL)
   mark_as_advanced(KOKKOS_MD5)
   ExternalProject_Add(kokkos_build
@@ -72,7 +73,7 @@ if(DOWNLOAD_KOKKOS)
   add_dependencies(LAMMPS::KOKKOSCORE kokkos_build)
   add_dependencies(LAMMPS::KOKKOSCONTAINERS kokkos_build)
 elseif(EXTERNAL_KOKKOS)
-  find_package(Kokkos 3.6.00 REQUIRED CONFIG)
+  find_package(Kokkos 3.6.01 REQUIRED CONFIG)
   target_link_libraries(lammps PRIVATE Kokkos::kokkos)
   target_link_libraries(lmp PRIVATE Kokkos::kokkos)
 else()

cmake/Modules/Packages/ML-HDNNP.cmake

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

cmake/Modules/Packages/ML-IAP.cmake

@@ -2,7 +2,13 @@
 set(MLIAP_ENABLE_PYTHON_DEFAULT OFF)
 if(PKG_PYTHON)
   find_package(Cythonize QUIET)
-  if(Cythonize_FOUND)
+  if (CMAKE_VERSION VERSION_GREATER_EQUAL 3.14)
+    find_package(Python COMPONENTS NumPy QUIET)
+  else()
+    # assume we have NumPy
+    set(Python_NumPy_FOUND ON)
+  endif()
+  if(Cythonize_FOUND AND Python_NumPy_FOUND)
     set(MLIAP_ENABLE_PYTHON_DEFAULT ON)
   endif()
 endif()
@@ -11,6 +17,9 @@ option(MLIAP_ENABLE_PYTHON "Build ML-IAP package with Python support" ${MLIAP_EN
 
 if(MLIAP_ENABLE_PYTHON)
   find_package(Cythonize REQUIRED)
+  if (CMAKE_VERSION VERSION_GREATER_EQUAL 3.14)
+    find_package(Python COMPONENTS NumPy REQUIRED)
+  endif()
   if(NOT PKG_PYTHON)
     message(FATAL_ERROR "Must enable PYTHON package for including Python support in ML-IAP")
   endif()

cmake/Modules/Packages/PLUMED.cmake

@@ -47,15 +47,15 @@ if(DOWNLOAD_PLUMED)
   endif()
   message(STATUS "PLUMED download requested - we will build our own")
   if(PLUMED_MODE STREQUAL "STATIC")
-    set(PLUMED_BUILD_BYPRODUCTS "<INSTALL_DIR>/lib/libplumed.a")
+    set(PLUMED_BUILD_BYPRODUCTS "<INSTALL_DIR>/lib/${CMAKE_STATIC_LIBRARY_PREFIX}plumed${CMAKE_STATIC_LIBRARY_SUFFIX}")
   elseif(PLUMED_MODE STREQUAL "SHARED")
-    set(PLUMED_BUILD_BYPRODUCTS "<INSTALL_DIR>/lib/libplumed${CMAKE_SHARED_LIBRARY_SUFFIX};<INSTALL_DIR>/lib/libplumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
+    set(PLUMED_BUILD_BYPRODUCTS "<INSTALL_DIR>/lib/${CMAKE_SHARED_LIBRARY_PREFIX}plumed${CMAKE_SHARED_LIBRARY_SUFFIX};<INSTALL_DIR>/lib/${CMAKE_SHARED_LIBRARY_PREFIX}plumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
   elseif(PLUMED_MODE STREQUAL "RUNTIME")
-    set(PLUMED_BUILD_BYPRODUCTS "<INSTALL_DIR>/lib/libplumedWrapper.a")
+    set(PLUMED_BUILD_BYPRODUCTS "<INSTALL_DIR>/lib/${CMAKE_STATIC_LIBRARY_PREFIX}plumedWrapper${CMAKE_STATIC_LIBRARY_PREFIX}")
   endif()
 
-  set(PLUMED_URL "https://github.com/plumed/plumed2/releases/download/v2.7.4/plumed-src-2.7.4.tgz" CACHE STRING "URL for PLUMED tarball")
-  set(PLUMED_MD5 "858e0b6aed173748fc85b6bc8a9dcb3e" CACHE STRING "MD5 checksum of PLUMED tarball")
+  set(PLUMED_URL "https://github.com/plumed/plumed2/releases/download/v2.8.1/plumed-src-2.8.1.tgz" CACHE STRING "URL for PLUMED tarball")
+  set(PLUMED_MD5 "6bfe72ebdae63dc38a9ca27d9b0e08f8" CACHE STRING "MD5 checksum of PLUMED tarball")
 
   mark_as_advanced(PLUMED_URL)
   mark_as_advanced(PLUMED_MD5)
@@ -78,12 +78,12 @@ if(DOWNLOAD_PLUMED)
   add_library(LAMMPS::PLUMED UNKNOWN IMPORTED)
   add_dependencies(LAMMPS::PLUMED plumed_build)
   if(PLUMED_MODE STREQUAL "STATIC")
-    set_target_properties(LAMMPS::PLUMED PROPERTIES IMPORTED_LOCATION ${INSTALL_DIR}/lib/libplumed.a INTERFACE_LINK_LIBRARIES "${PLUMED_LINK_LIBS};${CMAKE_DL_LIBS}")
+    set_target_properties(LAMMPS::PLUMED PROPERTIES IMPORTED_LOCATION ${INSTALL_DIR}/lib/${CMAKE_STATIC_LIBRARY_PREFIX}plumed${CMAKE_STATIC_LIBRARY_SUFFIX} INTERFACE_LINK_LIBRARIES "${PLUMED_LINK_LIBS};${CMAKE_DL_LIBS}")
   elseif(PLUMED_MODE STREQUAL "SHARED")
-    set_target_properties(LAMMPS::PLUMED PROPERTIES IMPORTED_LOCATION ${INSTALL_DIR}/lib/libplumed${CMAKE_SHARED_LIBRARY_SUFFIX} INTERFACE_LINK_LIBRARIES "${INSTALL_DIR}/lib/libplumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX};${CMAKE_DL_LIBS}")
+    set_target_properties(LAMMPS::PLUMED PROPERTIES IMPORTED_LOCATION ${INSTALL_DIR}/lib/${CMAKE_SHARED_LIBRARY_PREFIX}plumed${CMAKE_SHARED_LIBRARY_SUFFIX} INTERFACE_LINK_LIBRARIES "${INSTALL_DIR}/lib/${CMAKE_SHARED_LIBRARY_PREFIX}plumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX};${CMAKE_DL_LIBS}")
   elseif(PLUMED_MODE STREQUAL "RUNTIME")
-    set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_COMPILE_DEFINITIONS "__PLUMED_DEFAULT_KERNEL=${INSTALL_DIR}/lib/libplumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
-    set_target_properties(LAMMPS::PLUMED PROPERTIES IMPORTED_LOCATION ${INSTALL_DIR}/lib/libplumedWrapper.a INTERFACE_LINK_LIBRARIES "${CMAKE_DL_LIBS}")
+    set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_COMPILE_DEFINITIONS "__PLUMED_DEFAULT_KERNEL=${INSTALL_DIR}/lib/${CMAKE_SHARED_LIBRARY_PREFIX}plumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
+    set_target_properties(LAMMPS::PLUMED PROPERTIES IMPORTED_LOCATION ${INSTALL_DIR}/lib/${CMAKE_STATIC_LIBRARY_PREFIX}plumedWrapper${CMAKE_STATIC_LIBRARY_SUFFIX} INTERFACE_LINK_LIBRARIES "${CMAKE_DL_LIBS}")
   endif()
   set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_INCLUDE_DIRECTORIES ${INSTALL_DIR}/include)
   file(MAKE_DIRECTORY ${INSTALL_DIR}/include)
@@ -96,7 +96,7 @@ else()
   elseif(PLUMED_MODE STREQUAL "SHARED")
     include(${PLUMED_LIBDIR}/plumed/src/lib/Plumed.cmake.shared)
   elseif(PLUMED_MODE STREQUAL "RUNTIME")
-    set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_COMPILE_DEFINITIONS "__PLUMED_DEFAULT_KERNEL=${PLUMED_LIBDIR}/libplumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
+    set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_COMPILE_DEFINITIONS "__PLUMED_DEFAULT_KERNEL=${PLUMED_LIBDIR}/${CMAKE_SHARED_LIBRARY_PREFIX}plumedKernel${CMAKE_SHARED_LIBRARY_SUFFIX}")
     include(${PLUMED_LIBDIR}/plumed/src/lib/Plumed.cmake.runtime)
   endif()
   set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_LINK_LIBRARIES "${PLUMED_LOAD}")

cmake/Modules/Packages/VTK.cmake

@@ -1,4 +1,9 @@
 find_package(VTK REQUIRED NO_MODULE)
-include(${VTK_USE_FILE})
 target_compile_definitions(lammps PRIVATE -DLAMMPS_VTK)
-target_link_libraries(lammps PRIVATE ${VTK_LIBRARIES})
+if (VTK_MAJOR_VERSION VERSION_LESS 9.0)
+  include(${VTK_USE_FILE})
+  target_link_libraries(lammps PRIVATE ${VTK_LIBRARIES})
+else()
+  target_link_libraries(lammps PRIVATE VTK::CommonCore VTK::IOCore VTK::CommonDataModel VTK::IOXML VTK::IOLegacy VTK::IOParallelXML)
+  vtk_module_autoinit(TARGETS lammps MODULES VTK::CommonCore VTK::IOCore VTK::CommonDataModel VTK::IOXML VTK::IOLegacy VTK::IOParallelXML)
+endif()

cmake/presets/clang.cmake

@@ -21,10 +21,3 @@ set(MPI_CXX "clang++" CACHE STRING "" FORCE)
 set(MPI_CXX_COMPILER "mpicxx" CACHE STRING "" FORCE)
 
 unset(HAVE_OMP_H_INCLUDE CACHE)
-set(OpenMP_C "clang" CACHE STRING "" FORCE)
-set(OpenMP_C_FLAGS "-fopenmp" CACHE STRING "" FORCE)
-set(OpenMP_C_LIB_NAMES "omp" CACHE STRING "" FORCE)
-set(OpenMP_CXX "clang++" CACHE STRING "" FORCE)
-set(OpenMP_CXX_FLAGS "-fopenmp" CACHE STRING "" FORCE)
-set(OpenMP_CXX_LIB_NAMES "omp" CACHE STRING "" FORCE)
-set(OpenMP_omp_LIBRARY "libomp.so" CACHE PATH "" FORCE)

cmake/presets/gcc.cmake

@@ -19,11 +19,3 @@ set(CMAKE_Fortran_FLAGS_DEBUG "-Wall -Og -g -std=f2003" CACHE STRING "" FORCE)
 set(CMAKE_Fortran_FLAGS_RELWITHDEBINFO "-g -O2 -DNDEBUG -std=f2003" CACHE STRING "" FORCE)
 set(CMAKE_Fortran_FLAGS_RELEASE "-O3 -DNDEBUG -std=f2003" CACHE STRING "" FORCE)
 unset(HAVE_OMP_H_INCLUDE CACHE)
-
-set(OpenMP_C "gcc" CACHE STRING "" FORCE)
-set(OpenMP_C_FLAGS "-fopenmp" CACHE STRING "" FORCE)
-set(OpenMP_C_LIB_NAMES "gomp" CACHE STRING "" FORCE)
-set(OpenMP_CXX "g++" CACHE STRING "" FORCE)
-set(OpenMP_CXX_FLAGS "-fopenmp" CACHE STRING "" FORCE)
-set(OpenMP_CXX_LIB_NAMES "gomp" CACHE STRING "" FORCE)
-set(OpenMP_omp_LIBRARY "libgomp.so" CACHE PATH "" FORCE)

cmake/presets/nvhpc.cmake

@@ -0,0 +1,9 @@
+# preset that will enable Nvidia HPC SDK compilers with support for MPI and OpenMP (on Linux boxes)
+
+set(CMAKE_CXX_COMPILER "nvc++" CACHE STRING "" FORCE)
+set(CMAKE_C_COMPILER "nvc" CACHE STRING "" FORCE)
+set(CMAKE_Fortran_COMPILER "nvfortran" CACHE STRING "" FORCE)
+set(MPI_CXX "nvc++" CACHE STRING "" FORCE)
+set(MPI_CXX_COMPILER "mpicxx" CACHE STRING "" FORCE)
+unset(HAVE_OMP_H_INCLUDE CACHE)
+

cmake/presets/pedantic.cmake

@@ -1,4 +1,4 @@
-# preset that will restore gcc/g++ with support for MPI and OpenMP (on Linux boxes)
+# preset that will set gcc/g++ with extra warnings enabled and support for MPI and OpenMP (on Linux boxes)
 
 set(CMAKE_CXX_COMPILER "g++" CACHE STRING "" FORCE)
 set(CMAKE_C_COMPILER "gcc" CACHE STRING "" FORCE)
@@ -17,10 +17,3 @@ set(MPI_Fortran "gfortran" CACHE STRING "" FORCE)
 set(MPI_Fortran_COMPILER "mpifort" CACHE STRING "" FORCE)
 unset(HAVE_OMP_H_INCLUDE CACHE)
 
-set(OpenMP_C "gcc" CACHE STRING "" FORCE)
-set(OpenMP_C_FLAGS "-fopenmp" CACHE STRING "" FORCE)
-set(OpenMP_C_LIB_NAMES "gomp" CACHE STRING "" FORCE)
-set(OpenMP_CXX "g++" CACHE STRING "" FORCE)
-set(OpenMP_CXX_FLAGS "-fopenmp" CACHE STRING "" FORCE)
-set(OpenMP_CXX_LIB_NAMES "gomp" CACHE STRING "" FORCE)
-set(OpenMP_omp_LIBRARY "libgomp.so" CACHE PATH "" FORCE)

cmake/presets/pgi.cmake

@@ -7,10 +7,3 @@ set(MPI_CXX "pgc++" CACHE STRING "" FORCE)
 set(MPI_CXX_COMPILER "mpicxx" CACHE STRING "" FORCE)
 unset(HAVE_OMP_H_INCLUDE CACHE)
 
-set(OpenMP_C "pgcc" CACHE STRING "" FORCE)
-set(OpenMP_C_FLAGS "-mp" CACHE STRING "" FORCE)
-set(OpenMP_C_LIB_NAMES "omp" CACHE STRING "" FORCE)
-set(OpenMP_CXX "pgc++" CACHE STRING "" FORCE)
-set(OpenMP_CXX_FLAGS "-mp" CACHE STRING "" FORCE)
-set(OpenMP_CXX_LIB_NAMES "omp" CACHE STRING "" FORCE)
-set(OpenMP_omp_LIBRARY "libomp.so" CACHE PATH "" FORCE)

doc/Makefile

@@ -38,16 +38,14 @@ endif
 # override settings for PIP commands
 # PIP_OPTIONS = --cert /etc/pki/ca-trust/extracted/openssl/ca-bundle.trust.crt --proxy http://proxy.mydomain.org
 
-#SPHINXEXTRA = -j $(shell $(PYTHON) -c 'import multiprocessing;print(multiprocessing.cpu_count())') $(shell test -f $(BUILDDIR)/doxygen/xml/run.stamp && printf -- "-E")
-
 # temporarily disable caching so that the hack for the sphinx-tabs extensions to get proper non-html output works
-SPHINXEXTRA = -E -j $(shell $(PYTHON) -c 'import multiprocessing;print(multiprocessing.cpu_count())')
+SPHINXEXTRA = -j $(shell $(PYTHON) -c 'import multiprocessing;print(multiprocessing.cpu_count())')
 
 # grab list of sources from doxygen config file.
 # we only want to use explicitly listed files.
 DOXYFILES      = $(shell sed -n -e 's/\#.*$$//' -e '/^ *INPUT \+=/,/^[A-Z_]\+ \+=/p' doxygen/Doxyfile.in | sed -e 's/@LAMMPS_SOURCE_DIR@/..\/src/g' -e 's/\\//g' -e 's/ \+/ /' -e 's/[A-Z_]\+ \+= *\(YES\|NO\|\)//') 
 
-.PHONY: help clean-all clean clean-spelling epub mobi rst html pdf spelling anchor_check style_check char_check xmlgen fasthtml
+.PHONY: help clean-all clean clean-spelling epub mobi html pdf spelling anchor_check style_check char_check xmlgen fasthtml
 
 # ------------------------------------------
 
@@ -89,6 +87,8 @@ html: xmlgen $(VENV) $(SPHINXCONFIG)/conf.py $(ANCHORCHECK) $(MATHJAX)
 	@$(MAKE) $(MFLAGS) -C graphviz all
 	@(\
 		. $(VENV)/bin/activate ; env PYTHONWARNINGS= \
+		sphinx-build -E $(SPHINXEXTRA) -b html -c $(SPHINXCONFIG) -d $(BUILDDIR)/doctrees $(RSTDIR) html ;\
+		touch $(RSTDIR)/Fortran.rst ;\
 		sphinx-build $(SPHINXEXTRA) -b html -c $(SPHINXCONFIG) -d $(BUILDDIR)/doctrees $(RSTDIR) html ;\
 		ln -sf Manual.html html/index.html;\
 		rm -f $(BUILDDIR)/doxygen/xml/run.stamp;\
@@ -114,7 +114,9 @@ fasthtml: xmlgen $(VENV) $(SPHINXCONFIG)/conf.py $(ANCHORCHECK) $(MATHJAX)
 	@mkdir -p fasthtml
 	@(\
 		. $(VENV)/bin/activate ; env PYTHONWARNINGS= \
-		sphinx-build -j 4 -b html -c $(SPHINXCONFIG) -d $(BUILDDIR)/fasthtml/doctrees $(RSTDIR) fasthtml ;\
+		sphinx-build $(SPHINXEXTRA) -b html -c $(SPHINXCONFIG) -d $(BUILDDIR)/fasthtml/doctrees $(RSTDIR) fasthtml ;\
+		touch $(RSTDIR)/Fortran.rst ;\
+		sphinx-build $(SPHINXEXTRA) -b html -c $(SPHINXCONFIG) -d $(BUILDDIR)/fasthtml/doctrees $(RSTDIR) fasthtml ;\
 		deactivate ;\
 	)
 	@rm -rf fasthtml/_sources
@@ -144,6 +146,8 @@ epub: xmlgen $(VENV) $(SPHINXCONFIG)/conf.py $(ANCHORCHECK)
 	@cp src/JPG/*.* epub/JPG
 	@(\
 		. $(VENV)/bin/activate ;\
+		sphinx-build -E $(SPHINXEXTRA) -b epub -c $(SPHINXCONFIG) -d $(BUILDDIR)/doctrees $(RSTDIR) epub ;\
+		touch $(RSTDIR)/Fortran.rst ;\
 		sphinx-build $(SPHINXEXTRA) -b epub -c $(SPHINXCONFIG) -d $(BUILDDIR)/doctrees $(RSTDIR) epub ;\
 		rm -f $(BUILDDIR)/doxygen/xml/run.stamp;\
 		deactivate ;\
@@ -163,7 +167,9 @@ pdf: xmlgen $(VENV) $(SPHINXCONFIG)/conf.py $(ANCHORCHECK)
 	@if [ "$(HAS_PDFLATEX)" == "NO" ] ; then echo "PDFLaTeX or latexmk were not found! Please check README for further instructions" 1>&2; exit 1; fi
 	@(\
 		. $(VENV)/bin/activate ; env PYTHONWARNINGS= \
-		sphinx-build $(SPHINXEXTRA) -b latex -c $(SPHINXCONFIG) -d $(BUILDDIR)/doctrees $(RSTDIR) latex ;\
+		sphinx-build -E $(SPHINXEXTRA) -b latex -c $(SPHINXCONFIG) -d $(BUILDDIR)/doctrees $(RSTDIR) latex ;\
+		touch $(RSTDIR)/Fortran.rst ;\
+		sphinx-build  $(SPHINXEXTRA) -b latex -c $(SPHINXCONFIG) -d $(BUILDDIR)/doctrees $(RSTDIR) latex ;\
 		rm -f $(BUILDDIR)/doxygen/xml/run.stamp;\
 		echo "############################################" ;\
 		rst_anchor_check src/*.rst ;\
@@ -242,7 +248,6 @@ $(MATHJAX):
 $(ANCHORCHECK): $(VENV)
 	@( \
 		. $(VENV)/bin/activate; \
-		(cd utils/converters;\
-		python setup.py develop);\
+		pip $(PIP_OPTIONS) install -e utils/converters;\
 		deactivate;\
 	)

doc/lammps.1

@@ -161,7 +161,7 @@ list references for specific cite-able features used during a
 run.
 .TP
 \fB\-pk <style> [options]\fR or \fB\-package <style> [options]\fR
-Invoke the \fBpackage\R command with <style> and optional arguments.
+Invoke the \fBpackage\fR command with <style> and optional arguments.
 The syntax is the same as if the command appeared in an input script.
 For example "-pk gpu 2" is the same as "package gpu 2" in the input
 script. The possible styles and options are discussed in the

doc/src/Build_development.rst

@@ -147,6 +147,16 @@ compile and will download and compile a specific recent version of the
 `Googletest <https://github.com/google/googletest/>`_ C++ test framework
 for implementing the tests.
 
+.. admonition:: Software version requirements for testing
+   :class: note
+
+   The compiler and library version requirements for the testing
+   framework are more strict than for the main part of LAMMPS.  For
+   example the default GNU C++ and Fortran compilers of RHEL/CentOS 7.x
+   (version 4.8.x) are not sufficient.  The CMake configuration will try
+   to detect compatible versions and either skip incompatible tests or
+   stop with an error.
+
 After compilation is complete, the unit testing is started in the build
 folder using the ``ctest`` command, which is part of the CMake software.
 The output of this command will be looking something like this::

doc/src/Build_extras.rst

@@ -123,6 +123,7 @@ CMake build
    -D GPU_API=value             # value = opencl (default) or cuda or hip
    -D GPU_PREC=value            # precision setting
                                 # value = double or mixed (default) or single
+   -D HIP_PATH                  # path to HIP installation. Must be set if GPU_API=HIP
    -D GPU_ARCH=value            # primary GPU hardware choice for GPU_API=cuda
                                 # value = sm_XX, see below
                                 # default is sm_50
@@ -215,7 +216,7 @@ LAMMPS code.  This also applies to the ``-DLAMMPS_BIGBIG``\ ,
 Makefile you use.
 
 You can also build the library in one step from the ``lammps/src`` dir,
-using a command like these, which simply invoke the ``lib/gpu/Install.py``
+using a command like these, which simply invokes the ``lib/gpu/Install.py``
 script with the specified args:
 
 .. code-block:: bash
@@ -295,7 +296,7 @@ detailed information is available at:
 
 In addition to installing the KIM API, it is also necessary to install the
 library of KIM models (interatomic potentials).
-See `Obtaining KIM Models <http://openkim.org/doc/usage/obtaining-models>`_ to
+See `Obtaining KIM Models <https://openkim.org/doc/usage/obtaining-models>`_ to
 learn how to install a pre-build binary of the OpenKIM Repository of Models.
 See the list of all KIM models here: https://openkim.org/browse/models
 
@@ -331,7 +332,7 @@ minutes to hours) to build.  Of course you only need to do that once.)
       You can download and build the KIM library manually if you prefer;
       follow the instructions in ``lib/kim/README``.  You can also do
       this in one step from the lammps/src directory, using a command like
-      these, which simply invoke the ``lib/kim/Install.py`` script with
+      these, which simply invokes the ``lib/kim/Install.py`` script with
       the specified args.
 
       .. code-block:: bash
@@ -413,7 +414,7 @@ Enabling the extra unit tests have some requirements,
   ``EAM_Dynamo_MendelevAckland_2007v3_Zr__MO_004835508849_000``,
   ``EAM_Dynamo_ErcolessiAdams_1994_Al__MO_123629422045_005``, and
   ``LennardJones612_UniversalShifted__MO_959249795837_003`` KIM models.
-  See `Obtaining KIM Models <http://openkim.org/doc/usage/obtaining-models>`_
+  See `Obtaining KIM Models <https://openkim.org/doc/usage/obtaining-models>`_
   to learn how to install a pre-built binary of the OpenKIM Repository of
   Models or see
   `Installing KIM Models <https://openkim.org/doc/usage/obtaining-models/#installing_models>`_
@@ -905,7 +906,7 @@ more details.
       You can download and build the MS-CG library manually if you
       prefer; follow the instructions in ``lib/mscg/README``\ .  You can
       also do it in one step from the ``lammps/src`` dir, using a
-      command like these, which simply invoke the
+      command like these, which simply invokes the
       ``lib/mscg/Install.py`` script with the specified args:
 
       .. code-block:: bash
@@ -962,7 +963,7 @@ POEMS package
       ``lib/poems``\ .  You can do this manually if you prefer; follow
       the instructions in ``lib/poems/README``\ .  You can also do it in
       one step from the ``lammps/src`` dir, using a command like these,
-      which simply invoke the ``lib/poems/Install.py`` script with the
+      which simply invokes the ``lib/poems/Install.py`` script with the
       specified args:
 
       .. code-block:: bash
@@ -1025,7 +1026,7 @@ VORONOI package
 -----------------------------
 
 To build with this package, you must download and build the
-`Voro++ library <http://math.lbl.gov/voro++>`_ or install a
+`Voro++ library <https://math.lbl.gov/voro++>`_ or install a
 binary package provided by your operating system.
 
 .. tabs::
@@ -1051,7 +1052,7 @@ binary package provided by your operating system.
       You can download and build the Voro++ library manually if you
       prefer; follow the instructions in ``lib/voronoi/README``.  You
       can also do it in one step from the ``lammps/src`` dir, using a
-      command like these, which simply invoke the
+      command like these, which simply invokes the
       ``lib/voronoi/Install.py`` script with the specified args:
 
       .. code-block:: bash
@@ -1130,7 +1131,7 @@ The ATC package requires the MANYBODY package also be installed.
       ``lib/atc``.  You can do this manually if you prefer; follow the
       instructions in ``lib/atc/README``.  You can also do it in one
       step from the ``lammps/src`` dir, using a command like these,
-      which simply invoke the ``lib/atc/Install.py`` script with the
+      which simply invokes the ``lib/atc/Install.py`` script with the
       specified args:
 
       .. code-block:: bash
@@ -1181,7 +1182,7 @@ AWPMD package
       ``lib/awpmd``.  You can do this manually if you prefer; follow the
       instructions in ``lib/awpmd/README``.  You can also do it in one
       step from the ``lammps/src`` dir, using a command like these,
-      which simply invoke the ``lib/awpmd/Install.py`` script with the
+      which simply invokes the ``lib/awpmd/Install.py`` script with the
       specified args:
 
       .. code-block:: bash
@@ -1244,7 +1245,7 @@ be built for the most part with all major versions of the C++ language.
 
       In general, it is safer to use build setting consistent with the
       rest of LAMMPS.  This is best carried out from the LAMMPS src
-      directory using a command like these, which simply invoke the
+      directory using a command like these, which simply invokes the
       ``lib/colvars/Install.py`` script with the specified args:
 
       .. code-block:: bash
@@ -1285,20 +1286,30 @@ This package depends on the KSPACE package.
 
    .. tab:: CMake build
 
-      No additional settings are needed besides ``-D PKG_KSPACE=yes`` and ``-D
-      PKG_ELECTRODE=yes``.
+      No additional settings are needed besides ``-D PKG_KSPACE=yes`` and
+      ``-D PKG_ELECTRODE=yes``.
 
    .. tab:: Traditional make
 
-      The package is activated with ``make yes-KSPACE`` and ``make
-      yes-ELECTRODE``
+      Before building LAMMPS, you must configure the ELECTRODE support
+      libraries and settings in ``lib/electrode``.  You can do this
+      manually, if you prefer, or do it in one step from the
+      ``lammps/src`` dir, using a command like these, which simply
+      invokes the ``lib/electrode/Install.py`` script with the specified
+      args:
+
+      .. code-block:: bash
+
+         $ make lib-electrode                   # print help message
+         $ make lib-electrode args="-m serial"  # build with GNU g++ compiler and MPI STUBS (settings as with "make serial")
+         $ make lib-electrode args="-m mpi"     # build with default MPI compiler (settings as with "make mpi")
 
 
-      Note that the ``Makefile.lammps`` file has settings for the BLAS and
-      LAPACK linear algebra libraries.  As explained in ``lib/awpmd/README``
-      these can either exist on your system, or you can use the files provided
-      in ``lib/linalg``.  In the latter case you also need to build the library
-      in ``lib/linalg`` with a command like these:
+      Note that the ``Makefile.lammps`` file has settings for the BLAS
+      and LAPACK linear algebra libraries.  These can either exist on
+      your system, or you can use the files provided in ``lib/linalg``.
+      In the latter case you also need to build the library in
+      ``lib/linalg`` with a command like these:
 
       .. code-block:: bash
 
@@ -1307,6 +1318,9 @@ This package depends on the KSPACE package.
          $ make lib-linalg args="-m mpi"       # build with default MPI Fortran compiler (settings as with "make mpi")
          $ make lib-linalg args="-m gfortran"  # build with GNU Fortran compiler
 
+      The package itself is activated with ``make yes-KSPACE`` and
+      ``make yes-ELECTRODE``
+
 ----------
 
 .. _ml-pace:
@@ -1506,7 +1520,7 @@ the HDF5 library.
       ``lib/h5md``.  You can do this manually if you prefer; follow the
       instructions in ``lib/h5md/README``.  You can also do it in one
       step from the ``lammps/src`` dir, using a command like these,
-      which simply invoke the ``lib/h5md/Install.py`` script with the
+      which simply invokes the ``lib/h5md/Install.py`` script with the
       specified args:
 
       .. code-block:: bash
@@ -1562,7 +1576,7 @@ details please see ``lib/hdnnp/README`` and the `n2p2 build documentation
       You can download and build the *n2p2* library manually if you prefer;
       follow the instructions in ``lib/hdnnp/README``\ . You can also do it in
       one step from the ``lammps/src`` dir, using a command like these, which
-      simply invoke the ``lib/hdnnp/Install.py`` script with the specified args:
+      simply invokes the ``lib/hdnnp/Install.py`` script with the specified args:
 
       .. code-block:: bash
 
@@ -1699,7 +1713,7 @@ they will be downloaded the first time this package is installed.
       Before building LAMMPS, you must build the *mesont* library in
       ``lib/mesont``\ .  You can also do it in one step from the
       ``lammps/src`` dir, using a command like these, which simply
-      invoke the ``lib/mesont/Install.py`` script with the specified
+      invokes the ``lib/mesont/Install.py`` script with the specified
       args:
 
       .. code-block:: bash
@@ -1868,7 +1882,7 @@ verified to work in February 2020 with Quantum Espresso versions 6.3 to
       ``lib/qmmm``.  You can do this manually if you prefer; follow the
       first two steps explained in ``lib/qmmm/README``.  You can also do
       it in one step from the ``lammps/src`` dir, using a command like
-      these, which simply invoke the ``lib/qmmm/Install.py`` script with
+      these, which simply invokes the ``lib/qmmm/Install.py`` script with
       the specified args:
 
       .. code-block:: bash
@@ -1913,9 +1927,9 @@ within CMake will download the non-commercial use version.
 
       .. code-block:: bash
 
-         -D DOWNLOAD_QUIP=value   # download OpenKIM API v2 for build, value = no (default) or yes
-         -D QUIP_LIBRARY=path     # path to libquip.a (only needed if a custom location)
-
+         -D DOWNLOAD_QUIP=value       # download QUIP library for build, value = no (default) or yes
+         -D QUIP_LIBRARY=path         # path to libquip.a (only needed if a custom location)
+ 
       CMake will try to download and build the QUIP library from GitHub, if it is not
       found on the local machine. This requires to have git installed. It will use the same compilers
       and flags as used for compiling LAMMPS.  Currently this is only supported for the GNU and the
@@ -1965,7 +1979,7 @@ To build with this package, you must download and build the
       You can download and build the ScaFaCoS library manually if you
       prefer; follow the instructions in ``lib/scafacos/README``.  You
       can also do it in one step from the ``lammps/src`` dir, using a
-      command like these, which simply invoke the
+      command like these, which simply invokes the
       ``lib/scafacos/Install.py`` script with the specified args:
 
       .. code-block:: bash
@@ -2009,7 +2023,7 @@ Eigen3 is a template library, so you do not need to build it.
       You can download the Eigen3 library manually if you prefer; follow
       the instructions in ``lib/smd/README``.  You can also do it in one
       step from the ``lammps/src`` dir, using a command like these,
-      which simply invoke the ``lib/smd/Install.py`` script with the
+      which simply invokes the ``lib/smd/Install.py`` script with the
       specified args:
 
       .. code-block:: bash

doc/src/Commands_parse.rst

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

doc/src/Errors_messages.rst

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

doc/src/Errors_warnings.rst

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

doc/src/Howto_viscosity.rst

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

doc/src/Intro_citing.rst

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

doc/src/Manual.rst

@@ -49,12 +49,12 @@ descriptions of all commands included in the LAMMPS code.
 
 ----------
 
+.. _user_documentation:
 
 ************
 User Guide
 ************
 
-.. _user_documentation:
 .. toctree::
    :maxdepth: 2
    :numbered: 3
@@ -75,11 +75,12 @@ User Guide
    Errors
 
 
+.. _programmer_documentation:
+
 ******************
 Programmer Guide
 ******************
 
-.. _programmer_documentation:
 .. toctree::
    :maxdepth: 2
    :numbered: 3

doc/src/Modify_fix.rst

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

doc/src/bond_bpm_rotational.rst

@@ -38,8 +38,7 @@ Examples
    bond_style bpm/rotational
    bond_coeff 1 1.0 0.2 0.02 0.02 0.20 0.04 0.04 0.04 0.1 0.02 0.002 0.002
 
-   bond_style bpm/rotational myfix 1000 time id1 id2
-   fix myfix all store/local 1000 3
+   bond_style bpm/rotational store/local myfix 1000 time id1 id2
    dump 1 all local 1000 dump.broken f_myfix[1] f_myfix[2] f_myfix[3]
    dump_modify 1 write_header no
 

doc/src/compute_damage_atom.rst

@@ -6,7 +6,7 @@ compute damage/atom command
 Syntax
 """"""
 
-.. parsed-literal::
+.. code-block:: LAMMPS
 
    compute ID group-ID damage/atom
 
@@ -24,16 +24,18 @@ Description
 """""""""""
 
 Define a computation that calculates the per-atom damage for each atom
-in a group.  This is a quantity relevant for :doc:`Peridynamics models <pair_peri>`.  See `this document <PDF/PDLammps_overview.pdf>`_
-for an overview of LAMMPS commands for Peridynamics modeling.
+in a group.  This is a quantity relevant for :doc:`Peridynamics models
+<pair_peri>`.  See `this document <PDF/PDLammps_overview.pdf>`_ for an
+overview of LAMMPS commands for Peridynamics modeling.
 
 The "damage" of a Peridynamics particles is based on the bond breakage
 between the particle and its neighbors.  If all the bonds are broken
 the particle is considered to be fully damaged.
 
-See the `PDLAMMPS user guide <http://www.sandia.gov/~mlparks/papers/PDLAMMPS.pdf>`_ for a formal
-definition of "damage" and more details about Peridynamics as it is
-implemented in LAMMPS.
+See the `PDLAMMPS user guide
+<https://download.lammps.org/pdfs/PDLAMMPS_user_guide.pdf>`_ for a
+formal definition of "damage" and more details about Peridynamics as it
+is implemented in LAMMPS.
 
 This command can be used with all the Peridynamic pair styles.
 

doc/src/compute_dilatation_atom.rst

@@ -6,12 +6,12 @@ compute dilatation/atom command
 Syntax
 """"""
 
-.. parsed-literal::
+.. code-block:: LAMMPS
 
    compute ID group-ID dilatation/atom
 
 * ID, group-ID are documented in compute command
-* dilation/atom = style name of this compute command
+* dilatation/atom = style name of this compute command
 
 Examples
 """"""""
@@ -24,21 +24,23 @@ Description
 """""""""""
 
 Define a computation that calculates the per-atom dilatation for each
-atom in a group.  This is a quantity relevant for :doc:`Peridynamics models <pair_peri>`.  See `this document <PDF/PDLammps_overview.pdf>`_
+atom in a group.  This is a quantity relevant for :doc:`Peridynamics
+models <pair_peri>`.  See `this document <PDF/PDLammps_overview.pdf>`_
 for an overview of LAMMPS commands for Peridynamics modeling.
 
 For small deformation, dilatation of is the measure of the volumetric
 strain.
 
-The dilatation "theta" for each peridynamic particle I is calculated
-as a sum over its neighbors with unbroken bonds, where the
-contribution of the IJ pair is a function of the change in bond length
-(versus the initial length in the reference state), the volume
-fraction of the particles and an influence function.  See the
-`PDLAMMPS user guide <http://www.sandia.gov/~mlparks/papers/PDLAMMPS.pdf>`_ for a formal
-definition of dilatation.
+The dilatation :math:`\theta` for each peridynamic particle :math:`i` is
+calculated as a sum over its neighbors with unbroken bonds, where the
+contribution of the :math:`ij` pair is a function of the change in bond
+length (versus the initial length in the reference state), the volume
+fraction of the particles and an influence function.  See the `PDLAMMPS
+user guide <https://download.lammps.org/pdfs/PDLAMMPS_user_guide.pdf>`_
+for a formal definition of dilatation.
 
-This command can only be used with a subset of the Peridynamic :doc:`pair styles <pair_peri>`: peri/lps, peri/ves and peri/eps.
+This command can only be used with a subset of the Peridynamic
+:doc:`pair styles <pair_peri>`: peri/lps, peri/ves and peri/eps.
 
 The dilatation value will be 0.0 for atoms not in the specified
 compute group.

doc/src/dump.rst

@@ -27,6 +27,9 @@ dump command
 :doc:`dump custom/adios <dump_adios>` command
 =============================================
 
+:doc:`dump cfg/uef <dump_cfg_uef>` command
+==========================================
+
 Syntax
 """"""
 
@@ -36,7 +39,7 @@ Syntax
 
 * ID = user-assigned name for the dump
 * group-ID = ID of the group of atoms to be dumped
-* style = *atom* or *atom/gz* or *atom/zstd or *atom/mpiio* or *cfg* or *cfg/gz* or *cfg/zstd* or *cfg/mpiio* or *custom* or *custom/gz* or *custom/zstd* or *custom/mpiio* or *dcd* or *h5md* or *image* or *local* or *local/gz* or *local/zstd* or *molfile* or *movie* or *netcdf* or *netcdf/mpiio* or *vtk* or *xtc* or *xyz* or *xyz/gz* or *xyz/zstd* or *xyz/mpiio* or *yaml*
+* style = *atom* or *atom/gz* or *atom/zstd or *atom/mpiio* or *cfg* or *cfg/gz* or *cfg/zstd* or *cfg/mpiio* or *cfg/uef* or *custom* or *custom/gz* or *custom/zstd* or *custom/mpiio* or *dcd* or *h5md* or *image* or *local* or *local/gz* or *local/zstd* or *molfile* or *movie* or *netcdf* or *netcdf/mpiio* or *vtk* or *xtc* or *xyz* or *xyz/gz* or *xyz/zstd* or *xyz/mpiio* or *yaml*
 * N = dump every this many timesteps
 * file = name of file to write dump info to
 * args = list of arguments for a particular style
@@ -47,22 +50,23 @@ Syntax
        *atom/gz* args = none
        *atom/zstd* args = none
        *atom/mpiio* args = none
-       *atom/adios* args = none,  discussed on :doc:`dump atom/adios <dump_adios>` doc page
+       *atom/adios* args = none,  discussed on :doc:`dump atom/adios <dump_adios>` page
        *cfg* args = same as *custom* args, see below
        *cfg/gz* args = same as *custom* args, see below
        *cfg/zstd* args = same as *custom* args, see below
        *cfg/mpiio* args = same as *custom* args, see below
+       *cfg/uef* args = same as *custom* args, discussed on :doc:`dump cfg/uef <dump_cfg_uef>` page
        *custom*, *custom/gz*, *custom/zstd*, *custom/mpiio* args = see below
-       *custom/adios* args = same as *custom* args, discussed on :doc:`dump custom/adios <dump_adios>` doc page
+       *custom/adios* args = same as *custom* args, discussed on :doc:`dump custom/adios <dump_adios>` page
        *dcd* args = none
-       *h5md* args = discussed on :doc:`dump h5md <dump_h5md>` doc page
-       *image* args = discussed on :doc:`dump image <dump_image>` doc page
+       *h5md* args = discussed on :doc:`dump h5md <dump_h5md>` page
+       *image* args = discussed on :doc:`dump image <dump_image>` page
        *local*, *local/gz*, *local/zstd* args = see below
-       *molfile* args = discussed on :doc:`dump molfile <dump_molfile>` doc page
-       *movie* args = discussed on :doc:`dump image <dump_image>` doc page
-       *netcdf* args = discussed on :doc:`dump netcdf <dump_netcdf>` doc page
-       *netcdf/mpiio* args = discussed on :doc:`dump netcdf <dump_netcdf>` doc page
-       *vtk* args = same as *custom* args, see below, also :doc:`dump vtk <dump_vtk>` doc page
+       *molfile* args = discussed on :doc:`dump molfile <dump_molfile>` page
+       *movie* args = discussed on :doc:`dump image <dump_image>` page
+       *netcdf* args = discussed on :doc:`dump netcdf <dump_netcdf>` page
+       *netcdf/mpiio* args = discussed on :doc:`dump netcdf <dump_netcdf>` page
+       *vtk* args = same as *custom* args, see below, also :doc:`dump vtk <dump_vtk>` page
        *xtc* args = none
        *xyz* args = none
        *xyz/gz* args = none
@@ -155,7 +159,7 @@ timesteps in one of several styles.  The *image* and *movie* styles are
 the exception: the *image* style renders a JPG, PNG, or PPM image file
 of the atom configuration every N timesteps while the *movie* style
 combines and compresses them into a movie file; both are discussed in
-detail on the :doc:`dump image <dump_image>` doc page.  The timesteps on
+detail on the :doc:`dump image <dump_image>` page.  The timesteps on
 which dump output is written can also be controlled by a variable.
 See the :doc:`dump_modify every <dump_modify>` command.
 
@@ -194,7 +198,7 @@ or multiple smaller files).
 For the *atom*, *custom*, *cfg*, and *local* styles, sorting is off by
 default.  For the *dcd*, *xtc*, *xyz*, and *molfile* styles, sorting
 by atom ID is on by default. See the :doc:`dump_modify <dump_modify>`
-doc page for details.
+page for details.
 
 The *atom/gz*, *cfg/gz*, *custom/gz*, *local/gz*, and *xyz/gz* styles
 are identical in command syntax to the corresponding styles without
@@ -204,7 +208,7 @@ alternative approach to writing compressed files via a pipe, as done
 by the regular dump styles, which may be required on clusters where
 the interface to the high-speed network disallows using the fork()
 library call (which is needed for a pipe).  For the remainder of this
-doc page, you should thus consider the *atom* and *atom/gz* styles
+page, you should thus consider the *atom* and *atom/gz* styles
 (etc) to be inter-changeable, with the exception of the required
 filename suffix.
 
@@ -218,7 +222,7 @@ As explained below, the *atom/mpiio*, *cfg/mpiio*, *custom/mpiio*, and
 *xyz/mpiio* styles are identical in command syntax and in the format
 of the dump files they create, to the corresponding styles without
 "mpiio", except the single dump file they produce is written in
-parallel via the MPI-IO library.  For the remainder of this doc page,
+parallel via the MPI-IO library.  For the remainder of this page,
 you should thus consider the *atom* and *atom/mpiio* styles (etc) to
 be inter-changeable.  The one exception is how the filename is
 specified for the MPI-IO styles, as explained below.
@@ -664,7 +668,7 @@ so that each value is 0.0 to 1.0.  If the simulation box is triclinic
 (tilted), then all atom coords will still be between 0.0 and 1.0.
 I.e. actual unscaled (x,y,z) = xs\*A + ys\*B + zs\*C, where (A,B,C) are
 the non-orthogonal vectors of the simulation box edges, as discussed
-on the :doc:`Howto triclinic <Howto_triclinic>` doc page.
+on the :doc:`Howto triclinic <Howto_triclinic>` page.
 
 Use *xu*, *yu*, *zu* if you want the coordinates "unwrapped" by the
 image flags for each atom.  Unwrapped means that if the atom has
@@ -787,7 +791,7 @@ more info.
 The *atom/mpiio*, *cfg/mpiio*, *custom/mpiio*, and *xyz/mpiio* styles
 are part of the MPIIO package.  They are only enabled if LAMMPS was
 built with that package.  See the :doc:`Build package <Build_package>`
-doc page for more info.
+page for more info.
 
 The *xtc*, *dcd* and *yaml* styles are part of the EXTRA-DUMP package.
 They are only enabled if LAMMPS was built with that package.  See the
@@ -797,12 +801,12 @@ Related commands
 """"""""""""""""
 
 :doc:`dump atom/adios <dump_adios>`, :doc:`dump custom/adios <dump_adios>`,
-:doc:`dump h5md <dump_h5md>`, :doc:`dump image <dump_image>`,
-:doc:`dump molfile <dump_molfile>`, :doc:`dump_modify <dump_modify>`,
-:doc:`undump <undump>`
+:doc:`dump cfg/uef <dump_cfg_uef>`, :doc:`dump h5md <dump_h5md>`,
+:doc:`dump image <dump_image>`, :doc:`dump molfile <dump_molfile>`,
+:doc:`dump_modify <dump_modify>`, :doc:`undump <undump>`, :doc:`write_dump <write_dump>`
 
 Default
 """""""
 
 The defaults for the *image* and *movie* styles are listed on the
-:doc:`dump image <dump_image>` doc page.
+:doc:`dump image <dump_image>` page.

doc/src/fix_bond_react.rst

@@ -396,6 +396,12 @@ atoms are present. The velocity of each created atom is initialized in
 a random direction with a magnitude calculated from the instantaneous
 temperature of the reaction site.
 
+.. note::
+
+   The 'Coords' section must be included in the post-reaction template
+   when creating atoms because these coordinates are used to determine
+   where new atoms are inserted.
+
 The handedness of atoms that are chiral centers can be enforced by
 listing their IDs in the ChiralIDs section. A chiral atom must be
 bonded to four atoms with mutually different atom types. This feature

doc/src/fix_npt_sphere.rst

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

doc/src/fix_nve_sphere.rst

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

doc/src/fix_nvt_sphere.rst

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

doc/src/if.rst

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

doc/src/label.rst

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

doc/src/pair_cs.rst

@@ -114,7 +114,7 @@ Description
 """""""""""
 
 These pair styles are designed to be used with the adiabatic
-core/shell model of :ref:`(Mitchell and Finchham) <MitchellFinchham2>`.  See
+core/shell model of :ref:`(Mitchell and Fincham) <MitchellFincham3>`.  See
 the :doc:`Howto coreshell <Howto_coreshell>` page for an overview of
 the model as implemented in LAMMPS.
 
@@ -196,7 +196,7 @@ none
 
 ----------
 
-.. _MitchellFinchham2:
+.. _MitchellFincham3:
 
 **(Mitchell and Finchham)** Mitchell, Finchham, J Phys Condensed Matter,
 5, 1031-1038 (1993).

doc/src/pair_dipole.rst

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

doc/src/pair_meam.rst

@@ -413,10 +413,8 @@ none
 
 .. _Gullet:
 
-**(Gullet)** Gullet, Wagner, Slepoy, SANDIA Report 2003-8782 (2003).
-This report may be accessed on-line via `this link <sandreport_>`_.
-
-.. _sandreport: http://infoserve.sandia.gov/sand_doc/2003/038782.pdf
+**(Gullet)** Gullet, Wagner, Slepoy, SANDIA Report 2003-8782 (2003). DOI:10.2172/918395
+This report may be accessed on-line via `this link <https://download.lammps.org/pdfs/MEAM_report_2003.pdf>`_.
 
 .. _Lee:
 

doc/src/pair_peri.rst

@@ -51,8 +51,9 @@ Description
 """""""""""
 
 The peridynamic pair styles implement material models that can be used
-at the mesoscopic and macroscopic scales.  See `this document <PDF/PDLammps_overview.pdf>`_ for an overview of LAMMPS commands
-for Peridynamics modeling.
+at the mesoscopic and macroscopic scales.  See `this document
+<PDF/PDLammps_overview.pdf>`_ for an overview of LAMMPS commands for
+Peridynamics modeling.
 
 Style *peri/pmb* implements the Peridynamic bond-based prototype
 microelastic brittle (PMB) model.
@@ -66,10 +67,12 @@ peridynamic viscoelastic solid (VES) model.
 Style *peri/eps* implements the Peridynamic state-based elastic-plastic
 solid (EPS) model.
 
-The canonical papers on Peridynamics are :ref:`(Silling 2000) <Silling2000>`
-and :ref:`(Silling 2007) <Silling2007>`.  The implementation of Peridynamics
-in LAMMPS is described in :ref:`(Parks) <Parks>`.  Also see the `PDLAMMPS user guide <http://www.sandia.gov/~mlparks/papers/PDLAMMPS.pdf>`_ for
-more details about its implementation.
+The canonical papers on Peridynamics are :ref:`(Silling 2000)
+<Silling2000>` and :ref:`(Silling 2007) <Silling2007>`.  The
+implementation of Peridynamics in LAMMPS is described in :ref:`(Parks)
+<Parks>`.  Also see the `PDLAMMPS user guide
+<https://download.lammps.org/pdfs/PDLAMMPS_user_guide.pdf>`_ for more
+details about its implementation.
 
 The peridynamic VES and EPS models in PDLAMMPS were implemented by
 R. Rahman and J. T. Foster at University of Texas at San Antonio.  The
@@ -165,8 +168,9 @@ shift option.
 The :doc:`pair_modify <pair_modify>` table and tail options are not
 relevant for these pair styles.
 
-These pair styles write their information to :doc:`binary restart files <restart>`, so pair_style and pair_coeff commands do not need
-to be specified in an input script that reads a restart file.
+These pair styles write their information to :doc:`binary restart files
+<restart>`, so pair_style and pair_coeff commands do not need to be
+specified in an input script that reads a restart file.
 
 These pair styles can only be used via the *pair* keyword of the
 :doc:`run_style respa <run_style>` command.  They do not support the
@@ -177,8 +181,9 @@ These pair styles can only be used via the *pair* keyword of the
 Restrictions
 """"""""""""
 
-All of these styles are part of the PERI package. They are only
-enabled if LAMMPS was built with that package.  See the :doc:`Build package <Build_package>` page for more info.
+All of these styles are part of the PERI package. They are only enabled
+if LAMMPS was built with that package.  See the :doc:`Build package
+<Build_package>` page for more info.
 
 Related commands
 """"""""""""""""

doc/src/variable.rst

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

doc/utils/requirements.txt

@@ -1,8 +1,8 @@
-Sphinx
-sphinxcontrib-spelling
+Sphinx==5.2.3
+sphinxcontrib-spelling==7.6.2
 git+https://github.com/akohlmey/sphinx-fortran@parallel-read
-sphinx_tabs
-breathe
-Pygments
-six
-pyyaml
+sphinx_tabs==3.4.1
+breathe==4.34.0
+Pygments==2.13.0
+six==1.16.0
+pyyaml==6.0

doc/utils/sphinx-config/_static/css/lammps.css

@@ -39,6 +39,18 @@ hr {
      display: none;
 }
 
+#userdoc.toctree-wrapper.compound p {
+     display: none;
+}
+
+#progdoc.toctree-wrapper.compound p {
+     display: none;
+}
+
+#reference.toctree-wrapper.compound p {
+     display: none;
+}
+
 .ui.tabular.menu .item {
     padding-right: 1em;
     padding-left: 1em;

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

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

doc/utils/sphinx-config/_themes/lammps_theme/search.html

@@ -4,17 +4,8 @@
 {% endblock %}
 {% block body %}
   <div id="search-results">
-<script>
-  (function() {
-    var cx = '012685039201307511604:um7if1hinba';
-    var gcse = document.createElement('script');
-    gcse.type = 'text/javascript';
-    gcse.async = true;
-    gcse.src = 'https://cse.google.com/cse.js?cx=' + cx;
-    var s = document.getElementsByTagName('script')[0];
-    s.parentNode.insertBefore(gcse, s);
-  })();
-</script>
-<gcse:searchresults-only></gcse:searchresults-only>
+    <script async src="https://cse.google.com/cse.js?cx=000956471495417073164:5ggs_m6ymw0">
+    </script>
+    <div class="gcse-search"></div>
   </div>
 {% endblock %}

doc/utils/sphinx-config/conf.py.in

@@ -38,7 +38,7 @@ sys.path.append(os.path.join(LAMMPS_DOC_DIR, 'src', '_ext'))
 # -- General configuration ------------------------------------------------
 
 # If your documentation needs a minimal Sphinx version, state it here.
-#needs_sphinx = '1.0'
+needs_sphinx = '5.2.0'
 
 # Add any Sphinx extension module names here, as strings. They can be
 # extensions coming with Sphinx (named 'sphinx.ext.*') or your custom
@@ -53,10 +53,6 @@ extensions = [
     'tab_or_note',
     'breathe',
 ]
-# 2017-12-07: commented out, since this package is broken with Sphinx 16.x
-#             yet we can no longer use Sphinx 15.x, since that breaks with
-#             newer version of the multiprocessor module.
-#    'sphinxcontrib.images',
 
 images_config = {
   'default_image_width' : '25%',

examples/COUPLE/lammps_spparks/README

@@ -31,7 +31,7 @@ spkpath.h           contains path to SPPARKS home directory
 After editing the Makefile, lmppath.h, and spkpath.h to make them
 suitable for your box, type:
 
-g++ -f Makefile.g++
+make -f Makefile.g++
 
 and you should get the lmpspk executable.
 

examples/PACKAGES/pace/plugin/CMakeLists.txt

@@ -31,6 +31,23 @@ elseif(CMAKE_SYSTEM_NAME STREQUAL "Windows")
   if(CMAKE_CROSSCOMPILING)
     set_target_properties(paceplugin  PROPERTIES LINK_FLAGS "-Wl,--export-all-symbols")
   endif()
+
+  get_lammps_version(${LAMMPS_SOURCE_DIR}/version.h LAMMPS_VERSION)
+  find_program(MAKENSIS_PATH makensis)
+  if(MAKENSIS_PATH)
+    execute_process(COMMAND ${CMAKE_COMMAND} -E copy_if_different ${CMAKE_SOURCE_DIR}/lammps.ico
+      ${CMAKE_SOURCE_DIR}/lammps-text-logo-wide.bmp ${CMAKE_SOURCE_DIR}/paceplugin.nsis ${CMAKE_BINARY_DIR})
+    if(BUILD_MPI)
+      add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION}-MPI paceplugin.nsis
+        DEPENDS paceplugin
+        BYPRODUCTS LAMMPS-ML-PACE-plugin-${LAMMPS_VERSION}-MPI.exe)
+    else()
+      add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION} paceplugin.nsis
+        COMMAND ${CMAKE_COMMAND} -E echo ${PWD}
+        DEPENDS paceplugin lammps.ico lammps-text-logo-wide.bmp paceplugin.nsis
+        BYPRODUCTS LAMMPS-ML-PACE-plugin-${LAMMPS_VERSION}.exe)
+    endif()
+  endif()
 else()
   set_target_properties(paceplugin PROPERTIES LINK_FLAGS "-rdynamic")
 endif()

examples/plugins/LAMMPSInterfaceCXX.cmake

@@ -1,6 +1,9 @@
 # Cmake script code to define the LAMMPS C++ interface
 # settings required for building LAMMPS plugins
 
+set(LAMMPS_THIRDPARTY_URL "https://download.lammps.org/thirdparty"
+  CACHE STRING "URL for thirdparty package downloads")
+
 ################################################################################
 # helper function
 function(validate_option name values)
@@ -37,10 +40,56 @@ else()
 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)
+  # 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()

lib/atc/FE_Element.cpp

@@ -186,7 +186,7 @@ static const double localCoordinatesTolerance = 1.e-09;
       double c[3]={0,0,0};
       c[0] = y0*X[0] - y0*X[1] - y0*X[2] + y0*X[3] - x0*Y[0] + (X[1]*Y[0])*0.5 + (X[2]*Y[0])*0.5 + x0*Y[1] - (X[0]*Y[1])*0.5 - (X[3]*Y[1])*0.5 + x0*Y[2] - (X[0]*Y[2])*0.5 - (X[3]*Y[2])*0.5 - x0*Y[3] + (X[1]*Y[3])*0.5 + (X[2]*Y[3])*0.5;
       c[1] = -(y0*X[0]) + y0*X[1] - y0*X[2] + y0*X[3] + x0*Y[0] - X[1]*Y[0] - x0*Y[1] + X[0]*Y[1] + x0*Y[2] - X[3]*Y[2] - x0*Y[3] + X[2]*Y[3];
-      c[1] = (X[1]*Y[0])*0.5 - (X[2]*Y[0])*0.5 - (X[0]*Y[1])*0.5 + (X[3]*Y[1])*0.5 + (X[0]*Y[2])*0.5 - (X[3]*Y[2])*0.5 - (X[1]*Y[3])*0.5 + (X[2]*Y[3])*0.5;
+      c[2] = (X[1]*Y[0])*0.5 - (X[2]*Y[0])*0.5 - (X[0]*Y[1])*0.5 + (X[3]*Y[1])*0.5 + (X[0]*Y[2])*0.5 - (X[3]*Y[2])*0.5 - (X[1]*Y[3])*0.5 + (X[2]*Y[3])*0.5;
       double xi2[2]={0,0};
       int nroots = solve_quadratic(c,xi2);
       if (nroots == 0) throw ATC_Error("no real roots in 2D analytic projection guess");

lib/kokkos/CHANGELOG.md

@@ -1,5 +1,21 @@
 # Change Log
 
+## [3.6.01](https://github.com/kokkos/kokkos/tree/3.6.01) (2022-05-23)
+[Full Changelog](https://github.com/kokkos/kokkos/compare/3.6.00...3.6.01)
+
+### Bug Fixes:
+- Fix Threads: Fix serial resizing scratch space (3.6.01 cherry-pick) [\#5109](https://github.com/kokkos/kokkos/pull/5109)
+- Fix ScatterMin/ScatterMax to use proper atomics (3.6.01 cherry-pick) [\#5046](https://github.com/kokkos/kokkos/pull/5046)
+- Fix allocating large Views [\#4907](https://github.com/kokkos/kokkos/pull/4907)
+- Fix bounds errors with Kokkos::sort [\#4980](https://github.com/kokkos/kokkos/pull/4980)
+- Fix HIP version when printing the configuration [\#4872](https://github.com/kokkos/kokkos/pull/4872)
+- Fixed `_CUDA_ARCH__` to `__CUDA_ARCH__` for CUDA LDG [\#4893](https://github.com/kokkos/kokkos/pull/4893)
+- Fixed an incorrect struct initialization [\#5028](https://github.com/kokkos/kokkos/pull/5028)
+- Fix racing condition in `HIPParallelLaunch` [\#5008](https://github.com/kokkos/kokkos/pull/5008)
+- Avoid deprecation warnings with `OpenMPExec::validate_partition` [\#4982](https://github.com/kokkos/kokkos/pull/4982)
+- Make View self-assignment not produce double-free [\#5024](https://github.com/kokkos/kokkos/pull/5024)
+
+
 ## [3.6.00](https://github.com/kokkos/kokkos/tree/3.6.00) (2022-02-18)
 [Full Changelog](https://github.com/kokkos/kokkos/compare/3.5.00...3.6.00)
 

lib/kokkos/CMakeLists.txt

@@ -136,7 +136,7 @@ ENDIF()
 
 set(Kokkos_VERSION_MAJOR 3)
 set(Kokkos_VERSION_MINOR 6)
-set(Kokkos_VERSION_PATCH 00)
+set(Kokkos_VERSION_PATCH 01)
 set(Kokkos_VERSION "${Kokkos_VERSION_MAJOR}.${Kokkos_VERSION_MINOR}.${Kokkos_VERSION_PATCH}")
 math(EXPR KOKKOS_VERSION "${Kokkos_VERSION_MAJOR} * 10000 + ${Kokkos_VERSION_MINOR} * 100 + ${Kokkos_VERSION_PATCH}")
 

lib/kokkos/Makefile.kokkos

@@ -12,7 +12,7 @@ endif
 
 KOKKOS_VERSION_MAJOR = 3
 KOKKOS_VERSION_MINOR = 6
-KOKKOS_VERSION_PATCH = 00
+KOKKOS_VERSION_PATCH = 01
 KOKKOS_VERSION = $(shell echo $(KOKKOS_VERSION_MAJOR)*10000+$(KOKKOS_VERSION_MINOR)*100+$(KOKKOS_VERSION_PATCH) | bc)
 
 # Options: Cuda,HIP,SYCL,OpenMPTarget,OpenMP,Threads,Serial

lib/kokkos/algorithms/src/Kokkos_Sort.hpp

@@ -422,54 +422,34 @@ class BinSort {
 
 template <class KeyViewType>
 struct BinOp1D {
-  int max_bins_;
-  double mul_;
-  typename KeyViewType::const_value_type range_;
-  typename KeyViewType::const_value_type min_;
+  int max_bins_ = {};
+  double mul_   = {};
+  double min_   = {};
 
-  BinOp1D()
-      : max_bins_(0),
-        mul_(0.0),
-        range_(typename KeyViewType::const_value_type()),
-        min_(typename KeyViewType::const_value_type()) {}
+  BinOp1D() = default;
 
   // Construct BinOp with number of bins, minimum value and maxuimum value
   BinOp1D(int max_bins__, typename KeyViewType::const_value_type min,
           typename KeyViewType::const_value_type max)
       : max_bins_(max_bins__ + 1),
-        // Cast to int64_t to avoid possible overflow when using integer
-        mul_(std::is_integral<typename KeyViewType::const_value_type>::value
-                 ? 1.0 * max_bins__ / (int64_t(max) - int64_t(min))
-                 : 1.0 * max_bins__ / (max - min)),
-        range_(max - min),
-        min_(min) {
+        // Cast to double to avoid possible overflow when using integer
+        mul_(static_cast<double>(max_bins__) /
+             (static_cast<double>(max) - static_cast<double>(min))),
+        min_(static_cast<double>(min)) {
     // For integral types the number of bins may be larger than the range
     // in which case we can exactly have one unique value per bin
     // and then don't need to sort bins.
     if (std::is_integral<typename KeyViewType::const_value_type>::value &&
-        static_cast<uint64_t>(range_) <= static_cast<uint64_t>(max_bins__)) {
+        (static_cast<double>(max) - static_cast<double>(min)) <=
+            static_cast<double>(max_bins__)) {
       mul_ = 1.;
     }
   }
 
   // Determine bin index from key value
-  template <
-      class ViewType,
-      std::enable_if_t<!std::is_integral<typename ViewType::value_type>::value,
-                       bool> = true>
+  template <class ViewType>
   KOKKOS_INLINE_FUNCTION int bin(ViewType& keys, const int& i) const {
-    return int(mul_ * (keys(i) - min_));
-  }
-
-  // Determine bin index from key value
-  template <
-      class ViewType,
-      std::enable_if_t<std::is_integral<typename ViewType::value_type>::value,
-                       bool> = true>
-  KOKKOS_INLINE_FUNCTION int bin(ViewType& keys, const int& i) const {
-    // The cast to int64_t is necessary because otherwise HIP returns the wrong
-    // result.
-    return int(mul_ * (int64_t(keys(i)) - int64_t(min_)));
+    return static_cast<int>(mul_ * (static_cast<double>(keys(i)) - min_));
   }
 
   // Return maximum bin index + 1
@@ -486,10 +466,9 @@ struct BinOp1D {
 
 template <class KeyViewType>
 struct BinOp3D {
-  int max_bins_[3];
-  double mul_[3];
-  typename KeyViewType::non_const_value_type range_[3];
-  typename KeyViewType::non_const_value_type min_[3];
+  int max_bins_[3] = {};
+  double mul_[3]   = {};
+  double min_[3]   = {};
 
   BinOp3D() = default;
 
@@ -498,15 +477,15 @@ struct BinOp3D {
     max_bins_[0] = max_bins__[0];
     max_bins_[1] = max_bins__[1];
     max_bins_[2] = max_bins__[2];
-    mul_[0]      = 1.0 * max_bins__[0] / (max[0] - min[0]);
-    mul_[1]      = 1.0 * max_bins__[1] / (max[1] - min[1]);
-    mul_[2]      = 1.0 * max_bins__[2] / (max[2] - min[2]);
-    range_[0]    = max[0] - min[0];
-    range_[1]    = max[1] - min[1];
-    range_[2]    = max[2] - min[2];
-    min_[0]      = min[0];
-    min_[1]      = min[1];
-    min_[2]      = min[2];
+    mul_[0]      = static_cast<double>(max_bins__[0]) /
+              (static_cast<double>(max[0]) - static_cast<double>(min[0]));
+    mul_[1] = static_cast<double>(max_bins__[1]) /
+              (static_cast<double>(max[1]) - static_cast<double>(min[1]));
+    mul_[2] = static_cast<double>(max_bins__[2]) /
+              (static_cast<double>(max[2]) - static_cast<double>(min[2]));
+    min_[0] = static_cast<double>(min[0]);
+    min_[1] = static_cast<double>(min[1]);
+    min_[2] = static_cast<double>(min[2]);
   }
 
   template <class ViewType>
@@ -596,9 +575,9 @@ std::enable_if_t<Kokkos::is_execution_space<ExecutionSpace>::value> sort(
   // TODO: figure out better max_bins then this ...
   int64_t max_bins = view.extent(0) / 2;
   if (std::is_integral<typename ViewType::non_const_value_type>::value) {
-    // Cast to int64_t to avoid possible overflow when using integer
-    int64_t const max_val = result.max_val;
-    int64_t const min_val = result.min_val;
+    // Cast to double to avoid possible overflow when using integer
+    auto const max_val = static_cast<double>(result.max_val);
+    auto const min_val = static_cast<double>(result.min_val);
     // using 10M as the cutoff for special behavior (roughly 40MB for the count
     // array)
     if ((max_val - min_val) < 10000000) {
@@ -606,6 +585,10 @@ std::enable_if_t<Kokkos::is_execution_space<ExecutionSpace>::value> sort(
       sort_in_bins = false;
     }
   }
+  if (std::is_floating_point<typename ViewType::non_const_value_type>::value) {
+    KOKKOS_ASSERT(std::isfinite(static_cast<double>(result.max_val) -
+                                static_cast<double>(result.min_val)));
+  }
 
   BinSort<ViewType, CompType> bin_sort(
       view, CompType(max_bins, result.min_val, result.max_val), sort_in_bins);

lib/kokkos/algorithms/unit_tests/TestSort.hpp

@@ -353,6 +353,55 @@ void test_issue_1160_impl() {
   }
 }
 
+template <class ExecutionSpace>
+void test_issue_4978_impl() {
+  Kokkos::View<long long*, ExecutionSpace> element_("element", 9);
+
+  auto h_element = Kokkos::create_mirror_view(element_);
+
+  h_element(0) = LLONG_MIN;
+  h_element(1) = 0;
+  h_element(2) = 3;
+  h_element(3) = 2;
+  h_element(4) = 1;
+  h_element(5) = 3;
+  h_element(6) = 6;
+  h_element(7) = 4;
+  h_element(8) = 3;
+
+  ExecutionSpace exec;
+  Kokkos::deep_copy(exec, element_, h_element);
+
+  Kokkos::sort(exec, element_);
+
+  Kokkos::deep_copy(exec, h_element, element_);
+  exec.fence();
+
+  ASSERT_EQ(h_element(0), LLONG_MIN);
+  ASSERT_EQ(h_element(1), 0);
+  ASSERT_EQ(h_element(2), 1);
+  ASSERT_EQ(h_element(3), 2);
+  ASSERT_EQ(h_element(4), 3);
+  ASSERT_EQ(h_element(5), 3);
+  ASSERT_EQ(h_element(6), 3);
+  ASSERT_EQ(h_element(7), 4);
+  ASSERT_EQ(h_element(8), 6);
+}
+
+template <class ExecutionSpace, class T>
+void test_sort_integer_overflow() {
+  // array with two extrema in reverse order to expose integer overflow bug in
+  // bin calculation
+  T a[2]  = {Kokkos::Experimental::finite_max<T>::value,
+            Kokkos::Experimental::finite_min<T>::value};
+  auto vd = Kokkos::create_mirror_view_and_copy(
+      ExecutionSpace(), Kokkos::View<T[2], Kokkos::HostSpace>(a));
+  Kokkos::sort(vd, /*force using Kokkos bin sort*/ true);
+  auto vh = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), vd);
+  EXPECT_TRUE(std::is_sorted(vh.data(), vh.data() + 2))
+      << "view (" << vh[0] << ", " << vh[1] << ") is not sorted";
+}
+
 //----------------------------------------------------------------------------
 
 template <class ExecutionSpace, typename KeyType>
@@ -376,6 +425,11 @@ void test_issue_1160_sort() {
   test_issue_1160_impl<ExecutionSpace>();
 }
 
+template <class ExecutionSpace>
+void test_issue_4978_sort() {
+  test_issue_4978_impl<ExecutionSpace>();
+}
+
 template <class ExecutionSpace, typename KeyType>
 void test_sort(unsigned int N) {
   test_1D_sort<ExecutionSpace, KeyType>(N);
@@ -385,6 +439,10 @@ void test_sort(unsigned int N) {
   test_dynamic_view_sort<ExecutionSpace, KeyType>(N);
 #endif
   test_issue_1160_sort<ExecutionSpace>();
+  test_issue_4978_sort<ExecutionSpace>();
+  test_sort_integer_overflow<ExecutionSpace, long long>();
+  test_sort_integer_overflow<ExecutionSpace, unsigned long long>();
+  test_sort_integer_overflow<ExecutionSpace, int>();
 }
 }  // namespace Impl
 }  // namespace Test

lib/kokkos/bin/nvcc_wrapper

@@ -306,13 +306,13 @@ do
   -A)
     ;;
 
-  #strip of -std=c++98 due to nvcc warnings and Tribits will place both -std=c++11 and -std=c++98
+  #strip off -std=c++98 due to nvcc warnings and Tribits will place both -std=c++11 and -std=c++98
   -std=c++98|--std=c++98)
     ;;
-  #strip of pedantic because it produces endless warnings about #LINE added by the preprocessor
-  -pedantic|-Wpedantic|-ansi)
+  #strip off pedantic because it produces endless warnings about #LINE added by the preprocessor
+  -pedantic|-pedantic-errors|-Wpedantic|-ansi)
     ;;
-  #strip of -Woverloaded-virtual to avoid "cc1: warning: command line option ‘-Woverloaded-virtual’ is valid for C++/ObjC++ but not for C"
+  #strip off -Woverloaded-virtual to avoid "cc1: warning: command line option ‘-Woverloaded-virtual’ is valid for C++/ObjC++ but not for C"
   -Woverloaded-virtual)
     ;;
   #strip -Xcompiler because we add it

lib/kokkos/containers/src/Kokkos_ScatterView.hpp

@@ -369,18 +369,6 @@ struct ScatterValue<ValueType, Kokkos::Experimental::ScatterProd, DeviceType,
     Kokkos::atomic_div(&value, rhs);
   }
 
-  KOKKOS_FORCEINLINE_FUNCTION
-  void atomic_prod(ValueType& dest, const ValueType& src) const {
-    bool success = false;
-    while (!success) {
-      ValueType dest_old = dest;
-      ValueType dest_new = dest_old * src;
-      dest_new =
-          Kokkos::atomic_compare_exchange<ValueType>(&dest, dest_old, dest_new);
-      success = ((dest_new - dest_old) / dest_old <= 1e-15);
-    }
-  }
-
   KOKKOS_INLINE_FUNCTION
   void join(ValueType& dest, const ValueType& src) const {
     atomic_prod(&dest, src);
@@ -440,21 +428,9 @@ struct ScatterValue<ValueType, Kokkos::Experimental::ScatterMin, DeviceType,
   KOKKOS_FORCEINLINE_FUNCTION ScatterValue(ScatterValue&& other)
       : value(other.value) {}
 
-  KOKKOS_FORCEINLINE_FUNCTION
-  void atomic_min(ValueType& dest, const ValueType& src) const {
-    bool success = false;
-    while (!success) {
-      ValueType dest_old = dest;
-      ValueType dest_new = (dest_old > src) ? src : dest_old;
-      dest_new =
-          Kokkos::atomic_compare_exchange<ValueType>(&dest, dest_old, dest_new);
-      success = ((dest_new - dest_old) / dest_old <= 1e-15);
-    }
-  }
-
   KOKKOS_INLINE_FUNCTION
   void join(ValueType& dest, const ValueType& src) const {
-    atomic_min(dest, src);
+    atomic_min(&dest, src);
   }
 
   KOKKOS_INLINE_FUNCTION
@@ -511,21 +487,9 @@ struct ScatterValue<ValueType, Kokkos::Experimental::ScatterMax, DeviceType,
   KOKKOS_FORCEINLINE_FUNCTION ScatterValue(ScatterValue&& other)
       : value(other.value) {}
 
-  KOKKOS_FORCEINLINE_FUNCTION
-  void atomic_max(ValueType& dest, const ValueType& src) const {
-    bool success = false;
-    while (!success) {
-      ValueType dest_old = dest;
-      ValueType dest_new = (dest_old < src) ? src : dest_old;
-      dest_new =
-          Kokkos::atomic_compare_exchange<ValueType>(&dest, dest_old, dest_new);
-      success = ((dest_new - dest_old) / dest_old <= 1e-15);
-    }
-  }
-
   KOKKOS_INLINE_FUNCTION
   void join(ValueType& dest, const ValueType& src) const {
-    atomic_max(dest, src);
+    atomic_max(&dest, src);
   }
 
   KOKKOS_INLINE_FUNCTION

lib/kokkos/containers/src/Kokkos_Vector.hpp

@@ -162,7 +162,7 @@ class vector : public DualView<Scalar*, LayoutLeft, Arg1Type> {
     }
     DV::sync_host();
     DV::modify_host();
-    if (it < begin() || it > end())
+    if (std::less<>()(it, begin()) || std::less<>()(end(), it))
       Kokkos::abort("Kokkos::vector::insert : invalid insert iterator");
     if (count == 0) return it;
     ptrdiff_t start = std::distance(begin(), it);
@@ -189,27 +189,21 @@ class vector : public DualView<Scalar*, LayoutLeft, Arg1Type> {
                           iterator>::type
   insert(iterator it, InputIterator b, InputIterator e) {
     ptrdiff_t count = std::distance(b, e);
-    if (count == 0) return it;
 
     DV::sync_host();
     DV::modify_host();
-    if (it < begin() || it > end())
+    if (std::less<>()(it, begin()) || std::less<>()(end(), it))
       Kokkos::abort("Kokkos::vector::insert : invalid insert iterator");
 
-    bool resized = false;
-    if ((size() == 0) && (it == begin())) {
-      resize(count);
-      it      = begin();
-      resized = true;
-    }
     ptrdiff_t start = std::distance(begin(), it);
     auto org_size   = size();
-    if (!resized) resize(size() + count);
-    it = begin() + start;
+
+    // Note: resize(...) invalidates it; use begin() + start instead
+    resize(size() + count);
 
     std::copy_backward(begin() + start, begin() + org_size,
                        begin() + org_size + count);
-    std::copy(b, e, it);
+    std::copy(b, e, begin() + start);
 
     return begin() + start;
   }

lib/kokkos/containers/unit_tests/TestVector.hpp

@@ -172,6 +172,23 @@ struct test_vector_insert {
       run_test(a);
       check_test(a, size);
     }
+    { test_vector_insert_into_empty(size); }
+  }
+
+  void test_vector_insert_into_empty(const size_t size) {
+    using Vector = Kokkos::vector<Scalar, Device>;
+    {
+      Vector a;
+      Vector b(size);
+      a.insert(a.begin(), b.begin(), b.end());
+      ASSERT_EQ(a.size(), size);
+    }
+
+    {
+      Vector c;
+      c.insert(c.begin(), size, Scalar{});
+      ASSERT_EQ(c.size(), size);
+    }
   }
 };
 

lib/kokkos/core/src/CMakeLists.txt

@@ -8,6 +8,7 @@ KOKKOS_INCLUDE_DIRECTORIES(
 INSTALL (DIRECTORY
   "${CMAKE_CURRENT_SOURCE_DIR}/"
   DESTINATION ${KOKKOS_HEADER_DIR}
+  FILES_MATCHING
   PATTERN desul/src EXCLUDE
   PATTERN "*.inc"
   PATTERN "*.inc_*"

lib/kokkos/core/src/Cuda/Kokkos_Cuda_Instance.cpp

@@ -1007,6 +1007,15 @@ void CudaSpaceInitializer::print_configuration(std::ostream &msg,
 }
 
 }  // namespace Impl
+
+#ifdef KOKKOS_ENABLE_CXX14
+namespace Tools {
+namespace Experimental {
+constexpr DeviceType DeviceTypeTraits<Cuda>::id;
+}
+}  // namespace Tools
+#endif
+
 }  // namespace Kokkos
 
 #else

lib/kokkos/core/src/Cuda/Kokkos_Cuda_View.hpp

@@ -139,7 +139,7 @@ struct CudaLDGFetch {
 
   template <typename iType>
   KOKKOS_INLINE_FUNCTION ValueType operator[](const iType& i) const {
-#if defined(__CUDA_ARCH__) && (350 <= _CUDA_ARCH__)
+#if defined(__CUDA_ARCH__) && (350 <= __CUDA_ARCH__)
     AliasType v = __ldg(reinterpret_cast<const AliasType*>(&m_ptr[i]));
     return *(reinterpret_cast<ValueType*>(&v));
 #else

lib/kokkos/core/src/HIP/Kokkos_HIP_Instance.cpp

@@ -132,7 +132,8 @@ void HIPInternal::print_configuration(std::ostream &s) const {
   s << "macro  KOKKOS_ENABLE_HIP : defined" << '\n';
 #if defined(HIP_VERSION)
   s << "macro  HIP_VERSION = " << HIP_VERSION << " = version "
-    << HIP_VERSION / 100 << "." << HIP_VERSION % 100 << '\n';
+    << HIP_VERSION_MAJOR << '.' << HIP_VERSION_MINOR << '.' << HIP_VERSION_PATCH
+    << '\n';
 #endif
 
   for (int i = 0; i < dev_info.m_hipDevCount; ++i) {
@@ -467,7 +468,6 @@ void HIPInternal::finalize() {
 }
 
 char *HIPInternal::get_next_driver(size_t driverTypeSize) const {
-  std::lock_guard<std::mutex> const lock(m_mutexWorkArray);
   if (d_driverWorkArray == nullptr) {
     KOKKOS_IMPL_HIP_SAFE_CALL(
         hipHostMalloc(&d_driverWorkArray,

lib/kokkos/core/src/HIP/Kokkos_HIP_KernelLaunch.hpp

@@ -490,6 +490,8 @@ struct HIPParallelLaunch<
 
       KOKKOS_ENSURE_HIP_LOCK_ARRAYS_ON_DEVICE();
 
+      std::lock_guard<std::mutex> const lock(hip_instance->m_mutexWorkArray);
+
       // Invoke the driver function on the device
       DriverType *d_driver = reinterpret_cast<DriverType *>(
           hip_instance->get_next_driver(sizeof(DriverType)));

lib/kokkos/core/src/HIP/Kokkos_HIP_Locks.cpp

@@ -56,8 +56,7 @@ namespace Kokkos {
 
 #ifdef KOKKOS_ENABLE_HIP_RELOCATABLE_DEVICE_CODE
 namespace Impl {
-__device__ __constant__ HIPLockArrays g_device_hip_lock_arrays = {nullptr,
-                                                                  nullptr, 0};
+__device__ __constant__ HIPLockArrays g_device_hip_lock_arrays = {nullptr, 0};
 }
 #endif
 

lib/kokkos/core/src/HIP/Kokkos_HIP_Space.cpp

@@ -464,6 +464,15 @@ void HIPSpaceInitializer::print_configuration(std::ostream& msg,
 }
 
 }  // namespace Impl
+
+#ifdef KOKKOS_ENABLE_CXX14
+namespace Tools {
+namespace Experimental {
+constexpr DeviceType DeviceTypeTraits<Kokkos::Experimental::HIP>::id;
+}
+}  // namespace Tools
+#endif
+
 }  // namespace Kokkos
 
 //==============================================================================

lib/kokkos/core/src/HPX/Kokkos_HPX.cpp

@@ -199,6 +199,15 @@ void HPXSpaceInitializer::print_configuration(std::ostream &msg,
 }
 
 }  // namespace Impl
+
+#ifdef KOKKOS_ENABLE_CXX14
+namespace Tools {
+namespace Experimental {
+constexpr DeviceType DeviceTypeTraits<Kokkos::Experimental::HPX>::id;
+}
+}  // namespace Tools
+#endif
+
 }  // namespace Kokkos
 
 #else

lib/kokkos/core/src/Kokkos_Cuda.hpp

@@ -260,6 +260,7 @@ template <>
 struct DeviceTypeTraits<Cuda> {
   /// \brief An ID to differentiate (for example) Serial from OpenMP in Tooling
   static constexpr DeviceType id = DeviceType::Cuda;
+  static int device_id(const Cuda& exec) { return exec.cuda_device(); }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/Kokkos_HIP_Space.hpp

@@ -571,6 +571,9 @@ namespace Experimental {
 template <>
 struct DeviceTypeTraits<Kokkos::Experimental::HIP> {
   static constexpr DeviceType id = DeviceType::HIP;
+  static int device_id(const Kokkos::Experimental::HIP& exec) {
+    return exec.hip_device();
+  }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/Kokkos_HPX.hpp

@@ -500,6 +500,7 @@ namespace Experimental {
 template <>
 struct DeviceTypeTraits<Kokkos::Experimental::HPX> {
   static constexpr DeviceType id = DeviceType::HPX;
+  static int device_id(const Kokkos::Experimental::HPX &) { return 0; }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/Kokkos_OpenMP.hpp

@@ -179,6 +179,7 @@ namespace Experimental {
 template <>
 struct DeviceTypeTraits<OpenMP> {
   static constexpr DeviceType id = DeviceType::OpenMP;
+  static int device_id(const OpenMP&) { return 0; }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/Kokkos_OpenMPTarget.hpp

@@ -130,6 +130,9 @@ template <>
 struct DeviceTypeTraits<::Kokkos::Experimental::OpenMPTarget> {
   static constexpr DeviceType id =
       ::Kokkos::Profiling::Experimental::DeviceType::OpenMPTarget;
+  static int device_id(const Kokkos::Experimental::OpenMPTarget&) {
+    return omp_get_default_device();
+  }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/Kokkos_SYCL.hpp

@@ -182,6 +182,9 @@ template <>
 struct DeviceTypeTraits<Kokkos::Experimental::SYCL> {
   /// \brief An ID to differentiate (for example) Serial from OpenMP in Tooling
   static constexpr DeviceType id = DeviceType::SYCL;
+  static int device_id(const Kokkos::Experimental::SYCL& exec) {
+    return exec.sycl_device();
+  }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/Kokkos_Serial.hpp

@@ -226,6 +226,7 @@ namespace Experimental {
 template <>
 struct DeviceTypeTraits<Serial> {
   static constexpr DeviceType id = DeviceType::Serial;
+  static int device_id(const Serial&) { return 0; }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/Kokkos_Threads.hpp

@@ -175,6 +175,7 @@ namespace Experimental {
 template <>
 struct DeviceTypeTraits<Threads> {
   static constexpr DeviceType id = DeviceType::Threads;
+  static int device_id(const Threads&) { return 0; }
 };
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/OpenMP/Kokkos_OpenMP_Exec.cpp

@@ -67,8 +67,9 @@ __thread int t_openmp_hardware_id            = 0;
 __thread Impl::OpenMPExec *t_openmp_instance = nullptr;
 
 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE_3
-void OpenMPExec::validate_partition(const int nthreads, int &num_partitions,
-                                    int &partition_size) {
+void OpenMPExec::validate_partition_impl(const int nthreads,
+                                         int &num_partitions,
+                                         int &partition_size) {
   if (nthreads == 1) {
     num_partitions = 1;
     partition_size = 1;
@@ -506,6 +507,15 @@ void OpenMPSpaceInitializer::print_configuration(std::ostream &msg,
 }
 
 }  // namespace Impl
+
+#ifdef KOKKOS_ENABLE_CXX14
+namespace Tools {
+namespace Experimental {
+constexpr DeviceType DeviceTypeTraits<OpenMP>::id;
+}
+}  // namespace Tools
+#endif
+
 }  // namespace Kokkos
 
 #else

lib/kokkos/core/src/OpenMP/Kokkos_OpenMP_Exec.hpp

@@ -93,7 +93,11 @@ class OpenMPExec {
 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE_3
   KOKKOS_DEPRECATED static void validate_partition(const int nthreads,
                                                    int& num_partitions,
-                                                   int& partition_size);
+                                                   int& partition_size) {
+    validate_partition_impl(nthreads, num_partitions, partition_size);
+  }
+  static void validate_partition_impl(const int nthreads, int& num_partitions,
+                                      int& partition_size);
 #endif
 
  private:
@@ -179,8 +183,8 @@ KOKKOS_DEPRECATED void OpenMP::partition_master(F const& f, int num_partitions,
 
     Exec* prev_instance = Impl::t_openmp_instance;
 
-    Exec::validate_partition(prev_instance->m_pool_size, num_partitions,
-                             partition_size);
+    Exec::validate_partition_impl(prev_instance->m_pool_size, num_partitions,
+                                  partition_size);
 
     OpenMP::memory_space space;
 

lib/kokkos/core/src/SYCL/Kokkos_SYCL_Instance.hpp

@@ -72,7 +72,7 @@ class SYCLInternal {
                                   bool force_shrink = false);
 
   uint32_t impl_get_instance_id() const;
-  int m_syclDev = -1;
+  int m_syclDev = 0;
 
   size_t m_maxWorkgroupSize   = 0;
   uint32_t m_maxConcurrency   = 0;

lib/kokkos/core/src/Threads/Kokkos_ThreadsExec.cpp

@@ -399,27 +399,68 @@ bool ThreadsExec::wake() {
 
 //----------------------------------------------------------------------------
 
+void ThreadsExec::execute_resize_scratch_in_serial() {
+  const unsigned begin = s_threads_process.m_pool_base ? 1 : 0;
+
+  auto deallocate_scratch_memory = [](ThreadsExec &exec) {
+    if (exec.m_scratch) {
+      using Record =
+          Kokkos::Impl::SharedAllocationRecord<Kokkos::HostSpace, void>;
+      Record *const r = Record::get_record(exec.m_scratch);
+      exec.m_scratch  = nullptr;
+      Record::decrement(r);
+    }
+  };
+  if (s_threads_process.m_pool_base) {
+    for (unsigned i = s_thread_pool_size[0]; begin < i;) {
+      deallocate_scratch_memory(*s_threads_exec[--i]);
+    }
+  }
+
+  s_current_function     = &first_touch_allocate_thread_private_scratch;
+  s_current_function_arg = &s_threads_process;
+
+  // Make sure function and arguments are written before activating threads.
+  memory_fence();
+
+  for (unsigned i = s_thread_pool_size[0]; begin < i;) {
+    ThreadsExec &th = *s_threads_exec[--i];
+
+    th.m_pool_state = ThreadsExec::Active;
+
+    wait_yield(th.m_pool_state, ThreadsExec::Active);
+  }
+
+  if (s_threads_process.m_pool_base) {
+    deallocate_scratch_memory(s_threads_process);
+    s_threads_process.m_pool_state = ThreadsExec::Active;
+    first_touch_allocate_thread_private_scratch(s_threads_process, nullptr);
+    s_threads_process.m_pool_state = ThreadsExec::Inactive;
+  }
+
+  s_current_function_arg = nullptr;
+  s_current_function     = nullptr;
+
+  // Make sure function and arguments are cleared before proceeding.
+  memory_fence();
+}
+
+//----------------------------------------------------------------------------
+
 void *ThreadsExec::root_reduce_scratch() {
   return s_threads_process.reduce_memory();
 }
 
-void ThreadsExec::execute_resize_scratch(ThreadsExec &exec, const void *) {
-  using Record = Kokkos::Impl::SharedAllocationRecord<Kokkos::HostSpace, void>;
-
-  if (exec.m_scratch) {
-    Record *const r = Record::get_record(exec.m_scratch);
-
-    exec.m_scratch = nullptr;
-
-    Record::decrement(r);
-  }
-
+void ThreadsExec::first_touch_allocate_thread_private_scratch(ThreadsExec &exec,
+                                                              const void *) {
   exec.m_scratch_reduce_end = s_threads_process.m_scratch_reduce_end;
   exec.m_scratch_thread_end = s_threads_process.m_scratch_thread_end;
 
   if (s_threads_process.m_scratch_thread_end) {
     // Allocate tracked memory:
     {
+      using Record =
+          Kokkos::Impl::SharedAllocationRecord<Kokkos::HostSpace, void>;
       Record *const r =
           Record::allocate(Kokkos::HostSpace(), "Kokkos::thread_scratch",
                            s_threads_process.m_scratch_thread_end);
@@ -461,7 +502,7 @@ void *ThreadsExec::resize_scratch(size_t reduce_size, size_t thread_size) {
     s_threads_process.m_scratch_reduce_end = reduce_size;
     s_threads_process.m_scratch_thread_end = reduce_size + thread_size;
 
-    execute_resize_scratch(s_threads_process, nullptr);
+    execute_resize_scratch_in_serial();
 
     s_threads_process.m_scratch = s_threads_exec[0]->m_scratch;
   }
@@ -845,6 +886,15 @@ void ThreadsSpaceInitializer::print_configuration(std::ostream &msg,
 }
 
 }  // namespace Impl
+
+#ifdef KOKKOS_ENABLE_CXX14
+namespace Tools {
+namespace Experimental {
+constexpr DeviceType DeviceTypeTraits<Threads>::id;
+}
+}  // namespace Tools
+#endif
+
 } /* namespace Kokkos */
 //----------------------------------------------------------------------------
 //----------------------------------------------------------------------------

lib/kokkos/core/src/Threads/Kokkos_ThreadsExec.hpp

@@ -123,12 +123,15 @@ class ThreadsExec {
   static void global_unlock();
   static void spawn();
 
-  static void execute_resize_scratch(ThreadsExec &, const void *);
+  static void first_touch_allocate_thread_private_scratch(ThreadsExec &,
+                                                          const void *);
   static void execute_sleep(ThreadsExec &, const void *);
 
   ThreadsExec(const ThreadsExec &);
   ThreadsExec &operator=(const ThreadsExec &);
 
+  static void execute_resize_scratch_in_serial();
+
  public:
   KOKKOS_INLINE_FUNCTION int pool_size() const { return m_pool_size; }
   KOKKOS_INLINE_FUNCTION int pool_rank() const { return m_pool_rank; }

lib/kokkos/core/src/impl/Kokkos_Profiling_Interface.hpp

@@ -118,11 +118,14 @@ template <typename ExecutionSpace>
 constexpr uint32_t device_id_root() {
   constexpr auto device_id =
       static_cast<uint32_t>(DeviceTypeTraits<ExecutionSpace>::id);
-  return (device_id << num_instance_bits);
+  return (device_id << (num_instance_bits + num_device_bits));
 }
 template <typename ExecutionSpace>
 inline uint32_t device_id(ExecutionSpace const& space) noexcept {
-  return device_id_root<ExecutionSpace>() + space.impl_instance_id();
+  return device_id_root<ExecutionSpace>() +
+         (DeviceTypeTraits<ExecutionSpace>::device_id(space)
+          << num_instance_bits) +
+         space.impl_instance_id();
 }
 }  // namespace Experimental
 }  // namespace Tools

lib/kokkos/core/src/impl/Kokkos_Serial.cpp

@@ -233,6 +233,15 @@ void SerialSpaceInitializer::print_configuration(std::ostream& msg,
 }
 
 }  // namespace Impl
+
+#ifdef KOKKOS_ENABLE_CXX14
+namespace Tools {
+namespace Experimental {
+constexpr DeviceType DeviceTypeTraits<Serial>::id;
+}
+}  // namespace Tools
+#endif
+
 }  // namespace Kokkos
 
 #else

lib/kokkos/core/src/impl/Kokkos_ViewMapping.hpp

@@ -1005,15 +1005,15 @@ struct ViewOffset<
   /* Cardinality of the domain index space */
   KOKKOS_INLINE_FUNCTION
   constexpr size_type size() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5 *
-           m_dim.N6 * m_dim.N7;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5 * m_dim.N6 * m_dim.N7;
   }
 
   /* Span of the range space */
   KOKKOS_INLINE_FUNCTION
   constexpr size_type span() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5 *
-           m_dim.N6 * m_dim.N7;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5 * m_dim.N6 * m_dim.N7;
   }
 
   KOKKOS_INLINE_FUNCTION constexpr bool span_is_contiguous() const {
@@ -1026,23 +1026,24 @@ struct ViewOffset<
     return m_dim.N0;
   }
   KOKKOS_INLINE_FUNCTION constexpr size_type stride_2() const {
-    return m_dim.N0 * m_dim.N1;
+    return size_type(m_dim.N0) * m_dim.N1;
   }
   KOKKOS_INLINE_FUNCTION constexpr size_type stride_3() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2;
   }
   KOKKOS_INLINE_FUNCTION constexpr size_type stride_4() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3;
   }
   KOKKOS_INLINE_FUNCTION constexpr size_type stride_5() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4;
   }
   KOKKOS_INLINE_FUNCTION constexpr size_type stride_6() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5;
   }
   KOKKOS_INLINE_FUNCTION constexpr size_type stride_7() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5 *
-           m_dim.N6;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5 * m_dim.N6;
   }
 
   // Stride with [ rank ] value is the total length
@@ -1288,8 +1289,8 @@ struct ViewOffset<
   /* Cardinality of the domain index space */
   KOKKOS_INLINE_FUNCTION
   constexpr size_type size() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5 *
-           m_dim.N6 * m_dim.N7;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5 * m_dim.N6 * m_dim.N7;
   }
 
   /* Span of the range space */
@@ -1633,15 +1634,15 @@ struct ViewOffset<
   /* Cardinality of the domain index space */
   KOKKOS_INLINE_FUNCTION
   constexpr size_type size() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5 *
-           m_dim.N6 * m_dim.N7;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5 * m_dim.N6 * m_dim.N7;
   }
 
   /* Span of the range space */
   KOKKOS_INLINE_FUNCTION
   constexpr size_type span() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5 *
-           m_dim.N6 * m_dim.N7;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5 * m_dim.N6 * m_dim.N7;
   }
 
   KOKKOS_INLINE_FUNCTION constexpr bool span_is_contiguous() const {
@@ -1916,14 +1917,14 @@ struct ViewOffset<
   /* Cardinality of the domain index space */
   KOKKOS_INLINE_FUNCTION
   constexpr size_type size() const {
-    return m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 * m_dim.N5 *
-           m_dim.N6 * m_dim.N7;
+    return size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
+           m_dim.N5 * m_dim.N6 * m_dim.N7;
   }
 
   /* Span of the range space */
   KOKKOS_INLINE_FUNCTION
   constexpr size_type span() const {
-    return size() > 0 ? m_dim.N0 * m_stride : 0;
+    return size() > 0 ? size_type(m_dim.N0) * m_stride : 0;
   }
 
   KOKKOS_INLINE_FUNCTION constexpr bool span_is_contiguous() const {
@@ -2066,27 +2067,29 @@ struct ViewOffset<
                 stride(/* 2 <= rank */
                        m_dim.N1 *
                        (dimension_type::rank == 2
-                            ? 1
+                            ? size_t(1)
                             : m_dim.N2 *
                                   (dimension_type::rank == 3
-                                       ? 1
+                                       ? size_t(1)
                                        : m_dim.N3 *
                                              (dimension_type::rank == 4
-                                                  ? 1
+                                                  ? size_t(1)
                                                   : m_dim.N4 *
                                                         (dimension_type::rank ==
                                                                  5
-                                                             ? 1
+                                                             ? size_t(1)
                                                              : m_dim.N5 *
                                                                    (dimension_type::
                                                                                 rank ==
                                                                             6
-                                                                        ? 1
+                                                                        ? size_t(
+                                                                              1)
                                                                         : m_dim.N6 *
                                                                               (dimension_type::
                                                                                            rank ==
                                                                                        7
-                                                                                   ? 1
+                                                                                   ? size_t(
+                                                                                         1)
                                                                                    : m_dim
                                                                                          .N7)))))))) {
   }
@@ -2447,8 +2450,8 @@ struct ViewOffset<Dimension, Kokkos::LayoutStride, void> {
   constexpr size_type size() const {
     return dimension_type::rank == 0
                ? 1
-               : m_dim.N0 * m_dim.N1 * m_dim.N2 * m_dim.N3 * m_dim.N4 *
-                     m_dim.N5 * m_dim.N6 * m_dim.N7;
+               : size_type(m_dim.N0) * m_dim.N1 * m_dim.N2 * m_dim.N3 *
+                     m_dim.N4 * m_dim.N5 * m_dim.N6 * m_dim.N7;
   }
 
  private:

lib/kokkos/core/src/impl/Kokkos_ViewTracker.hpp

@@ -91,6 +91,7 @@ struct ViewTracker {
 
   template <class RT, class... RP>
   KOKKOS_INLINE_FUNCTION void assign(const View<RT, RP...>& vt) noexcept {
+    if (this == reinterpret_cast<const ViewTracker*>(&vt.m_track)) return;
     KOKKOS_IF_ON_HOST((
         if (view_traits::is_managed && Kokkos::Impl::SharedAllocationRecord<
                                            void, void>::tracking_enabled()) {
@@ -102,6 +103,7 @@ struct ViewTracker {
 
   KOKKOS_INLINE_FUNCTION ViewTracker& operator=(
       const ViewTracker& rhs) noexcept {
+    if (this == &rhs) return *this;
     KOKKOS_IF_ON_HOST((
         if (view_traits::is_managed && Kokkos::Impl::SharedAllocationRecord<
                                            void, void>::tracking_enabled()) {

lib/kokkos/core/unit_test/TestViewAPI.hpp

@@ -1087,6 +1087,20 @@ class TestViewAPI {
     dView4_unmanaged unmanaged_dx = dx;
     ASSERT_EQ(dx.use_count(), 1);
 
+    // Test self assignment
+#if defined(__clang__)
+#pragma GCC diagnostic push
+#pragma GCC diagnostic ignored "-Wself-assign-overloaded"
+#endif
+    dx = dx;  // copy-assignment operator
+#if defined(__clang__)
+#pragma GCC diagnostic pop
+#endif
+    ASSERT_EQ(dx.use_count(), 1);
+    dx = reinterpret_cast<typename dView4::uniform_type &>(
+        dx);  // conversion assignment operator
+    ASSERT_EQ(dx.use_count(), 1);
+
     dView4_unmanaged unmanaged_from_ptr_dx = dView4_unmanaged(
         dx.data(), dx.extent(0), dx.extent(1), dx.extent(2), dx.extent(3));
 

lib/kokkos/core/unit_test/TestViewAPI_e.hpp

@@ -240,6 +240,35 @@ struct TestViewOverloadResolution {
 TEST(TEST_CATEGORY, view_overload_resolution) {
   TestViewOverloadResolution<TEST_EXECSPACE>::test_function_overload();
 }
+
+template <typename MemorySpace>
+struct TestViewAllocationLargeRank {
+  using ViewType = Kokkos::View<char********, MemorySpace>;
+
+  KOKKOS_FUNCTION void operator()(int) const {
+    size_t idx = v.extent(0) - 1;
+    auto& lhs  = v(idx, idx, idx, idx, idx, idx, idx, idx);
+    lhs        = 42;  // This is where it segfaulted
+  }
+
+  ViewType v;
+};
+
+TEST(TEST_CATEGORY, view_allocation_large_rank) {
+  using ExecutionSpace = typename TEST_EXECSPACE::execution_space;
+  using MemorySpace    = typename TEST_EXECSPACE::memory_space;
+  constexpr int dim    = 16;
+  using FunctorType    = TestViewAllocationLargeRank<MemorySpace>;
+  typename FunctorType::ViewType v("v", dim, dim, dim, dim, dim, dim, dim, dim);
+
+  Kokkos::parallel_for(Kokkos::RangePolicy<ExecutionSpace>(0, 1),
+                       FunctorType{v});
+  typename FunctorType::ViewType v_single(v.data() + v.size() - 1, 1, 1, 1, 1,
+                                          1, 1, 1, 1);
+  auto result =
+      Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace{}, v_single);
+  ASSERT_EQ(result(0, 0, 0, 0, 0, 0, 0, 0), 42);
+}
 }  // namespace Test
 
 #include <TestViewIsAssignable.hpp>

lib/kokkos/core/unit_test/tools/TestEventCorrectness.hpp

@@ -238,13 +238,10 @@ TEST(kokkosp, test_id_gen) {
   using Kokkos::Tools::Experimental::DeviceTypeTraits;
   test_wrapper([&]() {
     Kokkos::DefaultExecutionSpace ex;
-    auto id      = device_id(ex);
-    auto id_ref  = identifier_from_devid(id);
-    auto success = (id_ref.instance_id == ex.impl_instance_id()) &&
-                   (id_ref.device_id ==
-                    static_cast<uint32_t>(
-                        DeviceTypeTraits<Kokkos::DefaultExecutionSpace>::id));
-    ASSERT_TRUE(success);
+    auto id     = device_id(ex);
+    auto id_ref = identifier_from_devid(id);
+    ASSERT_EQ(DeviceTypeTraits<decltype(ex)>::id, id_ref.type);
+    ASSERT_EQ(id_ref.instance_id, ex.impl_instance_id());
   });
 }
 
@@ -253,6 +250,7 @@ TEST(kokkosp, test_id_gen) {
  */
 TEST(kokkosp, test_kernel_sequence) {
   test_wrapper([&]() {
+    Kokkos::DefaultExecutionSpace ex;
     auto root = Kokkos::Tools::Experimental::device_id_root<
         Kokkos::DefaultExecutionSpace>();
     std::vector<FencePayload> expected{
@@ -260,11 +258,10 @@ TEST(kokkosp, test_kernel_sequence) {
         {"named_instance", FencePayload::distinguishable_devices::no,
          root + num_instances},
         {"test_kernel", FencePayload::distinguishable_devices::no,
-         root + num_instances}
+         Kokkos::Tools::Experimental::device_id(ex)}
 
     };
     expect_fence_events(expected, [=]() {
-      Kokkos::DefaultExecutionSpace ex;
       TestFunctor tf;
       ex.fence("named_instance");
       Kokkos::parallel_for(

lib/kokkos/master_history.txt

@@ -27,3 +27,4 @@ tag:  3.4.00     date: 04:26:2021    master: 1fb0c284    release: 5d7738d6
 tag:  3.4.01     date: 05:20:2021    master: 4b97a22f    release: 410b15c8
 tag:  3.5.00     date: 11:19:2021    master: c28a8b03    release: 21b879e4
 tag:  3.6.00     date: 04:14:2022    master: 2834f94a    release: 6ea708ff
+tag:  3.6.01     date: 06:16:2022    master: b52f8c83    release: afe9b404

lib/plumed/Install.py

@@ -17,7 +17,7 @@ parser = ArgumentParser(prog='Install.py',
 
 # settings
 
-version = "2.7.4"
+version = "2.8.1"
 mode = "static"
 
 # help message
@@ -58,7 +58,9 @@ checksums = { \
         '2.7.2' : 'cfa0b4dd90a81c25d3302e8d97bfeaea', \
         '2.7.3' : 'f00cc82edfefe6bb3df934911dbe32fb', \
         '2.7.4' : 'f858e0b6aed173748fc85b6bc8a9dcb3', \
+        '2.7.5' : '2aca1986d6c1ca3ba7e9eb51b1102792', \
         '2.8.0' : '489b23daba70da78cf0506cbc31689c6', \
+        '2.8.1' : '6bfe72ebdae63dc38a9ca27d9b0e08f8', \
         }
 
 # parse and process arguments

python/lammps/mliap/__init__.py

@@ -4,15 +4,27 @@
 # try to improperly start up a new interpreter.
 import sysconfig
 import ctypes
-library = sysconfig.get_config_vars('INSTSONAME')[0]
+import platform
+import warnings
+
+py_ver = sysconfig.get_config_vars('VERSION')[0]
+OS_name = platform.system()
+
+if OS_name == "Darwin":
+    SHLIB_SUFFIX = '.dylib'
+    library = 'libpython' + py_ver + SHLIB_SUFFIX
+elif OS_name == "Windows":
+    SHLIB_SUFFIX = '.dll'
+    library = 'python' + py_ver + SHLIB_SUFFIX
+else:
+    SHLIB_SUFFIX = '.so'
+    library = 'libpython' + py_ver + SHLIB_SUFFIX
+
 try:
     pylib = ctypes.CDLL(library)
-except OSError as e:
-    if pylib.endswith(".a"):
-        pylib.strip(".a") + ".so"
-        pylib = ctypes.CDLL(library)
-    else:
-        raise e
+except Exception as e:
+    raise OSError("Unable to locate python shared library") from e
+
 if not pylib.Py_IsInitialized():
     raise RuntimeError("This interpreter is not compatible with python-based mliap for LAMMPS.")
 del sysconfig, ctypes, library, pylib

src/CLASS2/pair_lj_class2.cpp

@@ -513,6 +513,7 @@ double PairLJClass2::init_one(int i, int j)
         pow(sigma[j][j], 3.0) / (pow(sigma[i][i], 6.0) + pow(sigma[j][j], 6.0));
     sigma[i][j] = pow((0.5 * (pow(sigma[i][i], 6.0) + pow(sigma[j][j], 6.0))), 1.0 / 6.0);
     cut[i][j] = mix_distance(cut[i][i], cut[j][j]);
+    did_mix = true;
   }
 
   lj1[i][j] = 18.0 * epsilon[i][j] * pow(sigma[i][j], 9.0);

src/CLASS2/pair_lj_class2_coul_cut.cpp

@@ -277,6 +277,7 @@ double PairLJClass2CoulCut::init_one(int i, int j)
     sigma[i][j] = pow((0.5 * (pow(sigma[i][i], 6.0) + pow(sigma[j][j], 6.0))), 1.0 / 6.0);
     cut_lj[i][j] = mix_distance(cut_lj[i][i], cut_lj[j][j]);
     cut_coul[i][j] = mix_distance(cut_coul[i][i], cut_coul[j][j]);
+    did_mix = true;
   }
 
   double cut = MAX(cut_lj[i][j], cut_coul[i][j]);

src/CLASS2/pair_lj_class2_coul_long.cpp

@@ -711,6 +711,7 @@ double PairLJClass2CoulLong::init_one(int i, int j)
         pow(sigma[j][j], 3.0) / (pow(sigma[i][i], 6.0) + pow(sigma[j][j], 6.0));
     sigma[i][j] = pow((0.5 * (pow(sigma[i][i], 6.0) + pow(sigma[j][j], 6.0))), 1.0 / 6.0);
     cut_lj[i][j] = mix_distance(cut_lj[i][i], cut_lj[j][j]);
+    did_mix = true;
   }
 
   double cut = MAX(cut_lj[i][j], cut_coul);

src/DIELECTRIC/atom_vec_dielectric.cpp

@@ -15,19 +15,23 @@
 
 #include "atom.h"
 #include "citeme.h"
+#include "error.h"
+#include "force.h"
+#include "pair.h"
+#include "pair_hybrid.h"
 
 #include <cmath>
 
 using namespace LAMMPS_NS;
 
 static const char cite_user_dielectric_package[] =
-    "DIELECTRIC package:\n\n"
+    "DIELECTRIC package: doi:10.1016/j.cpc.2019.03.006\n\n"
     "@Article{TrungCPC19,\n"
-    " author = {Trung Dac Nguyen, Honghao Li, Debarshee Bagchi,"
-    " Francisco J. Solis, Monica Olvera de la Cruz,\n"
-    " title = {Incorporating surface polarization effects into large-scale"
-    " coarse-grained Molecular Dynamics simulation},\n"
-    " journal = {Comp.~Phys.~Comm.},\n"
+    " author = {Trung Dac Nguyen and Honghao Li and Debarshee Bagchi and"
+    "   Francisco J. Solis and Olvera de la Cruz, Monica}\n"
+    " title = {Incorporating Surface Polarization Effects Into Large-Scale\n"
+    "   Coarse-Grained Molecular Dynamics Simulation},\n"
+    " journal = {Comput.\\ Phys.\\ Commun.},\n"
     " year =    2019,\n"
     " volume =  241,\n"
     " pages =   {80--91}\n"
@@ -87,6 +91,42 @@ AtomVecDielectric::AtomVecDielectric(LAMMPS *_lmp) : AtomVec(_lmp)
   bond_per_atom = angle_per_atom = dihedral_per_atom = improper_per_atom = 0;
 }
 
+/* ---------------------------------------------------------------------- */
+
+void AtomVecDielectric::init()
+{
+  AtomVec::init();
+
+  // since atom style dielectric modifies the charge q, it will produce incorrect results
+  // with pair styles using coulomb without dielectric support.
+
+  std::string pair_style(force->pair_style);
+  if ((pair_style != "none") && (pair_style != "zero") && !utils::strmatch(force->pair_style,"/dielectric")) {
+    bool mismatch = false;
+    if (utils::strmatch(force->pair_style,"^reaxff")) mismatch = true;
+    if (utils::strmatch(force->pair_style,"^comb")) mismatch = true;
+    if (utils::strmatch(force->pair_style,"coul")) mismatch = true;
+    if (utils::strmatch(force->pair_style,"tip4p")) mismatch = true;
+    if (utils::strmatch(force->pair_style,"dipole")) mismatch = true;
+
+    if (utils::strmatch(force->pair_style,"^hybrid")) {
+      auto hybrid = dynamic_cast<PairHybrid *>(force->pair);
+      if (hybrid) {
+        for (int i = 0; i < hybrid->nstyles; i++) {
+          if (utils::strmatch(hybrid->keywords[i],"^reaxff")) mismatch = true;
+          if (utils::strmatch(hybrid->keywords[i],"^comb")) mismatch = true;
+          if (utils::strmatch(hybrid->keywords[i],"coul")) mismatch = true;
+          if (utils::strmatch(hybrid->keywords[i],"tip4p")) mismatch = true;
+          if (utils::strmatch(hybrid->keywords[i],"dipole")) mismatch = true;
+        }
+      }
+    }
+    if (mismatch)
+      error->all(FLERR, "Pair style {} is not compatible with atom style {}",
+                 pair_style, atom->get_style());
+  }
+}
+
 /* ----------------------------------------------------------------------
    set local copies of all grow ptrs used by this class, except defaults
    needed in replicate when 2 atom classes exist and it calls pack_restart()
@@ -147,6 +187,24 @@ void AtomVecDielectric::data_atom_post(int ilocal)
   mu_one[3] = sqrt(mu_one[0] * mu_one[0] + mu_one[1] * mu_one[1] + mu_one[2] * mu_one[2]);
 }
 
+/* ----------------------------------------------------------------------
+   restore original data for writing the data file
+------------------------------------------------------------------------- */
+
+void AtomVecDielectric::pack_data_pre(int ilocal)
+{
+  atom->q[ilocal] = q_unscaled[ilocal];
+}
+
+/* ----------------------------------------------------------------------
+   undo restore and get back to post read data state
+------------------------------------------------------------------------- */
+
+void AtomVecDielectric::pack_data_post(int ilocal)
+{
+  atom->q[ilocal] /= epsilon[ilocal];
+}
+
 /* ----------------------------------------------------------------------
    initialize other atom quantities after AtomVec::unpack_restart()
 ------------------------------------------------------------------------- */

src/DIELECTRIC/atom_vec_dielectric.h

@@ -31,12 +31,15 @@ class AtomVecDielectric : public AtomVec {
  public:
   AtomVecDielectric(class LAMMPS *);
 
+  void init() override;
   void grow_pointers() override;
   void create_atom_post(int) override;
   void data_atom_post(int) override;
   void unpack_restart_init(int) override;
   int property_atom(const std::string &) override;
   void pack_property_atom(int, double *, int, int) override;
+  void pack_data_pre(int) override;
+  void pack_data_post(int) override;
 
  protected:
   int *num_bond, *num_angle, *num_dihedral, *num_improper;

src/DIELECTRIC/pair_coul_cut_dielectric.cpp

@@ -35,10 +35,10 @@ static constexpr double EPSILON = 1.0e-6;
 
 /* ---------------------------------------------------------------------- */
 
-PairCoulCutDielectric::PairCoulCutDielectric(LAMMPS *_lmp) : PairCoulCut(_lmp)
+PairCoulCutDielectric::PairCoulCutDielectric(LAMMPS *_lmp) : PairCoulCut(_lmp), efield(nullptr)
 {
-  efield = nullptr;
   nmax = 0;
+  no_virial_fdotr_compute = 1;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -54,7 +54,7 @@ void PairCoulCutDielectric::compute(int eflag, int vflag)
 {
   int i, j, ii, jj, inum, jnum, itype, jtype;
   double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, ecoul;
-  double fpair_i, fpair_j;
+  double fpair_i;
   double rsq, r2inv, rinv, forcecoul, factor_coul, efield_i;
   int *ilist, *jlist, *numneigh, **firstneigh;
 
@@ -75,9 +75,7 @@ void PairCoulCutDielectric::compute(int eflag, int vflag)
   double *curvature = atom->curvature;
   double *area = atom->area;
   int *type = atom->type;
-  int nlocal = atom->nlocal;
   double *special_coul = force->special_coul;
-  int newton_pair = force->newton_pair;
   double qqrd2e = force->qqrd2e;
 
   inum = list->inum;
@@ -137,16 +135,9 @@ void PairCoulCutDielectric::compute(int eflag, int vflag)
         efield[i][1] += dely * efield_i;
         efield[i][2] += delz * efield_i;
 
-        if (newton_pair && j >= nlocal) {
-          fpair_j = factor_coul * eps[j] * forcecoul * r2inv;
-          f[j][0] -= delx * fpair_j;
-          f[j][1] -= dely * fpair_j;
-          f[j][2] -= delz * fpair_j;
-        }
-
         if (eflag) {
-          ecoul = factor_coul * qqrd2e * scale[itype][jtype] * qtmp * q[j] * (etmp + eps[j]) * rinv;
-          ecoul *= 0.5;
+          ecoul = factor_coul * qqrd2e * scale[itype][jtype] * qtmp * q[j] * 0.5 * (etmp + eps[j]) *
+              rinv;
         }
         if (evflag) ev_tally_full(i, 0.0, ecoul, fpair_i, delx, dely, delz);
       }

src/DIELECTRIC/pair_coul_long_dielectric.cpp

@@ -36,10 +36,12 @@ using MathConst::MY_PIS;
 
 /* ---------------------------------------------------------------------- */
 
-PairCoulLongDielectric::PairCoulLongDielectric(LAMMPS *_lmp) : PairCoulLong(_lmp)
+PairCoulLongDielectric::PairCoulLongDielectric(LAMMPS *_lmp) :
+ PairCoulLong(_lmp), efield(nullptr)
 {
-  efield = nullptr;
   nmax = 0;
+  single_enable = 0;
+  no_virial_fdotr_compute = 1;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -55,15 +57,11 @@ void PairCoulLongDielectric::compute(int eflag, int vflag)
 {
   int i, j, ii, jj, inum, jnum, itable, itype, jtype;
   double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, ecoul;
-  double fpair_i, fpair_j;
+  double fpair_i;
   double fraction, table;
-  double r, r2inv, forcecoul, factor_coul;
+  double r, rsq, r2inv, forcecoul, factor_coul;
   double grij, expm2, prefactor, t, erfc, prefactorE, efield_i;
   int *ilist, *jlist, *numneigh, **firstneigh;
-  double rsq;
-
-  ecoul = 0.0;
-  ev_init(eflag, vflag);
 
   if (atom->nmax > nmax) {
     memory->destroy(efield);
@@ -71,17 +69,18 @@ void PairCoulLongDielectric::compute(int eflag, int vflag)
     memory->create(efield, nmax, 3, "pair:efield");
   }
 
+  ecoul = 0.0;
+  ev_init(eflag, vflag);
+
   double **x = atom->x;
   double **f = atom->f;
   double *q = atom->q;
+  double *eps = atom->epsilon;
   double **norm = atom->mu;
   double *curvature = atom->curvature;
   double *area = atom->area;
-  double *eps = atom->epsilon;
   int *type = atom->type;
-  int nlocal = atom->nlocal;
   double *special_coul = force->special_coul;
-  int newton_pair = force->newton_pair;
   double qqrd2e = force->qqrd2e;
 
   inum = list->inum;
@@ -103,6 +102,7 @@ void PairCoulLongDielectric::compute(int eflag, int vflag)
     jnum = numneigh[i];
 
     // self term Eq. (55) for I_{ii} and Eq. (52) and in Barros et al
+
     double curvature_threshold = sqrt(area[i]);
     if (curvature[i] < curvature_threshold) {
       double sf = curvature[i] / (4.0 * MY_PIS * curvature_threshold) * area[i] * q[i];
@@ -169,23 +169,15 @@ void PairCoulLongDielectric::compute(int eflag, int vflag)
         efield[i][1] += dely * efield_i;
         efield[i][2] += delz * efield_i;
 
-        if (newton_pair && j >= nlocal) {
-          fpair_j = eps[j] * forcecoul * r2inv;
-          f[j][0] -= delx * fpair_j;
-          f[j][1] -= dely * fpair_j;
-          f[j][2] -= delz * fpair_j;
-        }
-
         if (eflag) {
           if (!ncoultablebits || rsq <= tabinnersq)
-            ecoul = prefactor * (etmp + eps[j]) * erfc;
+            ecoul = prefactor * 0.5 * (etmp + eps[j]) * erfc;
           else {
             table = etable[itable] + fraction * detable[itable];
-            ecoul = scale[itype][jtype] * qtmp * q[j] * (etmp + eps[j]) * table;
+            ecoul = scale[itype][jtype] * qtmp * q[j] * 0.5 * (etmp + eps[j]) * table;
           }
-          ecoul *= 0.5;
           if (factor_coul < 1.0) ecoul -= (1.0 - factor_coul) * prefactor;
-        }
+        } else ecoul = 0.0;
 
         if (evflag) ev_tally_full(i, 0.0, ecoul, fpair_i, delx, dely, delz);
       }

src/DIELECTRIC/pair_lj_cut_coul_cut_dielectric.cpp

@@ -40,6 +40,7 @@ PairLJCutCoulCutDielectric::PairLJCutCoulCutDielectric(LAMMPS *_lmp) : PairLJCut
   efield = nullptr;
   epot = nullptr;
   nmax = 0;
+  no_virial_fdotr_compute = 1;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -56,8 +57,8 @@ void PairLJCutCoulCutDielectric::compute(int eflag, int vflag)
 {
   int i, j, ii, jj, inum, jnum, itype, jtype;
   double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, evdwl, ecoul;
-  double fpair_i, fpair_j;
-  double rsq, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj, efield_i, epot_i;
+  double fpair_i;
+  double rsq, rinv, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj, efield_i, epot_i;
   int *ilist, *jlist, *numneigh, **firstneigh;
 
   if (atom->nmax > nmax) {
@@ -79,10 +80,8 @@ void PairLJCutCoulCutDielectric::compute(int eflag, int vflag)
   double *curvature = atom->curvature;
   double *area = atom->area;
   int *type = atom->type;
-  int nlocal = atom->nlocal;
   double *special_coul = force->special_coul;
   double *special_lj = force->special_lj;
-  int newton_pair = force->newton_pair;
   double qqrd2e = force->qqrd2e;
 
   inum = list->inum;
@@ -95,15 +94,16 @@ void PairLJCutCoulCutDielectric::compute(int eflag, int vflag)
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
     qtmp = q[i];
+    etmp = eps[i];
     xtmp = x[i][0];
     ytmp = x[i][1];
     ztmp = x[i][2];
-    etmp = eps[i];
     itype = type[i];
     jlist = firstneigh[i];
     jnum = numneigh[i];
 
     // self term Eq. (55) for I_{ii} and Eq. (52) and in Barros et al
+
     double curvature_threshold = sqrt(area[i]);
     if (curvature[i] < curvature_threshold) {
       double sf = curvature[i] / (4.0 * MY_PIS * curvature_threshold) * area[i] * q[i];
@@ -114,7 +114,7 @@ void PairLJCutCoulCutDielectric::compute(int eflag, int vflag)
       efield[i][0] = efield[i][1] = efield[i][2] = 0;
     }
 
-    epot[i] = 0;
+    epot[i] = 0.0;
 
     for (jj = 0; jj < jnum; jj++) {
       j = jlist[jj];
@@ -130,9 +130,10 @@ void PairLJCutCoulCutDielectric::compute(int eflag, int vflag)
 
       if (rsq < cutsq[itype][jtype]) {
         r2inv = 1.0 / rsq;
+        rinv = sqrt(r2inv);
 
         if (rsq < cut_coulsq[itype][jtype] && rsq > EPSILON) {
-          efield_i = qqrd2e * q[j] * sqrt(r2inv);
+          efield_i = qqrd2e * q[j] * rinv;
           forcecoul = qtmp * efield_i;
           epot_i = efield_i;
         } else
@@ -156,19 +157,11 @@ void PairLJCutCoulCutDielectric::compute(int eflag, int vflag)
 
         epot[i] += epot_i;
 
-        if (newton_pair && j >= nlocal) {
-          fpair_j = (factor_coul * eps[j] * forcecoul + factor_lj * forcelj) * r2inv;
-          f[j][0] -= delx * fpair_j;
-          f[j][1] -= dely * fpair_j;
-          f[j][2] -= delz * fpair_j;
-        }
-
         if (eflag) {
           if (rsq < cut_coulsq[itype][jtype]) {
-            ecoul = factor_coul * qqrd2e * qtmp * q[j] * (etmp + eps[j]) * sqrt(r2inv);
+            ecoul = factor_coul * qqrd2e * qtmp * q[j] * 0.5 * (etmp + eps[j]) *rinv;
           } else
             ecoul = 0.0;
-          ecoul *= 0.5;
           if (rsq < cut_ljsq[itype][jtype]) {
             evdwl = r6inv * (lj3[itype][jtype] * r6inv - lj4[itype][jtype]) - offset[itype][jtype];
             evdwl *= factor_lj;

src/DIELECTRIC/pair_lj_cut_coul_debye_dielectric.cpp

@@ -40,6 +40,7 @@ PairLJCutCoulDebyeDielectric::PairLJCutCoulDebyeDielectric(LAMMPS *_lmp) : PairL
   efield = nullptr;
   epot = nullptr;
   nmax = 0;
+  no_virial_fdotr_compute = 1;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -56,7 +57,7 @@ void PairLJCutCoulDebyeDielectric::compute(int eflag, int vflag)
 {
   int i, j, ii, jj, inum, jnum, itype, jtype;
   double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, evdwl, ecoul;
-  double fpair_i, fpair_j;
+  double fpair_i;
   double rsq, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj, efield_i, epot_i;
   double r, rinv, screening;
   int *ilist, *jlist, *numneigh, **firstneigh;
@@ -80,10 +81,8 @@ void PairLJCutCoulDebyeDielectric::compute(int eflag, int vflag)
   double *curvature = atom->curvature;
   double *area = atom->area;
   int *type = atom->type;
-  int nlocal = atom->nlocal;
   double *special_coul = force->special_coul;
   double *special_lj = force->special_lj;
-  int newton_pair = force->newton_pair;
   double qqrd2e = force->qqrd2e;
 
   inum = list->inum;
@@ -159,19 +158,11 @@ void PairLJCutCoulDebyeDielectric::compute(int eflag, int vflag)
         efield[i][2] += delz * efield_i;
         epot[i] += epot_i;
 
-        if (newton_pair && j >= nlocal) {
-          fpair_j = (factor_coul * eps[j] * forcecoul + factor_lj * forcelj) * r2inv;
-          f[j][0] -= delx * fpair_j;
-          f[j][1] -= dely * fpair_j;
-          f[j][2] -= delz * fpair_j;
-        }
-
         if (eflag) {
           if (rsq < cut_coulsq[itype][jtype]) {
-            ecoul = factor_coul * qqrd2e * qtmp * q[j] * (etmp + eps[j]) * rinv * screening;
+            ecoul = factor_coul * qqrd2e * qtmp * q[j] * 0.5 * (etmp + eps[j]) * rinv * screening;
           } else
             ecoul = 0.0;
-          ecoul *= 0.5;
           if (rsq < cut_ljsq[itype][jtype]) {
             evdwl = r6inv * (lj3[itype][jtype] * r6inv - lj4[itype][jtype]) - offset[itype][jtype];
             evdwl *= factor_lj;

src/DIELECTRIC/pair_lj_cut_coul_long_dielectric.cpp

@@ -45,6 +45,7 @@ PairLJCutCoulLongDielectric::PairLJCutCoulLongDielectric(LAMMPS *_lmp) : PairLJC
   efield = nullptr;
   epot = nullptr;
   nmax = 0;
+  no_virial_fdotr_compute = 1;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -61,12 +62,11 @@ void PairLJCutCoulLongDielectric::compute(int eflag, int vflag)
 {
   int i, ii, j, jj, inum, jnum, itype, jtype, itable;
   double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, evdwl, ecoul;
-  double fpair_i, fpair_j;
+  double fpair_i;
   double fraction, table;
-  double r, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj;
+  double r, rsq, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj;
   double grij, expm2, prefactor, t, erfc, prefactorE, efield_i, epot_i;
   int *ilist, *jlist, *numneigh, **firstneigh;
-  double rsq;
 
   evdwl = ecoul = 0.0;
   ev_init(eflag, vflag);
@@ -87,10 +87,8 @@ void PairLJCutCoulLongDielectric::compute(int eflag, int vflag)
   double *curvature = atom->curvature;
   double *area = atom->area;
   int *type = atom->type;
-  int nlocal = atom->nlocal;
   double *special_coul = force->special_coul;
   double *special_lj = force->special_lj;
-  int newton_pair = force->newton_pair;
   double qqrd2e = force->qqrd2e;
 
   inum = list->inum;
@@ -123,7 +121,7 @@ void PairLJCutCoulLongDielectric::compute(int eflag, int vflag)
       efield[i][0] = efield[i][1] = efield[i][2] = 0;
     }
 
-    epot[i] = 0;
+    epot[i] = 0.0;
 
     for (jj = 0; jj < jnum; jj++) {
       j = jlist[jj];
@@ -196,23 +194,14 @@ void PairLJCutCoulLongDielectric::compute(int eflag, int vflag)
 
         epot[i] += epot_i;
 
-        if (newton_pair && j >= nlocal) {
-
-          fpair_j = (forcecoul * eps[j] + factor_lj * forcelj) * r2inv;
-          f[j][0] -= delx * fpair_j;
-          f[j][1] -= dely * fpair_j;
-          f[j][2] -= delz * fpair_j;
-        }
-
         if (eflag) {
           if (rsq < cut_coulsq) {
             if (!ncoultablebits || rsq <= tabinnersq)
-              ecoul = prefactor * (etmp + eps[j]) * erfc;
+              ecoul = prefactor * 0.5 * (etmp + eps[j]) * erfc;
             else {
               table = etable[itable] + fraction * detable[itable];
-              ecoul = qtmp * q[j] * (etmp + eps[j]) * table;
+              ecoul = qtmp * q[j] * 0.5 * (etmp + eps[j]) * table;
             }
-            ecoul *= 0.5;
             if (factor_coul < 1.0) ecoul -= (1.0 - factor_coul) * prefactor;
           } else
             ecoul = 0.0;

src/DIELECTRIC/pair_lj_cut_coul_msm_dielectric.cpp

@@ -47,6 +47,7 @@ PairLJCutCoulMSMDielectric::PairLJCutCoulMSMDielectric(LAMMPS *_lmp) : PairLJCut
   nmax = 0;
   ftmp = nullptr;
   efield = nullptr;
+  no_virial_fdotr_compute = 1;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -63,7 +64,7 @@ void PairLJCutCoulMSMDielectric::compute(int eflag, int vflag)
 {
   int i, ii, j, jj, inum, jnum, itype, jtype, itable;
   double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, evdwl, ecoul, fpair;
-  double fpair_i, fpair_j;
+  double fpair_i;
   double fraction, table;
   double r, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj;
   double egamma, fgamma, prefactor, prefactorE, efield_i;
@@ -106,10 +107,8 @@ void PairLJCutCoulMSMDielectric::compute(int eflag, int vflag)
   double *curvature = atom->curvature;
   double *area = atom->area;
   int *type = atom->type;
-  int nlocal = atom->nlocal;
   double *special_coul = force->special_coul;
   double *special_lj = force->special_lj;
-  int newton_pair = force->newton_pair;
   double qqrd2e = force->qqrd2e;
 
   inum = list->inum;
@@ -208,12 +207,6 @@ void PairLJCutCoulMSMDielectric::compute(int eflag, int vflag)
           efield[i][1] += dely * efield_i;
           efield[i][2] += delz * efield_i;
 
-          if (newton_pair && j >= nlocal) {
-            fpair_j = (forcecoul * eps[j] + factor_lj * forcelj) * r2inv;
-            f[j][0] -= delx * fpair_j;
-            f[j][1] -= dely * fpair_j;
-            f[j][2] -= delz * fpair_j;
-          }
         } else {
 
           // separate LJ and Coulombic forces
@@ -223,11 +216,6 @@ void PairLJCutCoulMSMDielectric::compute(int eflag, int vflag)
           f[i][0] += delx * fpair;
           f[i][1] += dely * fpair;
           f[i][2] += delz * fpair;
-          if (newton_pair) {
-            f[j][0] -= delx * fpair;
-            f[j][1] -= dely * fpair;
-            f[j][2] -= delz * fpair;
-          }
 
           fpair_i = (forcecoul * etmp) * r2inv;
           ftmp[i][0] += delx * fpair_i;
@@ -238,22 +226,15 @@ void PairLJCutCoulMSMDielectric::compute(int eflag, int vflag)
           efield[i][0] += delx * efield_i;
           efield[i][1] += dely * efield_i;
           efield[i][2] += delz * efield_i;
-
-          if (newton_pair && j >= nlocal) {
-            fpair_j = (forcecoul * eps[j]) * r2inv;
-            ftmp[j][0] -= delx * fpair_j;
-            ftmp[j][1] -= dely * fpair_j;
-            ftmp[j][2] -= delz * fpair_j;
-          }
         }
 
         if (eflag) {
           if (rsq < cut_coulsq) {
             if (!ncoultablebits || rsq <= tabinnersq)
-              ecoul = prefactor * (etmp + eps[j]) * egamma;
+              ecoul = prefactor * 0.5 * (etmp + eps[j]) * egamma;
             else {
               table = etable[itable] + fraction * detable[itable];
-              ecoul = qtmp * q[j] * (etmp + eps[j]) * table;
+              ecoul = qtmp * q[j] * 0.5 * (etmp + eps[j]) * table;
             }
             if (factor_coul < 1.0) ecoul -= (1.0 - factor_coul) * prefactor;
           } else