Merge pull request #4177 from lammps/maintenance

Fourth Set of Collected Bug Fixes and Maintenance Updates for 2 August 2023 Stable release

.github/workflows/unittest-macos.yml

@@ -15,7 +15,7 @@ jobs:
   build:
     name: MacOS Unit Test
     if: ${{ github.repository == 'lammps/lammps' }}
-    runs-on: macos-latest
+    runs-on: macos-13
     env:
       CCACHE_DIR: ${{ github.workspace }}/.ccache
 
@@ -43,6 +43,8 @@ jobs:
       working-directory: build
       run: |
         ccache -z
+        python3 -m venv macosenv
+        source macosenv/bin/activate
         python3 -m pip install numpy
         python3 -m pip install pyyaml
         cmake -C ../cmake/presets/clang.cmake \

cmake/CMakeLists.txt

@@ -12,6 +12,11 @@ endif()
 if(POLICY CMP0075)
   cmake_policy(SET CMP0075 NEW)
 endif()
+# set policy to silence warnings about requiring execute permission for find_program
+# we use OLD because the python-config script for the Fedora MinGW cross-compiler requires it currently
+if(POLICY CMP0109)
+  cmake_policy(SET CMP0109 OLD)
+endif()
 # set policy to silence warnings about timestamps of downloaded files. review occasionally if it may be set to NEW
 if(POLICY CMP0135)
   cmake_policy(SET CMP0135 OLD)
@@ -29,6 +34,7 @@ project(lammps CXX)
 set(SOVERSION 0)
 get_property(BUILD_IS_MULTI_CONFIG GLOBAL PROPERTY GENERATOR_IS_MULTI_CONFIG)
 
+include(GNUInstallDirs)
 get_filename_component(LAMMPS_DIR ${CMAKE_CURRENT_SOURCE_DIR}/.. ABSOLUTE)
 get_filename_component(LAMMPS_LIB_BINARY_DIR ${CMAKE_BINARY_DIR}/lib ABSOLUTE)
 # collect all executables and shared libs in the top level build folder
@@ -126,6 +132,19 @@ if((CMAKE_CXX_COMPILER_ID STREQUAL "NVHPC") OR (CMAKE_CXX_COMPILER_ID STREQUAL "
   set(CMAKE_TUNE_DEFAULT "-Minform=severe")
 endif()
 
+# this hack is required to compile fmt lib with CrayClang version 15.0.2
+# CrayClang is only directly recognized by version 3.28 and later
+if(CMAKE_VERSION VERSION_LESS 3.28)
+  get_filename_component(_exe "${CMAKE_CXX_COMPILER}" NAME)
+  if((CMAKE_CXX_COMPILER_ID STREQUAL "Clang") AND (_exe STREQUAL "crayCC"))
+    set(CMAKE_TUNE_DEFAULT "-DFMT_STATIC_THOUSANDS_SEPARATOR")
+  endif()
+else()
+  if(CMAKE_CXX_COMPILER_ID STREQUAL "CrayClang")
+    set(CMAKE_TUNE_DEFAULT "-DFMT_STATIC_THOUSANDS_SEPARATOR")
+  endif()
+endif()
+
 # silence nvcc warnings
 if((PKG_KOKKOS) AND (Kokkos_ENABLE_CUDA) AND NOT (CMAKE_CXX_COMPILER_ID STREQUAL "Clang"))
   set(CMAKE_TUNE_DEFAULT "${CMAKE_TUNE_DEFAULT} -Xcudafe --diag_suppress=unrecognized_pragma")
@@ -156,6 +175,22 @@ if(MSVC)
   add_compile_definitions(_CRT_SECURE_NO_WARNINGS)
 endif()
 
+# warn about potentially problematic GCC compiler versions
+if(CMAKE_CXX_COMPILER_ID STREQUAL "GNU")
+  if (CMAKE_CXX_STANDARD GREATER_EQUAL 17)
+    if (CMAKE_CXX_COMPILER_VERSION VERSION_LESS 9.0)
+      message(WARNING "Using ${CMAKE_CXX_COMPILER_ID} compiler version ${CMAKE_CXX_COMPILER_VERSION} "
+        "with C++17 is not recommended. Please use ${CMAKE_CXX_COMPILER_ID} compiler version 9.x or later")
+    endif()
+  endif()
+  if (CMAKE_CXX_STANDARD GREATER_EQUAL 11)
+    if (CMAKE_CXX_COMPILER_VERSION VERSION_LESS 5.0)
+      message(WARNING "Using ${CMAKE_CXX_COMPILER_ID} compiler version ${CMAKE_CXX_COMPILER_VERSION} "
+        "with C++11 is not recommended. Please use ${CMAKE_CXX_COMPILER_ID} compiler version 5.x or later")
+    endif()
+  endif()
+endif()
+
 # export all symbols when building a .dll file on windows
 if((CMAKE_SYSTEM_NAME STREQUAL "Windows") AND BUILD_SHARED_LIBS)
   set(CMAKE_WINDOWS_EXPORT_ALL_SYMBOLS ON)
@@ -198,7 +233,6 @@ else()
   unset(CMAKE_CXX_CLANG_TIDY CACHE)
 endif()
 
-include(GNUInstallDirs)
 file(GLOB ALL_SOURCES ${CONFIGURE_DEPENDS} ${LAMMPS_SOURCE_DIR}/[^.]*.cpp)
 file(GLOB MAIN_SOURCES ${CONFIGURE_DEPENDS} ${LAMMPS_SOURCE_DIR}/main.cpp)
 list(REMOVE_ITEM ALL_SOURCES ${MAIN_SOURCES})
@@ -212,6 +246,10 @@ endif()
 add_executable(lmp ${MAIN_SOURCES})
 target_link_libraries(lmp PRIVATE lammps)
 set_target_properties(lmp PROPERTIES OUTPUT_NAME ${LAMMPS_BINARY})
+# re-export all symbols for plugins
+if(PKG_PLUGIN AND (NOT ((CMAKE_SYSTEM_NAME STREQUAL "Windows"))))
+  set_target_properties(lmp PROPERTIES ENABLE_EXPORTS TRUE)
+endif()
 install(TARGETS lmp EXPORT LAMMPS_Targets DESTINATION ${CMAKE_INSTALL_BINDIR})
 
 option(CMAKE_VERBOSE_MAKEFILE "Generate verbose Makefiles" OFF)
@@ -426,6 +464,7 @@ if(BUILD_OMP)
       (CMAKE_CXX_COMPILER_ID STREQUAL "IntelLLVM") OR (CMAKE_CXX_COMPILER_ID STREQUAL "XLClang") OR
       ((CMAKE_CXX_COMPILER_ID STREQUAL "AppleClang") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 10.0)) OR
       ((CMAKE_CXX_COMPILER_ID STREQUAL "Clang") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 10.0)) OR
+      ((CMAKE_CXX_COMPILER_ID STREQUAL "CrayClang") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 10.0)) OR
       ((CMAKE_CXX_COMPILER_ID STREQUAL "Intel") AND (CMAKE_CXX_COMPILER_VERSION VERSION_GREATER_EQUAL 19.0)))
     # GCC 9.x and later plus Clang 10.x and later implement strict OpenMP 4.0 semantics for consts.
     # Intel 18.0 was tested to support both, so we switch to OpenMP 4+ from 19.x onward to be safe.

cmake/Modules/ExternalCMakeProject.cmake

@@ -43,5 +43,5 @@ function(ExternalCMakeProject target url hash basedir cmakedir cmakefile)
       "${CMAKE_BINARY_DIR}/_deps/${target}-src/${cmakedir}/CMakeLists.txt")
   endif()
   add_subdirectory("${CMAKE_BINARY_DIR}/_deps/${target}-src/${cmakedir}"
-    "${CMAKE_BINARY_DIR}/_deps/${target}-build")
+    "${CMAKE_BINARY_DIR}/_deps/${target}-build" EXCLUDE_FROM_ALL)
 endfunction(ExternalCMakeProject)

cmake/Modules/Packages/KOKKOS.cmake

@@ -132,8 +132,12 @@ if(PKG_KSPACE)
                                  ${KOKKOS_PKG_SOURCES_DIR}/remap_kokkos.cpp)
   if(Kokkos_ENABLE_CUDA)
     if(NOT (FFT STREQUAL "KISS"))
+      find_library(CUFFT_LIBRARY cufft)
+      if (CUFFT_LIBRARY STREQUAL "CUFFT_LIBRARY-NOTFOUND")
+        message(FATAL_ERROR "Required cuFFT library not found. Check your environment or set CUFFT_LIBRARY to its location")
+      endif()
       target_compile_definitions(lammps PRIVATE -DFFT_CUFFT)
-      target_link_libraries(lammps PRIVATE cufft)
+      target_link_libraries(lammps PRIVATE ${CUFFT_LIBRARY})
     endif()
   elseif(Kokkos_ENABLE_HIP)
     if(NOT (FFT STREQUAL "KISS"))

cmake/Modules/Packages/MDI.cmake

@@ -8,8 +8,8 @@ option(DOWNLOAD_MDI "Download and compile the MDI library instead of using an al
 
 if(DOWNLOAD_MDI)
   message(STATUS "MDI download requested - we will build our own")
-  set(MDI_URL "https://github.com/MolSSI-MDI/MDI_Library/archive/v1.4.16.tar.gz" CACHE STRING "URL for MDI tarball")
-  set(MDI_MD5 "407db44e2d79447ab5c1233af1965f65" CACHE STRING "MD5 checksum for MDI tarball")
+  set(MDI_URL "https://github.com/MolSSI-MDI/MDI_Library/archive/v1.4.26.tar.gz" CACHE STRING "URL for MDI tarball")
+  set(MDI_MD5 "3124bb85259471e2a53a891f04bf697a" CACHE STRING "MD5 checksum for MDI tarball")
   mark_as_advanced(MDI_URL)
   mark_as_advanced(MDI_MD5)
   GetFallbackURL(MDI_URL MDI_FALLBACK)

cmake/Modules/Packages/ML-PACE.cmake

@@ -1,6 +1,12 @@
-set(PACELIB_URL "https://github.com/ICAMS/lammps-user-pace/archive/refs/tags/v.2023.10.04.tar.gz" CACHE STRING "URL for PACE evaluator library sources")
+# PACE library support for ML-PACE package
 
-set(PACELIB_MD5 "70ff79f4e59af175e55d24f3243ad1ff" CACHE STRING "MD5 checksum of PACE evaluator library tarball")
+# set policy to silence warnings about timestamps of downloaded files. review occasionally if it may be set to NEW
+if(POLICY CMP0135)
+  cmake_policy(SET CMP0135 OLD)
+endif()
+
+set(PACELIB_URL "https://github.com/ICAMS/lammps-user-pace/archive/refs/tags/v.2023.11.25.fix.tar.gz" CACHE STRING "URL for PACE evaluator library sources")
+set(PACELIB_MD5 "b45de9a633f42ed65422567e3ce56f9f" CACHE STRING "MD5 checksum of PACE evaluator library tarball")
 mark_as_advanced(PACELIB_URL)
 mark_as_advanced(PACELIB_MD5)
 GetFallbackURL(PACELIB_URL PACELIB_FALLBACK)

cmake/Modules/Packages/ML-QUIP.cmake

@@ -18,14 +18,14 @@ if(DOWNLOAD_QUIP)
     set(temp "${temp}F77FLAGS += -fpp -fixed -fPIC\n")
     set(temp "${temp}F95_PRE_FILENAME_FLAG = -Tf\n")
   elseif(CMAKE_Fortran_COMPILER_ID STREQUAL GNU)
-    set(temp "${temp}FPP=${CMAKE_Fortran_COMPILER} -E -x f95-cpp-input\nOPTIM=${CMAKE_Fortran_FLAGS_${BTYPE}}\n")
+    set(temp "${temp}FPP=${CMAKE_Fortran_COMPILER} -E -x f95-cpp-input\nOPTIM=${CMAKE_Fortran_FLAGS_${BTYPE}} -fmax-stack-var-size=6553600\n")
     set(temp "${temp}DEFINES += -DGETARG_F2003 -DGETENV_F2003 -DGFORTRAN -DFORTRAN_UNDERSCORE\n")
     set(temp "${temp}F95FLAGS += -x f95-cpp-input -ffree-line-length-none -ffree-form -fno-second-underscore -fPIC\n")
     set(temp "${temp}F77FLAGS += -x f77-cpp-input -fno-second-underscore -fPIC\n")
   else()
     message(FATAL_ERROR "The ${CMAKE_Fortran_COMPILER_ID} Fortran compiler is not (yet) supported for building QUIP")
   endif()
-  set(temp "${temp}CFLAGS += -fPIC \nCPLUSPLUSFLAGS += -fPIC\nAR_ADD=src\n")
+  set(temp "${temp}CFLAGS += -fPIC -Wno-return-mismatch \nCPLUSPLUSFLAGS += -fPIC -Wno-return-mismatch\nAR_ADD=src\n")
   set(temp "${temp}MATH_LINKOPTS=")
   foreach(flag ${BLAS_LIBRARIES})
     set(temp "${temp} ${flag}")
@@ -41,6 +41,11 @@ if(DOWNLOAD_QUIP)
   set(temp "${temp}HAVE_TURBOGAP=0\nHAVE_QR=1\nHAVE_THIRDPARTY=0\nHAVE_FX=0\nHAVE_SCME=0\nHAVE_MTP=0\n")
   set(temp "${temp}HAVE_MBD=0\nHAVE_TTM_NF=0\nHAVE_CH4=0\nHAVE_NETCDF4=0\nHAVE_MDCORE=0\nHAVE_ASAP=0\n")
   set(temp "${temp}HAVE_CGAL=0\nHAVE_METIS=0\nHAVE_LMTO_TBE=0\nHAVE_SCALAPACK=0\n")
+  # for gfortran, the -std= flag, if present, *must* be -std=gnu or else the compilation will fail.
+  if(CMAKE_Fortran_COMPILER_ID STREQUAL GNU)
+    string(REGEX REPLACE -std=f[0-9]+ -std=gnu newtemp "${temp}")
+    set(temp "${newtemp}")
+  endif()
   file(WRITE ${CMAKE_BINARY_DIR}/quip.config "${temp}")
 
   message(STATUS "QUIP download via git requested - we will build our own")
@@ -56,7 +61,7 @@ if(DOWNLOAD_QUIP)
     GIT_SUBMODULES "src/fox;src/GAP"
     PATCH_COMMAND ${CMAKE_COMMAND} -E copy_if_different ${CMAKE_BINARY_DIR}/quip.config <SOURCE_DIR>/arch/Makefile.lammps
     CONFIGURE_COMMAND env QUIP_ARCH=lammps make config
-    BUILD_COMMAND env QUIP_ARCH=lammps make libquip
+    BUILD_COMMAND env QUIP_ARCH=lammps make -j1 libquip
     INSTALL_COMMAND ""
     BUILD_IN_SOURCE YES
     BUILD_BYPRODUCTS <SOURCE_DIR>/build/lammps/${CMAKE_STATIC_LIBRARY_PREFIX}quip${CMAKE_STATIC_LIBRARY_SUFFIX}

cmake/Modules/Packages/PLUMED.cmake

@@ -1,5 +1,10 @@
 # Plumed2 support for PLUMED package
 
+# set policy to silence warnings about timestamps of downloaded files. review occasionally if it may be set to NEW
+if(POLICY CMP0135)
+  cmake_policy(SET CMP0135 OLD)
+endif()
+
 if(BUILD_MPI)
   set(PLUMED_CONFIG_MPI "--enable-mpi")
   set(PLUMED_CONFIG_CC  ${CMAKE_MPI_C_COMPILER})
@@ -21,9 +26,11 @@ else()
   set(PLUMED_CONFIG_OMP "--disable-openmp")
 endif()
 
-set(PLUMED_URL "https://github.com/plumed/plumed2/releases/download/v2.8.2/plumed-src-2.8.2.tgz"
+# Note: must also adjust check for supported API versions in
+# fix_plumed.cpp when version changes from v2.n.x to v2.n+1.y
+set(PLUMED_URL "https://github.com/plumed/plumed2/releases/download/v2.9.1/plumed-src-2.9.1.tgz"
   CACHE STRING "URL for PLUMED tarball")
-set(PLUMED_MD5 "599092b6a0aa6fff992612537ad98994" CACHE STRING "MD5 checksum of PLUMED tarball")
+set(PLUMED_MD5 "c3b2d31479c1e9ce211719d40e9efbd7" CACHE STRING "MD5 checksum of PLUMED tarball")
 
 mark_as_advanced(PLUMED_URL)
 mark_as_advanced(PLUMED_MD5)
@@ -75,6 +82,9 @@ if((CMAKE_SYSTEM_NAME STREQUAL "Windows") AND (CMAKE_CROSSCOMPILING))
     DEPENDS plumed_build
     COMMENT "Copying Plumed files"
   )
+  if(CMAKE_PROJECT_NAME STREQUAL "lammps")
+    target_link_libraries(lammps INTERFACE LAMMPS::PLUMED)
+  endif()
 
 else()
 
@@ -149,6 +159,9 @@ else()
     endif()
     set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_INCLUDE_DIRECTORIES ${INSTALL_DIR}/include)
     file(MAKE_DIRECTORY ${INSTALL_DIR}/include)
+    if(CMAKE_PROJECT_NAME STREQUAL "lammps")
+      target_link_libraries(lammps PRIVATE LAMMPS::PLUMED)
+    endif()
   else()
     find_package(PkgConfig REQUIRED)
     pkg_check_modules(PLUMED REQUIRED plumed)
@@ -163,7 +176,9 @@ else()
     endif()
     set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_LINK_LIBRARIES "${PLUMED_LOAD}")
     set_target_properties(LAMMPS::PLUMED PROPERTIES INTERFACE_INCLUDE_DIRECTORIES "${PLUMED_INCLUDE_DIRS}")
+    if(CMAKE_PROJECT_NAME STREQUAL "lammps")
+      target_link_libraries(lammps PUBLIC LAMMPS::PLUMED)
+    endif()
   endif()
 endif()
 
-target_link_libraries(lammps PRIVATE LAMMPS::PLUMED)

cmake/presets/mingw-cross.cmake

@@ -32,7 +32,6 @@ set(WIN_PACKAGES
   FEP
   GPU
   GRANULAR
-  INTEL
   INTERLAYER
   KSPACE
   LEPTON

cmake/presets/oneapi.cmake

@@ -18,11 +18,11 @@ set(MPI_CXX_COMPILER "mpicxx" CACHE STRING "" FORCE)
 
 unset(HAVE_OMP_H_INCLUDE CACHE)
 set(OpenMP_C "icx" CACHE STRING "" FORCE)
-set(OpenMP_C_FLAGS "-qopenmp -qopenmp-simd" CACHE STRING "" FORCE)
+set(OpenMP_C_FLAGS "-qopenmp;-qopenmp-simd" CACHE STRING "" FORCE)
 set(OpenMP_C_LIB_NAMES "omp" CACHE STRING "" FORCE)
 set(OpenMP_CXX "icpx" CACHE STRING "" FORCE)
-set(OpenMP_CXX_FLAGS "-qopenmp -qopenmp-simd" CACHE STRING "" FORCE)
+set(OpenMP_CXX_FLAGS "-qopenmp;-qopenmp-simd" CACHE STRING "" FORCE)
 set(OpenMP_CXX_LIB_NAMES "omp" CACHE STRING "" FORCE)
-set(OpenMP_Fortran_FLAGS "-qopenmp -qopenmp-simd" CACHE STRING "" FORCE)
+set(OpenMP_Fortran_FLAGS "-qopenmp;-qopenmp-simd" CACHE STRING "" FORCE)
 set(OpenMP_omp_LIBRARY "libiomp5.so" CACHE PATH "" FORCE)
 

doc/src/Build_extras.rst

@@ -1450,6 +1450,11 @@ in lib/pace or somewhere else, which must be done before building
 LAMMPS with this package. The code for the library can be found
 at: `https://github.com/ICAMS/lammps-user-pace/ <https://github.com/ICAMS/lammps-user-pace/>`_
 
+Instead of including the ML-PACE package directly into LAMMPS, it
+is also possible to skip this step and build the ML-PACE package as
+a plugin using the CMake script files in the ``examples/PACKAGE/pace/plugin``
+folder and then load this plugin at runtime with the :doc:`plugin command <plugin>`.
+
 .. tabs::
 
    .. tab:: CMake build
@@ -1614,6 +1619,11 @@ try a different one, switch to a different build system, consider a
 global PLUMED installation or consider downloading PLUMED during the
 LAMMPS build.
 
+Instead of including the PLUMED package directly into LAMMPS, it
+is also possible to skip this step and build the PLUMED package as
+a plugin using the CMake script files in the ``examples/PACKAGE/plumed/plugin``
+folder and then load this plugin at runtime with the :doc:`plugin command <plugin>`.
+
 .. tabs::
 
    .. tab:: CMake build

doc/src/Developer_plugins.rst

@@ -283,7 +283,7 @@ in the ``examples/kim/plugin`` folder.  No changes to the sources of the
 KIM package themselves are needed; only the plugin interface and loader
 code needs to be added.  This example only supports building with CMake,
 but is probably a more typical example. To compile you need to run CMake
-with -DLAMMPS_SOURCE_DIR=<path/to/lammps/src/folder>.  Other
+with ``-DLAMMPS_SOURCE_DIR=<path/to/lammps/src/folder>``.  Other
 configuration setting are identical to those for compiling LAMMPS.
 
 A second example for a plugin from a package is in the

doc/src/Manual.rst

@@ -144,26 +144,28 @@ Indices and tables
    * :ref:`genindex`
    * :ref:`search`
 
-.. _webbrowser:
-.. admonition:: Web Browser Compatibility
-   :class: note
+.. only:: html
 
-   The HTML version of the manual makes use of advanced features present
-   in "modern" web browsers.  This leads to incompatibilities with older
-   web browsers and specific vendor browsers (e.g. Internet Explorer on Windows)
-   where parts of the pages are not rendered as expected (e.g. the layout is
-   broken or mathematical expressions not typeset).  In that case we
-   recommend to install/use a different/newer web browser or use
-   the `PDF version of the manual <https://docs.lammps.org/Manual.pdf>`_.
+  .. _webbrowser:
+  .. admonition:: Web Browser Compatibility
+     :class: note
 
-   The following web browser versions have been verified to work as
-   expected on Linux, macOS, and Windows where available:
+     The HTML version of the manual makes use of advanced features present
+     in "modern" web browsers.  This leads to incompatibilities with older
+     web browsers and specific vendor browsers (e.g. Internet Explorer on Windows)
+     where parts of the pages are not rendered as expected (e.g. the layout is
+     broken or mathematical expressions not typeset).  In that case we
+     recommend to install/use a different/newer web browser or use
+     the `PDF version of the manual <https://docs.lammps.org/Manual.pdf>`_.
 
-   - Safari version 11.1 and later
-   - Firefox version 54 and later
-   - Chrome version 54 and later
-   - Opera version 41 and later
-   - Edge version 80 and later
+     The following web browser versions have been verified to work as
+     expected on Linux, macOS, and Windows where available:
 
-   Also Android version 7.1 and later and iOS version 11 and later have
-   been verified to render this website as expected.
+     - Safari version 11.1 and later
+     - Firefox version 54 and later
+     - Chrome version 54 and later
+     - Opera version 41 and later
+     - Edge version 80 and later
+
+     Also Android version 7.1 and later and iOS version 11 and later have
+     been verified to render this website as expected.

doc/src/Packages_details.rst

@@ -1797,7 +1797,8 @@ Aidan Thompson^3, Gabor Csanyi^2, Christoph Ortner^4, Ralf Drautz^1.
 **Install:**
 
 This package has :ref:`specific installation instructions <ml-pace>` on the
-:doc:`Build extras <Build_extras>` page.
+:doc:`Build extras <Build_extras>` page.  This package may also be compiled
+as a plugin to avoid licensing conflicts when distributing binaries.
 
 **Supporting info:**
 
@@ -2357,7 +2358,9 @@ and Gareth Tribello.
 
 **Install:**
 
-This package has :ref:`specific installation instructions <plumed>` on the :doc:`Build extras <Build_extras>` page.
+This package has :ref:`specific installation instructions <plumed>` on the
+:doc:`Build extras <Build_extras>` page.  This package may also be compiled
+as a plugin to avoid licensing conflicts when distributing binaries.
 
 **Supporting info:**
 

doc/src/Tools.rst

@@ -802,7 +802,7 @@ and LAMMPS GUI can be launched from anywhere from the command line.
 
 The standard CMake build procedure can be applied and the
 ``mingw-cross.cmake`` preset used. By using ``mingw64-cmake`` the CMake
-command will automatically include a suitable CMake toolset file (the
+command will automatically include a suitable CMake toolchain file (the
 regular cmake command can be used after that).  After building the
 libraries and executables, you can build the target 'zip'
 (i.e. ``cmake --build <build dir> --target zip`` or ``make zip``

doc/src/commands_list.rst

@@ -31,15 +31,6 @@ Commands
    dihedral_write
    dimension
    displace_atoms
-   dump
-   dump_adios
-   dump_cfg_uef
-   dump_h5md
-   dump_image
-   dump_modify
-   dump_molfile
-   dump_netcdf
-   dump_vtk
    dynamical_matrix
    echo
    fix

doc/src/compute_smd_rho.rst

@@ -45,7 +45,8 @@ Restrictions
 """"""""""""
 
 This compute is part of the MACHDYN package. It is only enabled if
-LAMMPS was built with that package. See the :doc:`Build package <Build_package>` page for more info.
+LAMMPS was built with that package. See the :doc:`Build package
+<Build_package>` page for more info.
 
 Related commands
 """"""""""""""""

doc/src/compute_smd_vol.rst

@@ -24,7 +24,7 @@ Description
 """""""""""
 
 Define a computation that provides the per-particle volume and the sum
-of the per-particle volumes of the group for which the fix is defined.
+of the per-particle volumes of the group for which the compute is defined.
 
 See `this PDF guide <PDF/MACHDYN_LAMMPS_userguide.pdf>`_ to using Smooth
 Mach Dynamics in LAMMPS.
@@ -41,13 +41,14 @@ The per-particle vector values will be given in :doc:`units <units>` of
 volume.
 
 Additionally, the compute returns a scalar, which is the sum of the
-per-particle volumes of the group for which the fix is defined.
+per-particle volumes of the group for which the compute is defined.
 
 Restrictions
 """"""""""""
 
 This compute is part of the MACHDYN package.  It is only enabled if
-LAMMPS was built with that package. See the :doc:`Build package <Build_package>` page for more info.
+LAMMPS was built with that package. See the :doc:`Build package
+<Build_package>` page for more info.
 
 Related commands
 """"""""""""""""

doc/src/dumps.rst

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

doc/src/fix_ave_correlate.rst

@@ -65,7 +65,6 @@ Examples
    fix 1 all ave/correlate 1 50 10000 &
              c_thermo_press[1] c_thermo_press[2] c_thermo_press[3] &
              type upper ave running title1 "My correlation data"
-
    fix 1 all ave/correlate 1 50 10000 c_thermo_press[*]
 
 Description

doc/src/fix_ave_correlate_long.rst

@@ -20,11 +20,11 @@ Syntax
   .. parsed-literal::
 
        c_ID = global scalar calculated by a compute with ID
-       c_ID[I] = Ith component of global vector calculated by a compute with ID
+       c_ID[I] = Ith component of global vector calculated by a compute with ID, I can include wildcard (see below)
        f_ID = global scalar calculated by a fix with ID
-       f_ID[I] = Ith component of global vector calculated by a fix with ID
+       f_ID[I] = Ith component of global vector calculated by a fix with ID, I can include wildcard (see below)
        v_name = global value calculated by an equal-style variable with name
-       v_name[I] = Ith component of global vector calculated by a vector-style variable with name
+       v_name[I] = Ith component of a vector-style variable with name, I can include wildcard (see below)
 
 * zero or more keyword/arg pairs may be appended
 * keyword = *type* or *start* or *file* or *overwrite* or *title1* or *title2* or *ncorr* or *nlen* or *ncount*
@@ -63,6 +63,7 @@ Examples
    fix 1 all ave/correlate/long 1 10000 &
              c_thermo_press[1] c_thermo_press[2] c_thermo_press[3] &
              type upper title1 "My correlation data" nlen 15 ncount 3
+   fix 1 all ave/correlate/long 1 10000 c_thermo_press[*]
 
 Description
 """""""""""
@@ -80,8 +81,10 @@ specified values may represent calculations performed by computes and
 fixes which store their own "group" definitions.
 
 Each listed value can be the result of a compute or fix or the
-evaluation of an equal-style variable. See the
-:doc:`fix ave/correlate <fix_ave_correlate>` page for details.
+evaluation of an equal-style or vector-style variable.  For
+vector-style variables, the specified indices can include a wildcard
+character.  See the :doc:`fix ave/correlate <fix_ave_correlate>` page
+for details.
 
 The *Nevery* and *Nfreq* arguments specify on what time steps the input
 values will be used to calculate correlation data and the frequency

doc/src/fix_charge_regulation.rst

@@ -253,11 +253,11 @@ built with that package.  See the :doc:`Build package <Build_package>`
 page for more info.
 
 The :doc:`atom_style <atom_style>`, used must contain the charge
-property, for example, the style could be *charge* or *full*. Only
-usable for 3D simulations. Atoms specified as free ions cannot be part
-of rigid bodies or molecules and cannot have bonding interactions. The
-scheme is limited to integer charges, any atoms with non-integer charges
-will not be considered by the fix.
+property and have per atom type masses, for example, the style could be
+*charge* or *full*. Only usable for 3D simulations.  Atoms specified as
+free ions cannot be part of rigid bodies or molecules and cannot have
+bonding interactions.  The scheme is limited to integer charges, any
+atoms with non-integer charges will not be considered by the fix.
 
 All interaction potentials used must be continuous, otherwise the MD
 integration and the particle exchange MC moves do not correspond to the

doc/src/fix_gcmc.rst

@@ -440,8 +440,11 @@ This fix is part of the MC package.  It is only enabled if LAMMPS was
 built with that package.  See the :doc:`Build package <Build_package>`
 doc page for more info.
 
+This fix style requires an :doc:`atom style <atom_style>` with per atom
+type masses.
+
 Do not set "neigh_modify once yes" or else this fix will never be
-called.  Reneighboring is required.
+called.  Reneighboring is **required**.
 
 Only usable for 3D simulations.
 

doc/src/fix_qeq.rst

@@ -232,8 +232,6 @@ These fixes are part of the QEQ package.  They are only enabled if
 LAMMPS was built with that package.  See the :doc:`Build package
 <Build_package>` page for more info.
 
-These qeq fixes are not compatible with the GPU and USER-INTEL packages.
-
 These qeq fixes will ignore electric field contributions from
 :doc:`fix efield <fix_efield>`.
 

doc/src/fix_sgcmc.rst

@@ -155,6 +155,9 @@ This fix is part of the MC package. It is only enabled if LAMMPS was
 built with that package.  See the :doc:`Build package <Build_package>`
 page for more info.
 
+This fix style requires an :doc:`atom style <atom_style>` with per atom
+type masses.
+
 At present the fix provides optimized subroutines for EAM type
 potentials (see above) that calculate potential energy changes due to
 *local* atom type swaps very efficiently.  Other potentials are

doc/src/fix_widom.rst

@@ -195,8 +195,11 @@ doc page for more info.
 Do not set "neigh_modify once yes" or else this fix will never be
 called.  Reneighboring is **required**.
 
-Can be run in parallel, but aspects of the GCMC part will not scale well
-in parallel. Only usable for 3D simulations.
+This fix style requires an :doc:`atom style <atom_style>` with per atom
+type masses.
+
+Can be run in parallel, but some aspects of the insertion procedure
+will not scale well in parallel. Only usable for 3D simulations.
 
 
 Related commands

doc/src/neigh_modify.rst

@@ -267,6 +267,8 @@ The value of the *page* setting must be at least 10x larger than the
 *one* setting.  This ensures neighbor pages are not mostly empty
 space.
 
+The *exclude group* setting is currently not compatible with dynamic groups.
+
 Related commands
 """"""""""""""""
 

doc/src/pair_dsmc.rst

@@ -138,8 +138,12 @@ This pair style can only be used via the *pair* keyword of the
 Restrictions
 """"""""""""
 
-This style is part of the MC package.  It is only enabled if LAMMPS
-was built with that package.  See the :doc:`Build package <Build_package>` page for more info.
+This pair style is part of the MC package.  It is only enabled if LAMMPS
+was built with that package.  See the :doc:`Build package
+<Build_package>` page for more info.
+
+This pair style requires an :doc:`atom style <atom_style>` with per
+atom type masses.
 
 Related commands
 """"""""""""""""

doc/src/pair_meam.rst

@@ -427,7 +427,7 @@ package. They are only enabled if LAMMPS was built with that package.
 See the :doc:`Build package <Build_package>` page for more info.
 
 The maximum number of elements, that can be read from the MEAM library
-file, is determined at compile time. The default is 5.  If you need
+file, is determined at compile time. The default is 8.  If you need
 support for more elements, you have to change the the constant 'maxelt'
 at the beginning of the file ``src/MEAM/meam.h`` and update/recompile
 LAMMPS.  There is no limit on the number of atoms types.

doc/utils/requirements.txt

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

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

@@ -68,7 +68,7 @@ images_config = {
 templates_path = ['_templates']
 
 # The suffix of source filenames.
-source_suffix = '.rst'
+source_suffix = {'.rst': 'restructuredtext'}
 
 # The encoding of source files.
 #source_encoding = 'utf-8-sig'
@@ -296,6 +296,7 @@ latex_elements = {
 \setcounter{tocdepth}{2}
 \renewcommand{\sfdefault}{ptm}  % Use Times New Roman font for \textrm
 \renewcommand{\sfdefault}{phv}  % Use Helvetica font for \textsf
+\usepackage[columns=1]{idxlayout} % create index with only one column
 % Set up math fonts to match text fonts
 \DeclareSymbolFont{operators}   {OT1}{ptm}{m}{n}
 \DeclareSymbolFont{bold}        {OT1}{ptm}{bx}{n}
@@ -340,6 +341,11 @@ latex_elements = {
 '''
 }
 
+# copy custom style file for tweaking index layout
+latex_additional_files = [
+    'latex/idxlayout.sty', 'latex/ellipse.sty',
+    'latex/pict2e.sty', 'latex/p2e-pdftex.def',
+]
 # Grouping the document tree into LaTeX files. List of tuples
 # (source start file, target name, title,
 #  author, documentclass [howto, manual, or own class]).
@@ -382,6 +388,12 @@ man_pages = [
 #man_show_urls = False
 
 
+# strip off LAMMPS_NS:: from index entries
+cpp_index_common_prefix = [
+    'LAMMPS_NS::',
+    '_LMP_STYLE_CONST::', '_LMP_TYPE_CONST::', '_LMP_VAR_CONST::',
+]
+
 # -- Options for Texinfo output -------------------------------------------
 
 # Grouping the document tree into Texinfo files. List of tuples

doc/utils/sphinx-config/false_positives.txt

@@ -790,6 +790,7 @@ dispersionflag
 dissipative
 Dissipative
 distharm
+distutils
 dl
 dlabel
 dlambda

doc/utils/sphinx-config/latex/ellipse.sty

@@ -0,0 +1,318 @@
+%%
+%% This is file `ellipse.sty',
+%% generated with the docstrip utility.
+%%
+%% The original source files were:
+%%
+%% ellipse.dtx  (with options: `package')
+%% 
+%% Copyright (C) 2015
+%% Daan Leijen
+%% 
+%% This work may be distributed and/or modified under the
+%% conditions of the LaTeX Project Public License, either version 1.3
+%% of this license or (at your option) any later version.
+%% The latest version of this license is in
+%%   http://www.latex-project.org/lppl.txt
+%% and version 1.3 or later is part of all distributions of LaTeX
+%% version 2003/12/01 or later.
+%% 
+%% This work has the LPPL maintenance status "author-maintained".
+%% 
+\NeedsTeXFormat{LaTeX2e}[1999/12/01]
+\ProvidesPackage{ellipse}
+    [2004/11/05 v1.0 .dtx ellipse file]
+\RequirePackage{pict2e}
+
+\providecommand*\pIIe@csedef[1]{\expandafter\edef\csname #1\endcsname}
+\newcommand*\pIIe@ellip@csqrt[3]{%
+  \@ovxx=#1\relax
+  \ifdim\@ovxx<\z@\@ovxx-\@ovxx\fi
+  \@ovyy=#2\relax
+  \ifdim\@ovyy<\z@\@ovyy-\@ovyy\fi
+  \edef\pIIe@csname{@csqrt(\number\@ovxx,\number\@ovyy)}%
+  \expandafter\ifx\csname\pIIe@csname\endcsname\relax
+    \pIIe@ellip@csqrt@%
+    \pIIe@csedef{\pIIe@csname}{\the\dimen@}%
+    #3\dimen@
+  \else
+    #3\dimexpr\csname\pIIe@csname\endcsname\relax
+  \fi
+}
+\newcommand*\pIIe@ellip@csqrt@{%
+  \@ovdx\@ovxx
+  \advance\@ovdx by \@ovyy
+  \dimen@0.7071067\@ovdx
+  \ifdim\dimen@<\@ovyy\dimen@\@ovyy\fi
+  \ifdim\dimen@<\@ovxx\dimen@\@ovxx\fi
+  \ifdim\@ovdx<128\p@
+    \edef\@tempa{\strip@pt\@ovxx}%
+    \@ovxx\@tempa\@ovxx
+    \edef\@tempa{\strip@pt\@ovyy}%
+    \@ovyy\@tempa\@ovyy
+    \advance\@ovxx by \@ovyy
+    \advance\dimen@ by \dimexpr1pt * \@ovxx/\dimen@\relax
+    \divide\dimen@ by 2%
+    \advance\dimen@ by \dimexpr1pt * \@ovxx/\dimen@\relax
+    \divide\dimen@ by 2%
+  \fi
+}
+
+\newcommand*\pIIe@atan@{%
+  \@tempdima\dimen@
+  \@tempdimb\@tempdima
+  \ifdim\@tempdimb<\z@\@tempdimb-\@tempdimb\fi
+  \dimen@0.0663\@tempdimb
+  \advance\dimen@ 0.2447pt\relax
+  \advance\@tempdimb -1pt\relax
+  \edef\@tempa{\strip@pt\@tempdimb}%
+  \dimen@\@tempa\dimen@
+  \edef\@tempa{\strip@pt\@tempdima}%
+  \dimen@\@tempa\dimen@
+  \dimen@-\dimen@
+  \advance\dimen@ 0.7853\@tempdima
+}
+
+\newcommand*\pIIe@atantwo[3]{%
+  \edef\pIIe@csname{@atan2(\number\dimexpr#1\relax,\number\dimexpr#2\relax)}%
+  \expandafter\ifx\csname\pIIe@csname\endcsname\relax
+    \pIIe@atantwo@{#1}{#2}{#3}%
+    \pIIe@csedef{\pIIe@csname}{\the\dimexpr#3\relax}%
+  \else
+    #3\dimexpr\csname\pIIe@csname\endcsname\relax
+  \fi
+}
+
+\newcommand*\pIIe@atantwo@[3]{%
+  \@tempdima\dimexpr#2\relax
+  \@tempdimb\dimexpr#1\relax
+  \ifdim\@tempdima=\z@\relax
+    \ifdim\@tempdimb>\z@\relax\dimen@90\p@
+    \else\ifdim\@tempdimb<\z@\relax\dimen@-90\p@
+    \else\dimen@0\p@
+    \fi\fi
+  \else
+    \@tempdimd\z@
+    \ifdim\@tempdima<\z@\relax
+      \ifdim\@tempdimb<\z@\relax\@tempdimd-180\p@
+      \else\@tempdimd180\p@
+      \fi
+    \fi
+    \dimen@\dimexpr1pt * \@tempdimb/\@tempdima\relax
+    \@tempdimc\dimen@
+    \ifdim\@tempdimc<\z@\relax\@tempdimc-\@tempdimc\fi
+    \ifdim\@tempdimc>\p@\relax
+      \dimen@\dimexpr1pt * \@tempdima/\@tempdimb\relax
+      \ifdim\dimen@<\z@\relax\def\@tempsign{-}\else\def\@tempsign{}\fi
+      \pIIe@atan@
+      \dimen@-\dimen@
+      \advance\dimen@ by \@tempsign1.5707pt\relax
+    \else
+      \pIIe@atan@
+    \fi
+    \dimen@57.29578\dimen@
+    \advance\dimen@ by \@tempdimd
+  \fi
+  #3\dimen@%
+}
+
+\newcommand*\pIIe@noneto[2]{}
+\newcommand*\pIIe@ellip@sincost@[2]{%
+  \CalculateSin{#1}%
+  \CalculateCos{#1}%
+  \@tempdima\UseSin{#1}\p@
+  \@tempdimb\UseCos{#1}\p@
+  \ifdim\@tempdima=\p@\relax
+    \pIIe@csedef{@ellipsin#2}{1}%
+    \pIIe@csedef{@ellipcos#2}{0}%
+  \else\ifdim\@tempdima=-\p@\relax
+    \pIIe@csedef{@ellipsin#2}{-1}%
+    \pIIe@csedef{@ellipcos#2}{0}%
+  \else
+    \@tempdimc\@ellipratio\dimexpr1pt * \@tempdima/\@tempdimb\relax
+    %\typeout{ i#2=\the\@tempdimc, sin(#1)=\the\@tempdima}%
+    \pIIe@ellip@csqrt{\p@}{\@tempdimc}\@tempdimd
+    \ifdim\@tempdimb<\z@\relax\@tempdimd-\@tempdimd\fi
+    \pIIe@csedef{@ellipsin#2}{\strip@pt\dimexpr1pt * \@tempdimc/\@tempdimd\relax}%
+    \pIIe@csedef{@ellipcos#2}{\strip@pt\dimexpr1pt * \p@/\@tempdimd\relax}%
+  \fi\fi
+}
+
+\newcommand*\pIIe@ellip@sincost[2]{%
+  %\typeout{ calc sin cos: angles (#1,#2), radii: (\the\@ovro,\the\@ovri)}%
+  \edef\@ellipratio{\strip@pt\dimexpr1pt * \@ovro/\@ovri\relax}%
+  \pIIe@ellip@sincost@{#1}{one}%
+  \pIIe@ellip@sincost@{#2}{two}%
+  %\typeout{ sincos(a=#1)=(\@ellipsinone,\@ellipcosone), sincos(a=#2)=(\@ellipsintwo,\@ellipcostwo), }%
+}
+\newcommand*\pIIe@omega[3]{%
+  \@tempdima\csname @ellipcos#3\endcsname\@ovro
+  \advance\@tempdima by #1\relax
+  \@tempdimb\csname @ellipsin#3\endcsname\@ovri
+  \advance\@tempdimb by #2\relax
+}
+
+\newcommand*\pIIe@omegai[1]{%
+  \@tempdimc\csname @ellipsin#1\endcsname\@ovro
+  \@tempdimc-\@tempdimc
+  \@tempdimd\csname @ellipcos#1\endcsname\@ovri
+}
+
+\newcommand*\pIIe@ellip@kappa{%
+  \@ovyy\@ellipsinone\p@
+  \@ovxx\@ellipcosone\p@
+  \@tempdima\@ellipcostwo\@ovyy
+  \@tempdima-\@tempdima
+  \advance\@tempdima by \@ellipsintwo\@ovxx
+  \@tempdimb\@ellipcostwo\@ovxx
+  \advance\@tempdimb by \@ellipsintwo\@ovyy
+  \ifdim\@tempdima=\z@\relax
+    \edef\@ellipkappa{0}%
+  \else
+    \dimen@\dimexpr1pt - \@tempdimb\relax
+    \dimen@\dimexpr1pt * \dimen@/\@tempdima\relax
+    \pIIe@ellip@csqrt{2\p@}{1.73205\dimen@}{\dimen@}%
+    \advance\dimen@ by -\p@
+    \divide\dimen@ by 3%
+    \edef\@tempa{\strip@pt\@tempdima}%
+    \dimen@\@tempa\dimen@
+    \edef\@ellipkappa{\strip@pt\dimen@}%
+  \fi
+  %\typeout{ calculated kappa: \@ellipkappa}%
+}
+
+\newcommand*\pIIe@elliparc@[5]{%
+  %\typeout{elliparc: #1, center: (#2, #3), radius (\the\@ovro, \the\@ovri),angle (#4, #5)}%
+  \ifcase #1\relax
+      \let\@ellip@startto\pIIe@lineto
+  \or \let\@ellip@startto\pIIe@moveto
+  \or \let\@ellip@startto\pIIe@noneto%
+  \else\PackageWarning{ellipse}{Illegal initial action in \protect\elliparc: %
+          must be one of 0 (lineto), 1 (moveto) or 2 (do nothing) but I got: #1}%
+  \fi
+  \ifdim\@ovro=\z@\relax\@ovri\z@\fi
+  \ifdim\@ovri=\z@\relax
+    \@ellip@startto{#2}{#3}%
+  \else
+    \pIIe@ellip@sincost{#4}{#5}%
+    \pIIe@elliparc@draw{#2}{#3}%
+  \fi
+}
+\newcommand*\pIIe@elliparc@t[5]{%
+  \ifcase #1\relax
+      \let\@ellip@startto\pIIe@lineto
+  \or \let\@ellip@startto\pIIe@moveto
+  \or \let\@ellip@startto\pIIe@noneto%
+  \else\PackageWarning{ellipse}{Illegal initial action in \protect\elliparc: %
+          must be one of 0 (lineto), 1 (moveto) or 2 (do nothing) but I got: #1}%
+  \fi
+  \ifdim\@ovro=\z@\relax\@ovri\z@\fi
+  \ifdim\@ovri=\z@\relax
+    \@ellip@startto{#2}{#3}%
+  \else
+    \CalculateSin{#4}\CalculateCos{#4}%
+    \edef\@ellipsinone{\UseSin{#4}}%
+    \edef\@ellipcosone{\UseCos{#4}}%
+    \CalculateSin{#5}\CalculateCos{#5}%
+    \edef\@ellipsintwo{\UseSin{#5}}%
+    \edef\@ellipcostwo{\UseCos{#5}}%
+    \pIIe@elliparc@draw{#2}{#3}%
+  \fi
+}
+\newcommand*\pIIe@elliparc@draw[2]{%
+  \pIIe@ellip@kappa%
+  \pIIe@omega{#1}{#2}{one}%
+  %\typeout{ point one: (\the\@tempdima,\the\@tempdimb)}%
+  \@ellip@startto\@tempdima\@tempdimb
+  \pIIe@omegai{one}%
+  \advance\@tempdima by \@ellipkappa\@tempdimc
+  \advance\@tempdimb by \@ellipkappa\@tempdimd
+  \pIIe@add@nums\@tempdima\@tempdimb
+  %\typeout{ control one: (\the\@tempdima,\the\@tempdimb)}%
+  \pIIe@omega{#1}{#2}{two}%
+  \pIIe@omegai{two}%
+  \@tempdimc\@ellipkappa\@tempdimc
+  \@tempdimd\@ellipkappa\@tempdimd
+  \@tempdimc-\@tempdimc
+  \@tempdimd-\@tempdimd
+  \advance\@tempdimc by \@tempdima
+  \advance\@tempdimd by \@tempdimb
+  \pIIe@add@nums\@tempdimc\@tempdimd
+  %\typeout{ control two: (\the\@tempdimc,\the\@tempdimd)}%
+  \pIIe@add@CP\@tempdima\@tempdimb
+  %\typeout{ point two: (\the\@tempdima,\the\@tempdimb)}%
+  \pIIe@addtoGraph\pIIe@curveto@op
+}
+\newcommand*\pIIe@elliparc[7][0]{%
+  \@ovro #4\relax
+  \@ovri #5\relax
+  \iffalse%dim\@ovro=\@ovri
+    \pIIe@arc[#1]{#2}{#3}{#4}{#6}{#7}
+  \else
+    \ifdim \@ovro<\z@ \pIIe@badcircarg\else
+      \ifdim \@ovri<\z@ \pIIe@badcircarg\else
+        \@arclen #7\p@ \advance\@arclen -#6\p@
+        \ifdim \@arclen<\z@ \def\@tempsign{-}\else\def\@tempsign{}\fi
+        \ifdim \@tempsign\@arclen>720\p@
+          \PackageWarning {ellipse}{The arc angle is reduced to -720..720}%
+          \@whiledim \@tempsign\@arclen>720\p@ \do {\advance\@arclen-\@tempsign360\p@}%
+          \@tempdima #6\p@ \advance\@tempdima \@arclen
+          \edef\@angleend{\strip@pt\@tempdima}%
+          \pIIe@@elliparc{#1}{#2}{#3}{#6}{\@angleend}%
+        \else
+          \pIIe@@elliparc{#1}{#2}{#3}{#6}{#7}%
+        \fi
+      \fi
+    \fi
+  \fi
+}
+\newcommand*\pIIe@@elliparc[5]{%
+  \begingroup
+  \ifdim \@tempsign\@arclen>90\p@
+    \divide\@arclen 2%
+    \@tempdima #4\p@\advance\@tempdima by \@arclen
+    \edef\@anglemid{\strip@pt\@tempdima}%
+    \def\@tempa{\pIIe@@elliparc{#1}{#2}{#3}{#4}}%
+    \expandafter\@tempa\expandafter{\@anglemid}%
+    \def\@tempa{\pIIe@@elliparc{2}{#2}{#3}}%
+    \expandafter\@tempa\expandafter{\@anglemid}{#5}%
+  \else
+    \pIIe@elliparc@{#1}{#2}{#3}{#4}{#5}%
+  \fi
+  \endgroup
+}%
+
+\newcommand*\pIIeelliparc[7][0]{%
+  \@killglue
+  \pIIe@elliparc[#1]{#2\unitlength}{#3\unitlength}{#4\unitlength}{#5\unitlength}{#6}{#7}%
+  \ignorespaces%
+}
+\ifx\undefined\elliparc\else
+  \PackageWarning{ellipse}{\protect\elliparc\space is redefined}%
+\fi
+\let\elliparc\pIIeelliparc
+
+\newcommand*\pIIeearc
+  {\@ifstar{\@tempswatrue\pIIe@earc@}{\@tempswafalse\pIIe@earc@}}
+\newcommand*\pIIe@earc@[3][0,360]{\pIIe@earc@@(#1){#2}{#3}}
+\def\pIIe@earc@@(#1,#2)#3#4{%
+  \if@tempswa
+    \pIIe@moveto\z@\z@
+    \pIIe@elliparc{\z@}{\z@}{#3\unitlength}{#4\unitlength}{#1}{#2}%
+    \pIIe@closepath\pIIe@fillGraph
+  \else
+    \pIIe@elliparc[1]{\z@}{\z@}{#3\unitlength}{#4\unitlength}{#1}{#2}%
+    \pIIe@strokeGraph
+  \fi}
+\ifx\undefined\earc\else
+  \PackageWarning{ellipse}{\protect\earc\space is redefined}%
+\fi
+\let\earc\pIIeearc
+
+\newcommand*\pIIeellipse
+  {\@ifstar{\@tempswatrue\pIIe@earc@}{\@tempswafalse\pIIe@earc@}}
+\let\ellipse\pIIeellipse
+
+\endinput
+%%
+%% End of file `ellipse.sty'.

doc/utils/sphinx-config/latex/idxlayout.sty

@@ -0,0 +1,279 @@
+%%
+%% This is file `idxlayout.sty',
+%% generated with the docstrip utility.
+%%
+%% The original source files were:
+%%
+%% idxlayout.dtx  (with options: `package')
+%% 
+%% Copyright (C) 2010--2012 by Thomas Titz <thomas.titz@chello.at>
+%% 
+%% Permission is granted to distribute and/or modify this work under the
+%% terms of the LaTeX Project Public License (LPPL), version 1.3c or
+%% later.
+%% 
+%% The LPPL maintenance status of this work is "maintained".
+%% 
+%% This work consists of the files idxlayout.dtx, idxlayout.ins and
+%% README and the derived files idxlayout.pdf and idxlayout.sty.
+%% 
+\NeedsTeXFormat{LaTeX2e}
+\ProvidesPackage{idxlayout}[2012/03/30 v0.4d Configurable index layout]
+\RequirePackage{etoolbox,kvoptions,multicol}
+\SetupKeyvalOptions{family=ila,prefix=ila@}
+\newcommand*{\ila@defradiokey}[3][]{%
+  \define@key{ila}{#2}[#1]{%
+    \ifcsname ila@#3@##1\endcsname
+      \csname ila@#3@##1\expandafter\endcsname
+    \else
+      \PackageError{idxlayout}{Unknown value ##1 for option #2}%
+    \fi
+  }%
+}
+\newcounter{idxcols}
+\define@key{ila}{columns}{\setcounter{idxcols}{#1}}
+\setkeys{ila}{columns=2}
+\@ifclassloaded{memoir}{%
+  \ifbool{onecolindex}{%
+    \setkeys{ila}{columns=1}%
+  }{%
+  }%
+  \appto{\onecolindex}{\setkeys{ila}{columns=1}}%
+  \appto{\twocolindex}{\setkeys{ila}{columns=2}}%
+}{%
+  \newlength{\indexcolsep}%
+  \setlength{\indexcolsep}{35\p@}%
+  \newlength{\indexrule}%
+  \setlength{\indexrule}{\z@}%
+}
+\define@key{ila}{columnsep}{\setlength{\indexcolsep}{#1}}
+\define@key{ila}{rule}{\setlength{\indexrule}{#1}}
+\DeclareBoolOption{unbalanced}
+\newlength{\ila@indentunit}
+\define@key{ila}{indentunit}{\setlength{\ila@indentunit}{#1}}
+\setkeys{ila}{indentunit=20\p@}
+\DeclareStringOption{hangindent}
+\DeclareStringOption{subindent}
+\DeclareStringOption{subsubindent}
+\renewcommand{\@idxitem}{\par\setlength{\hangindent}{\ila@hangindent}}
+\def\ila@it@abshang{%
+  \renewcommand*{\ila@hangindent}{2\ila@indentunit}%
+  \renewcommand*{\ila@subindent}{\ila@indentunit}%
+  \renewcommand*{\ila@subsubindent}{1.5\ila@indentunit}%
+  \renewcommand{\subitem}{\@idxitem\hspace*{\ila@subindent}}%
+  \renewcommand{\subsubitem}{\@idxitem\hspace*{\ila@subsubindent}}%
+}
+\def\ila@it@relhang{%
+  \renewcommand*{\ila@hangindent}{1.5\ila@indentunit}%
+  \renewcommand*{\ila@subindent}{\ila@indentunit}%
+  \renewcommand*{\ila@subsubindent}{2\ila@indentunit}%
+  \renewcommand{\subitem}{%
+    \par
+    \deflength{\hangindent}{\ila@hangindent + \ila@subindent}%
+    \hspace*{\ila@subindent}%
+  }%
+  \renewcommand{\subsubitem}{%
+    \par
+    \deflength{\hangindent}{\ila@hangindent + \ila@subsubindent}%
+    \hspace*{\ila@subsubindent}%
+  }%
+}
+\newcommand*{\indexsubsdelim}{; }
+\def\ila@it@singlepar{%
+  \renewcommand*{\ila@hangindent}{\ila@indentunit}%
+  \renewcommand{\subitem}{\unskip\indexsubsdelim}%
+  \renewcommand{\subsubitem}{\unskip\indexsubsdelim}%
+}
+\ila@defradiokey{itemlayout}{it}
+\setkeys{ila}{itemlayout=abshang}
+\newlength{\ila@initsep}
+\define@key{ila}{initsep}{\setlength{\ila@initsep}{#1}}
+\setkeys{ila}{initsep=10\p@ \@plus 5\p@ \@minus 3\p@}
+\DeclareStringOption[\ila@initsep]{notesep}
+\DeclareStringOption[80\p@]{minspace}
+\renewcommand{\indexspace}{\par\vspace{\ila@initsep}}
+\DeclareBoolOption{columnnote}
+\newcommand*{\indexfont}{}
+\def\ila@fo@current{\renewcommand*{\indexfont}{}}
+\def\ila@fo@normalsize{\renewcommand*{\indexfont}{\normalsize}}
+\def\ila@fo@small{\renewcommand*{\indexfont}{\small}}
+\def\ila@fo@footnotesize{\renewcommand*{\indexfont}{\footnotesize}}
+\ila@defradiokey{font}{fo}
+\setkeys{ila}{font=current}
+\newcommand*{\indexjustific}{}
+\def\ila@ju@standard{%
+  \renewcommand*{\indexjustific}{%
+    \setlength{\parindent}{\z@}%
+    \setlength{\parfillskip}{\z@ \@plus 1fil}%
+  }%
+}
+\def\ila@ju@raggedright{\renewcommand*{\indexjustific}{\raggedright}}
+\newcommand*{\ila@RaggedRight}{}
+\def\ila@ju@RaggedRight{%
+  \renewcommand*{\indexjustific}{\ila@RaggedRight}%
+}
+\AtEndPreamble{%
+  \IfFileExists{ragged2e.sty}{%
+    \RequirePackage{ragged2e}%
+    \renewcommand*{\ila@RaggedRight}{\RaggedRight}%
+  }{%
+    \PackageWarning{idxlayout}{%
+      Package ragged2e not available, therefore\MessageBreak
+      substituting command raggedright for RaggedRight\MessageBreak
+    }%
+    \renewcommand*{\ila@RaggedRight}{\raggedright}%
+  }%
+}
+\ila@defradiokey{justific}{ju}
+\setkeys{ila}{justific=standard}
+\newcommand{\ila@prenote}{}
+\newcommand{\setindexprenote}[1]{%
+  \def\ila@prenote{%
+    \begingroup#1\par\nobreak\endgroup
+    \vspace{\ila@notesep}%
+  }%
+}
+\newcommand*{\noindexprenote}{\let\ila@prenote\relax}
+\noindexprenote
+\newcommand*{\indexstheadcase}{\MakeUppercase}
+\newcommand*{\ila@classtype}{0}
+\@ifclassloaded{memoir}{%
+  \def\ila@classtype{2}%
+}{%
+  \ifundef{\KOMAClassName}{%
+  }{%
+    \def\ila@classtype{1}%
+  }%
+}
+\ifcase\ila@classtype\relax
+  \DeclareBoolOption{totoc}%
+  \def\ila@prologue{%
+    \ifundef{\chapter}{%
+      \section*{\indexname}%
+      \ifbool{ila@totoc}{%
+        \addcontentsline{toc}{section}{\indexname}%
+      }{%
+      }%
+    }{%
+      \chapter*{\indexname}%
+      \ifbool{ila@totoc}{%
+        \addcontentsline{toc}{chapter}{\indexname}%
+      }{%
+      }%
+    }%
+    \@mkboth{\indexstheadcase\indexname}{\indexstheadcase\indexname}%
+  }
+\or
+  \def\ila@tc@true{\KOMAoptions{index=totoc}}%
+  \def\ila@tc@false{\KOMAoptions{index=nottotoc}}%
+  \ila@defradiokey[true]{totoc}{tc}%
+  \providecommand{\MakeMarkcase}[1]{#1}%
+  \def\ila@prologue{%
+    \ifundef{\chapter}{%
+    }{%
+      \ifundef{\index@preamble}{%
+      }{%
+        \setchapterpreamble{\index@preamble}%
+      }%
+    }%
+    \idx@@heading{\indexname}%
+    \@mkboth{\MakeMarkcase{\indexname}}{\MakeMarkcase{\indexname}}%
+    \ifundef{\chapter}{%
+    }{%
+      \thispagestyle{\indexpagestyle}%
+    }%
+  }
+\or
+  \def\ila@tc@true{\boolfalse{noindexintoc}}%
+  \def\ila@tc@false{\booltrue{noindexintoc}}%
+  \ila@defradiokey[true]{totoc}{tc}%
+  \def\ila@prologue{%
+    \chapter*{\indexname}%
+    \ifbool{noindexintoc}{%
+    }{%
+      \addcontentsline{toc}{chapter}{\indexname}%
+    }%
+    \ifbool{artopt}{%
+    }{%
+      \thispagestyle{indextitlepagestyle}%
+    }%
+    \indexmark
+    \preindexhook
+  }
+\fi
+\@ifpackageloaded{index}{%
+  \@ifpackagelater{index}{2004/01/21}{%
+    \let\ila@packindadjust\relax
+  }{%
+    \def\ila@packindadjust{%
+      \edef\indexname{\the\@nameuse{idxtitle@\@indextype}}%
+      \ifdefempty{\index@prologue}{%
+      }{%
+        \setindexprenote{\index@prologue}%
+      }%
+    }%
+  }%
+}{%
+  \let\ila@packindadjust\relax
+}
+\ProcessKeyvalOptions*
+\newcommand*{\idxlayout}[1]{\setkeys{ila}{#1}}
+\renewenvironment{theindex}{%
+  \ifbool{@twocolumn}{%
+    \boolfalse{@restonecol}%
+    \onecolumn
+  }{%
+    \booltrue{@restonecol}%
+  }%
+  \setlength{\multicolsep}{\z@}%
+  \setlength{\columnsep}{\indexcolsep}%
+  \setlength{\columnseprule}{\indexrule}%
+  \ila@packindadjust
+  \def\ila@prologueplus{%
+    \ila@prologue
+    \indexfont
+    \ifbool{ila@columnnote}{%
+    }{%
+      \ila@prenote
+    }%
+  }%
+  \ifnumcomp{\theidxcols}{<}{\tw@}{%
+    \ila@prologueplus
+  }{%
+    \ifbool{ila@unbalanced}{%
+      \begin{multicols*}{\theidxcols}[\ila@prologueplus][\ila@minspace]%
+    }{%
+      \begin{multicols}{\theidxcols}[\ila@prologueplus][\ila@minspace]%
+    }%
+  }%
+  \setlength{\parskip}{\z@ \@plus 0.3\p@}%
+  \indexjustific
+  \ifundef{\chapter}{%
+    \ifundef{\index@preamble}{%
+    }{%
+      \index@preamble\par\nobreak
+    }%
+  }{%
+  }%
+  \ifbool{ila@columnnote}{%
+    \ila@prenote
+  }{%
+  }%
+  \let\item\@idxitem
+}{%
+  \ifnumcomp{\theidxcols}{<}{\tw@}{%
+  }{%
+    \ifbool{ila@unbalanced}{%
+      \end{multicols*}%
+    }{%
+      \end{multicols}%
+    }%
+  }%
+  \ifbool{@restonecol}{%
+  }{%
+    \twocolumn
+  }%
+}
+\endinput
+%%
+%% End of file `idxlayout.sty'.

doc/utils/sphinx-config/latex/p2e-pdftex.def

@@ -0,0 +1,41 @@
+%%
+%% This is file `p2e-pdftex.def',
+%% generated with the docstrip utility.
+%%
+%% The original source files were:
+%%
+%% p2e-drivers.dtx  (with options: `pdftex')
+%% 
+%% Copyright (C) 2003-2016
+%% Rolf Niepraschk, Rolf.Niepraschk@gmx.de
+%% Hubert Gaesslein,
+%% Josef Tkadlec, j.tkadlec@email.cz
+%% 
+%% This work may be distributed and/or modified under the
+%% conditions of the LaTeX Project Public License, either version 1.3
+%% of this license or (at your option) any later version.
+%% The latest version of this license is in
+%%   http://www.latex-project.org/lppl.txt
+%% and version 1.3 or later is part of all distributions of LaTeX
+%% version 2003/12/01 or later.
+%% 
+%% This work has the LPPL maintenance status "author-maintained".
+%% 
+%% This work consists of all files listed in `manifest.txt'.
+%% 
+\ProvidesFile{p2e-pdftex.def}
+        [2016/02/05 v0.1u
+                          Driver-dependant file (RN,HjG,JT)]
+\begingroup
+  \@ifundefined{pdfoutput}{}{%
+    \ifnum\pdfoutput>0\relax
+      \gdef\pIIe@mode{2}
+    \fi
+  }
+\endgroup
+\ifcase\pIIe@mode\relax \or\or
+  \def\pIIe@code#1{\pdfliteral{ q #1 Q }}
+\fi
+\endinput
+%%
+%% End of file `p2e-pdftex.def'.

doc/utils/sphinx-config/latex/pict2e.sty

@@ -0,0 +1,820 @@
+%%
+%% This is file `pict2e.sty',
+%% generated with the docstrip utility.
+%%
+%% The original source files were:
+%%
+%% pict2e.dtx  (with options: `package')
+%% 
+%% Copyright (C) 2003-2016
+%% Rolf Niepraschk, Rolf.Niepraschk@gmx.de
+%% Hubert Gaesslein,
+%% Josef Tkadlec, j.tkadlec@email.cz
+%% 
+%% This work may be distributed and/or modified under the
+%% conditions of the LaTeX Project Public License, either version 1.3
+%% of this license or (at your option) any later version.
+%% The latest version of this license is in
+%%   http://www.latex-project.org/lppl.txt
+%% and version 1.3 or later is part of all distributions of LaTeX
+%% version 2003/12/01 or later.
+%% 
+%% This work has the LPPL maintenance status "author-maintained".
+%% 
+%% This work consists of all files listed in `manifest.txt'.
+%% 
+\NeedsTeXFormat{LaTeX2e}[1999/12/01]
+\ProvidesPackage{pict2e}[2020/09/30 v0.4b Improved picture commands (HjG,RN,JT)]
+\edef\Gin@codes{%
+ \catcode`\noexpand\^^A\the\catcode`\^^A\relax
+ \catcode`\noexpand\"\the\catcode`\"\relax
+ \catcode`\noexpand\!\the\catcode`\!\relax
+ \catcode`\noexpand\:\the\catcode`\:\relax}
+\catcode`\^^A=\catcode`\%
+\@makeother\"%
+\@makeother\!%
+\@makeother\:%
+\def\@defaultunitsset#1#2#3{%
+  \@defaultunits#1\dimexpr#2#3\relax\relax\@nnil}
+\newcommand*\pIIe@mode{-1}
+\newcommand*\pIIe@code[1]{}
+\providecommand*\Gin@driver{}
+\newcommand*\pIIe@tempa{}
+\newcommand*\pIIe@tempb{}
+\newcommand*\pIIe@tempc{}
+\DeclareOption{dvips}{\def\Gin@driver{dvips.def}}
+\DeclareOption{xdvi}{\ExecuteOptions{dvips}}
+\DeclareOption{dvipdf}{\def\Gin@driver{dvipdf.def}}
+\DeclareOption{dvipdfm}{\def\Gin@driver{dvipdfm.def}}
+\DeclareOption{dvipdfmx}{\def\Gin@driver{dvipdfmx.def}}
+\DeclareOption{pdftex}{\def\Gin@driver{pdftex.def}}
+\DeclareOption{luatex}{\def\Gin@driver{luatex.def}}
+\DeclareOption{xetex}{\def\Gin@driver{xetex.def}}
+\DeclareOption{dvipsone}{\def\Gin@driver{dvipsone.def}}
+\DeclareOption{dviwindo}{\ExecuteOptions{dvipsone}}
+\DeclareOption{oztex}{\ExecuteOptions{dvips}}
+\DeclareOption{textures}{\def\Gin@driver{textures.def}}
+\DeclareOption{pctexps}{\def\Gin@driver{pctexps.def}}
+\DeclareOption{pctex32}{\def\Gin@driver{pctex32.def}}
+\DeclareOption{vtex}{\def\Gin@driver{vtex.def}}
+\DeclareOption{original}{\def\pIIe@mode{0}}
+\newif\ifpIIe@pdfliteral@ok
+\pIIe@pdfliteral@oktrue
+\ifx\pIIe@pdfliteral\@undefined
+  \ifx\pdfliteral\@undefined
+    \pIIe@pdfliteral@okfalse
+    \def\pIIe@pdfliteral#1{%
+      \PackageWarning{pict2e}{pdfliteral not supported}%
+    }%
+  \else
+    \let\pIIe@pdfliteral\pdfliteral
+  \fi
+\fi
+\def\pIIe@buttcap{%
+  \ifpIIe@pdfliteral@ok
+    \buttcap
+  \fi
+}
+\newcommand*\pIIe@FAL{1.52}%
+\newcommand*\pIIe@FAW{3.2}%
+\newcommand*\pIIe@CAW{1.5pt}%
+\newcommand*\pIIe@FAI{0.25}%
+\newcommand*\ltxarrows{%
+  \let\pIIe@vector=\pIIe@vector@ltx
+}
+\newcommand*\pstarrows{%
+  \let\pIIe@vector=\pIIe@vector@pst
+}
+\DeclareOption{ltxarrows}{\AtEndOfPackage{\ltxarrows}}
+\DeclareOption{pstarrows}{\AtEndOfPackage{\pstarrows}}
+\newcommand*\pIIe@debug@comment{}
+\DeclareOption{debug}{%
+  \def\pIIe@debug@comment{^^J^^J\@percentchar\space >>> pict2e <<<^^J}%
+  \begingroup
+    \@ifundefined{pdfcompresslevel}{}{\global\pdfcompresslevel\z@}%
+  \endgroup}
+\DeclareOption{hide}{\AtEndOfPackage{%
+  \let\pIIe@code\@gobble
+}}
+\DeclareOption*{\ExecuteOptions{original}}
+\ExecuteOptions{ltxarrows}
+\InputIfFileExists{pict2e.cfg}{}{}
+\ProcessOptions\relax
+\ifnum\pIIe@mode=\z@
+  \PackageInfo{pict2e}{Package option `original' requested}
+\else
+  \if!\Gin@driver!
+    \PackageError{pict2e}
+      {No driver specified at all}
+      {You should make a default driver option in a file\MessageBreak
+       pict2e.cfg\MessageBreak eg: \protect\ExecuteOptions{dvips}}%
+  \else
+    \PackageInfo{pict2e}{Driver file: \Gin@driver}
+    \@ifundefined{ver@\Gin@driver}{\input{\Gin@driver}}{}
+    \PackageInfo{pict2e}{Driver file for pict2e: p2e-\Gin@driver}
+    \InputIfFileExists{p2e-\Gin@driver}{}{%
+      \PackageError{pict2e}%
+        {Driver file ``p2e-\Gin@driver'' not found}%
+        {Q: Is the file properly installed? A: No!}}
+  \fi
+\fi
+\ifnum\pIIe@mode>\z@
+  \ifnum\pIIe@mode<\thr@@
+    \RequirePackage{trig}
+    \let\pIIe@oldline\line
+    \let\pIIe@old@sline\@sline
+    \let\pIIe@oldvector\vector
+    \let\pIIe@old@circle\@circle
+    \let\pIIe@old@dot\@dot
+    \let\pIIe@old@bezier\@bezier
+    \AtBeginDocument{%
+      \@ifundefined{@cbezier}{%
+        \def\pIIe@old@cbezier[#1](#2,#3)(#4,#5)(#6,#7)(#8,#9){}%
+        }{\let\pIIe@old@cbezier\@cbezier}}
+    \let\pIIe@oldoval\oval
+    \let\pIIe@old@oval\@oval
+    \newcommand*\OriginalPictureCmds{%
+      \let\@sline\pIIe@old@sline
+      \let\line\pIIe@oldline
+      \let\vector\pIIe@oldvector
+      \let\@circle\pIIe@old@circle
+      \let\@dot\pIIe@old@dot
+      \let\@bezier\pIIe@old@bezier
+      \let\@cbezier\pIIe@old@cbezier
+      \renewcommand*\oval[1][]{\pIIe@oldoval}%
+      \let\@oval\pIIe@old@oval
+    }
+  \else
+    \PackageError{pict2e}
+      {Unsupported mode (\pIIe@mode) specified}
+      {The driver you specified requested a mode\MessageBreak
+       not supported by this version of this package}
+  \fi
+\else
+  \ifnum\pIIe@mode<\z@
+    \PackageError{pict2e}
+      {No suitable driver specified}
+      {You should make a default driver option in a file\MessageBreak
+       pict2e.cfg\MessageBreak eg: \protect\ExecuteOptions{dvips}}
+  \fi
+\fi
+\ifnum\pIIe@mode>\z@
+  \ifcase\pIIe@mode\relax
+  \or
+    \newcommand*\pIIe@moveto@op{moveto}
+    \newcommand*\pIIe@lineto@op{lineto}
+    \newcommand*\pIIe@setlinewidth@op{setlinewidth}
+    \newcommand*\pIIe@stroke@op{stroke}
+    \newcommand*\pIIe@fill@op{fill}
+    \newcommand*\pIIe@curveto@op{curveto}
+    \newcommand*\pIIe@concat@op{concat}
+    \newcommand*\pIIe@closepath@op{closepath}
+  \or
+    \newcommand*\pIIe@moveto@op{m}
+    \newcommand*\pIIe@lineto@op{l}
+    \newcommand*\pIIe@setlinewidth@op{w}
+    \newcommand*\pIIe@stroke@op{S}
+    \newcommand*\pIIe@fill@op{f}
+    \newcommand*\pIIe@curveto@op{c}
+    \newcommand*\pIIe@concat@op{cm}
+    \newcommand*\pIIe@closepath@op{h}
+  \fi
+  \@ifdefinable\pIIe@GRAPH{\newtoks\pIIe@GRAPH}
+  \newcommand*\pIIe@addtoGraph[1]{%
+    \begingroup
+      \edef\x{\the\pIIe@GRAPH\space#1}%
+      \global\pIIe@GRAPH\expandafter{\x}%
+    \endgroup}
+  \newcommand*\pIIe@fillGraph{\begingroup \@tempswatrue\pIIe@drawGraph}
+  \newcommand*\pIIe@strokeGraph{\begingroup \@tempswafalse\pIIe@drawGraph}
+  \newcommand*\pIIe@drawGraph{%
+      \edef\x{\pIIe@debug@comment\space
+              \pIIe@scale@PTtoBP}%
+      \if@tempswa
+        \edef\y{\pIIe@fill@op}%
+      \else
+        \edef\x{\x\space
+          \strip@pt\@wholewidth\space\pIIe@setlinewidth@op
+          \pIIe@linecap\pIIe@linejoin\space}%
+        \edef\y{\pIIe@stroke@op}%
+      \fi
+      \expandafter\pIIe@code\expandafter{%
+        \expandafter\x\the\pIIe@GRAPH\space\y}%
+      \global\pIIe@GRAPH{}\xdef\pIIe@CPx{}\xdef\pIIe@CPy{}%
+    \endgroup}
+  \newcommand*\pIIe@CPx{} \newcommand*\pIIe@CPy{}
+  \newcommand*\pIIe@add@CP[2]{%
+    \begingroup
+      \@tempdima#1\xdef\pIIe@CPx{\the\@tempdima}%
+      \@tempdimb#2\xdef\pIIe@CPy{\the\@tempdimb}%
+      \pIIe@addtoGraph{\strip@pt\@tempdima\space\strip@pt\@tempdimb}%
+    \endgroup}
+  \newcommand*\pIIe@add@nums[2]{%
+    \begingroup
+      \@tempdima#1\relax
+      \@tempdimb#2\relax
+      \pIIe@addtoGraph{\strip@pt\@tempdima\space\strip@pt\@tempdimb}%
+    \endgroup}
+  \newcommand*\pIIe@add@num[1]{%
+    \begingroup
+      \@tempdima#1\relax
+      \pIIe@addtoGraph{\strip@pt\@tempdima}%
+    \endgroup}
+  \newcommand*\pIIe@PTtoBP{0.99626401 }
+  \ifcase\pIIe@mode\relax
+  \or
+    \newcommand*\pIIe@concat[6]{%
+      \begingroup
+        \pIIe@addtoGraph{[}%
+        \@tempdima#1\relax \@tempdimb#2\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#3\relax \@tempdimb#4\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#5\relax \@tempdimb#6\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \pIIe@addtoGraph{] \pIIe@concat@op}%
+      \endgroup}
+    \newcommand*\pIIe@translate[2]{\pIIe@add@nums{#1}{#2}\pIIe@addtoGraph{translate}}
+    \newcommand*\pIIe@rotate[1]{\pIIe@add@num{#1}\pIIe@addtoGraph{rotate}}
+    \newcommand*\pIIe@scale[2]{\pIIe@add@nums{#1}{#2}\pIIe@addtoGraph{scale}}
+    \newcommand*\pIIe@scale@PTtoBP{\pIIe@PTtoBP \pIIe@PTtoBP scale}
+  \or
+    \newcommand*\pIIe@concat[6]{%
+      \begingroup
+        \@tempdima#1\relax \@tempdimb#2\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#3\relax \@tempdimb#4\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#5\relax \@tempdimb#6\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \pIIe@addtoGraph\pIIe@concat@op
+      \endgroup}
+    \newcommand*\pIIe@translate[2]{\pIIe@concat\p@\z@\z@\p@{#1}{#2}}
+    \newcommand*\pIIe@rotate[1]{%
+      \begingroup
+        \@tempdima#1\relax
+        \edef\pIIe@tempa{\strip@pt\@tempdima}%
+        \CalculateSin\pIIe@tempa
+        \CalculateCos\pIIe@tempa
+        \edef\pIIe@tempb{\UseSin\pIIe@tempa}%
+        \edef\pIIe@tempc{\UseCos\pIIe@tempa}%
+        \pIIe@concat{\pIIe@tempc\p@}{\pIIe@tempb\p@}%
+          {-\pIIe@tempb\p@}{\pIIe@tempc\p@}\z@\z@
+      \endgroup}
+    \newcommand*\pIIe@scale[2]{\pIIe@concat{#1}\z@\z@{#2}\z@\z@}
+    \newcommand*\pIIe@scale@PTtoBP{\pIIe@PTtoBP 0 0 \pIIe@PTtoBP 0 0 \pIIe@concat@op}
+  \fi
+  \newcommand*\pIIe@moveto[2]{%
+    \pIIe@add@CP{#1}{#2}\pIIe@addtoGraph\pIIe@moveto@op}
+  \newcommand*\pIIe@lineto[2]{%
+    \pIIe@add@CP{#1}{#2}\pIIe@addtoGraph\pIIe@lineto@op}
+  \ifcase\pIIe@mode\relax
+  \or
+    \newcommand*\pIIe@rcurveto[6]{%
+      \begingroup
+        \@tempdima#1\relax \@tempdimb#2\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#3\relax \@tempdimb#4\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#5\relax \@tempdimb#6\relax
+        \pIIe@add@CP\@tempdima\@tempdimb
+        \pIIe@addtoGraph{rcurveto}%
+      \endgroup}
+  \or
+    \newcommand*\pIIe@rcurveto[6]{%
+      \begingroup
+        \@tempdima#1\advance\@tempdima\pIIe@CPx\relax
+        \@tempdimb#2\advance\@tempdimb\pIIe@CPy\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#3\advance\@tempdima\pIIe@CPx\relax
+        \@tempdimb#4\advance\@tempdimb\pIIe@CPy\relax
+        \pIIe@add@nums\@tempdima\@tempdimb
+        \@tempdima#5\advance\@tempdima\pIIe@CPx\relax
+        \@tempdimb#6\advance\@tempdimb\pIIe@CPy\relax
+        \pIIe@add@CP\@tempdima\@tempdimb
+        \pIIe@addtoGraph\pIIe@curveto@op
+      \endgroup}
+  \fi
+  \newcommand*\pIIe@curveto[6]{%
+    \begingroup
+      \@tempdima#1\relax \@tempdimb#2\relax
+      \pIIe@add@nums\@tempdima\@tempdimb
+      \@tempdima#3\relax \@tempdimb#4\relax
+      \pIIe@add@nums\@tempdima\@tempdimb
+      \@tempdima#5\relax \@tempdimb#6\relax
+      \pIIe@add@CP\@tempdima\@tempdimb
+      \pIIe@addtoGraph\pIIe@curveto@op
+    \endgroup}
+  \newcommand*\pIIe@closepath{\pIIe@addtoGraph\pIIe@closepath@op}
+  \newcommand*\pIIe@pyth[3]{%
+    \begingroup
+      \@tempdima=#1\relax
+      \ifnum\@tempdima<\z@\@tempdima=-\@tempdima\fi
+      \@tempdimb=#2\relax
+      \ifnum\@tempdimb<\z@\@tempdimb=-\@tempdimb\fi
+      \advance\@tempdimb\@tempdima
+      \ifnum\@tempdimb=\z@
+        \@tempdimc=\z@
+      \else
+        \multiply\@tempdima 8\relax
+        \pIIe@divide\@tempdima\@tempdimb\@tempdimc
+        \advance\@tempdimc -4pt
+        \multiply\@tempdimc 2
+        \edef\pIIe@tempa{\strip@pt\@tempdimc}%
+        \@tempdima=\pIIe@tempa\@tempdimc
+        \advance\@tempdima 64pt
+        \divide\@tempdima 2\relax
+        \@dashdim=7pt
+        \pIIe@@pyth\pIIe@@pyth\pIIe@@pyth
+        \edef\pIIe@tempa{\strip@pt\@dashdim}%
+        \@tempdimc=\pIIe@tempa\@tempdimb
+        \global\divide\@tempdimc 8
+      \fi
+      \edef\x{\endgroup#3=\the\@tempdimc}%
+    \x}
+  \newcommand*\pIIe@@pyth{%
+    \pIIe@divide\@tempdima\@dashdim\@tempdimc
+    \advance\@dashdim\@tempdimc
+    \divide\@dashdim\tw@}
+  \newcommand*\pIIe@divide[3]{%
+    \begingroup
+    \dimendef\Numer=254\relax \dimendef\Denom=252\relax
+    \countdef\Num=254\relax   \countdef\Den=252\relax
+    \countdef\I=250\relax     \countdef\Numb=248\relax
+    \Numer #1\relax \Denom #2\relax
+    \ifdim\Denom<\z@ \Denom -\Denom \Numer=-\Numer \fi
+    \ifdim\Numer<\z@ \def\sign{-}\Numer=-\Numer \else \def\sign{}\fi
+    \ifdim\Denom=\z@
+      \edef\Q{\strip@pt\maxdimen}%
+      \PackageWarning{pict2e}%
+        {Division by 0, \sign\strip@pt\maxdimen\space used}{}%
+    \else
+      \Num=\Numer \Den=\Denom
+      \Numb=\Num \divide\Numb\Den
+      \ifnum\Numb>16383
+        \edef\Q{\strip@pt\maxdimen}%
+        \PackageWarning{pict2e}%
+          {Division overflow, \sign\strip@pt\maxdimen\space used}{}%
+      \else
+        \edef\Q{\number\Numb.}%
+        \multiply \Numb\Den \advance\Num -\Numb
+        \I=6\relax
+        \@whilenum \I>\z@ \do{\pIIe@@divide\advance\I\m@ne}%
+      \fi
+    \fi
+    \edef\tempend{\noexpand\endgroup\noexpand#3=\sign\Q\p@}%
+    \tempend}
+  \def\pIIe@@divide{%
+    \@whilenum \Num>214748364 \do{\divide\Num\tw@ \divide\Den\tw@}%
+    \multiply \Num 10
+    \Numb=\Num \divide\Numb\Den
+    \edef\Q{\Q\number\Numb}%
+    \multiply \Numb\Den \advance \Num -\Numb
+    \ifnum\Num>\z@\else\I=0\fi}
+  \newcommand*\pIIe@checkslopeargsline[2]{%
+    \pIIe@checkslopeargs{#1}{#2}{16383}}
+  \newcommand*\pIIe@checkslopeargsvector[2]{%
+    \pIIe@checkslopeargs{#1}{#2}{1000}}
+  \newcommand*\pIIe@checkslopeargs[3]{%
+    \edef\@tempa{#1}\expandafter\pIIe@checkslopearg\@tempa.:{#3}%
+    \edef\@tempa{#2}\expandafter\pIIe@checkslopearg\@tempa.:{#3}%
+    \ifdim #1\p@=\z@ \ifdim #2\p@=\z@ \@badlinearg \fi\fi}
+  \def\pIIe@checkslopearg #1.#2:#3{%
+    \def\@tempa{#1}%
+    \ifx\@tempa\empty\def\@tempa{0}\fi
+    \ifx\@tempa\space\def\@tempa{0}\fi
+    \ifnum\ifnum\@tempa<\z@-\fi\@tempa>#3 \@badlinearg \fi}
+  \def\@badlinearg{\PackageError
+    {pict2e}{Bad \protect\line\space or \protect\vector\space argument}{}}
+  \def\line(#1,#2)#3{%
+    \begingroup
+    \pIIe@checkslopeargsline{#1}{#2}%
+    \@tempdima=#1pt\relax \@tempdimb=#2pt\relax
+    \@defaultunitsset\@linelen{#3}\unitlength
+    \ifdim\@linelen<\z@ \@badlinearg \else
+      \pIIe@sline
+      \pIIe@moveto\z@\z@
+      \pIIe@lineto\@xdim\@ydim
+      \pIIe@strokeGraph
+      \box\@tempboxa
+    \fi
+    \endgroup}
+  \newcommand*\pIIe@sline{%
+    \ifdim\@tempdima=\z@
+      \ifdim\@tempdimb<\z@\@linelen-\@linelen\fi
+      \@ydim=\@linelen
+      \@xdim=\z@
+    \else
+      \ifdim\@tempdima<\z@\@linelen-\@linelen\fi
+      \ifdim\@tempdimb=\z@
+        \@xdim=\@linelen
+        \@ydim=\z@
+      \else
+        \pIIe@divide\@tempdimb\@tempdima\dimen@
+        \@ydim=\strip@pt\dimen@\@linelen
+        \@xdim=\@linelen
+      \fi
+    \fi
+    \@ovxx=\ifnum\@xdim=\z@ \z@\else\@linelen\fi
+    \@ovyy=\ifnum\@ydim<\z@ \z@\else\@ydim\fi
+    \@ovdy=\ifnum\@ydim<\z@ -\@ydim\else\z@\fi
+    \setbox\@tempboxa\hbox{%
+      \vrule\@height \@ovyy \@depth \@ovdy \@width \z@
+      \vrule\@height \z@ \@depth \z@ \@width \@ovxx}}
+  \def\vector(#1,#2)#3{%
+    \begingroup
+    \pIIe@checkslopeargsvector{#1}{#2}%
+    \@tempdima=#1pt\relax \@tempdimb=#2pt\relax
+    \@defaultunitsset\@linelen{#3}\unitlength
+    \ifdim\@linelen<\z@ \@badlinearg \else
+      \pIIe@sline
+      \@defaultunitsset\@linelen{#3}\unitlength
+      \pIIe@pyth{\@tempdima}{\@tempdimb}\dimen@
+      \ifdim\@tempdima=\z@ \else
+        \ifdim\@tempdimb=\z@ \else
+          \pIIe@divide\dimen@{\@tempdima}\@xdim
+          \@linelen\strip@pt\@xdim\@linelen
+          \ifdim\@linelen<\z@\@linelen-\@linelen\fi
+        \fi
+      \fi
+      \pIIe@divide{\@tempdimb}\dimen@\@ydim
+      \pIIe@divide{\@tempdima}\dimen@\@xdim
+      \pIIe@concat\@xdim\@ydim{-\@ydim}\@xdim\z@\z@
+      \pIIe@vector
+      \pIIe@fillGraph
+      \box\@tempboxa
+    \fi
+    \endgroup}
+  \newcommand*\pIIe@vector{}
+  \newcommand*\pIIe@vector@ltx{%
+    \@ydim\pIIe@FAW\@wholewidth \advance\@ydim\pIIe@CAW\relax
+    \@ovxx\pIIe@FAL\@ydim
+    \@xdim\@linelen \advance\@xdim-\@ovxx
+    \divide\@ydim\tw@
+    \divide\@ovxx\tw@ \advance\@ovxx\@xdim
+    \@ovyy\@ydim
+    \divide\@ovyy\tw@ \advance\@ovyy-\pIIe@FAI\@ydim
+    \pIIe@bezier@QtoC\@linelen\@ovxx\@ovro
+    \pIIe@bezier@QtoC\z@\@ovyy\@ovri
+    \pIIe@bezier@QtoC\@xdim\@ovxx\@clnwd
+    \pIIe@bezier@QtoC\@ydim\@ovyy\@clnht
+    \pIIe@moveto\@linelen\z@
+    \pIIe@curveto\@ovro{-\@ovri}\@clnwd{-\@clnht}\@xdim{-\@ydim}%
+    \ifdim\@xdim>\z@
+      \pIIe@lineto\@xdim{-\@halfwidth}%
+      \pIIe@lineto\z@{-\@halfwidth}%
+      \pIIe@lineto\z@{\@halfwidth}%
+      \pIIe@lineto\@xdim{\@halfwidth}%
+    \fi
+    \pIIe@lineto\@xdim\@ydim
+    \pIIe@curveto\@clnwd\@clnht\@ovro\@ovri\@linelen\z@}
+  \newcommand*\pIIe@vector@pst{%
+    \@ydim\pIIe@FAW\@wholewidth \advance\@ydim\pIIe@CAW\relax
+    \@ovxx\pIIe@FAL\@ydim
+    \@xdim\@linelen \advance\@xdim-\@ovxx
+    \divide\@ydim\tw@
+    \@ovyy\@ydim \advance\@ovyy-\@halfwidth
+    \@ovdx\pIIe@FAI\@ovxx
+    \pIIe@divide\@ovdx\@ydim\@tempdimc
+    \@ovxx\strip@pt\@ovyy\@tempdimc
+    \advance\@ovxx\@xdim
+    \advance\@ovdx\@xdim
+    \pIIe@moveto\@linelen\z@
+    \pIIe@lineto\@xdim{-\@ydim}%
+    \ifdim\@xdim>\z@
+      \pIIe@lineto\@ovxx{-\@halfwidth}%
+      \pIIe@lineto\z@{-\@halfwidth}%
+      \pIIe@lineto\z@{\@halfwidth}%
+      \pIIe@lineto\@ovxx{\@halfwidth}%
+    \else
+      \pIIe@lineto\@ovdx\z@
+    \fi
+    \pIIe@lineto\@xdim\@ydim
+    \pIIe@lineto\@linelen\z@}
+    \def\@circle#1{\begingroup \@tempswafalse\pIIe@circ{#1}}
+    \def\@dot#1{\begingroup \@tempswatrue\pIIe@circ{#1}}
+  \newcommand*\pIIe@circ[1]{%
+      \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+      \ifdim\pIIe@tempdima<\z@ \pIIe@badcircarg \fi
+      \divide\pIIe@tempdima\tw@
+      \pIIe@circle\pIIe@tempdima
+      \if@tempswa \pIIe@fillGraph \else \buttcap \pIIe@strokeGraph \fi
+    \endgroup}
+  \newcommand*\pIIe@circle[1]{%
+    \pIIe@qcircle[1]\z@{#1}\pIIe@qcircle  \@ne{#1}%
+    \pIIe@qcircle  \tw@{#1}\pIIe@qcircle\thr@@{#1}\pIIe@closepath}
+  \newcommand*\pIIe@qcircle[3][0]{%
+    \begingroup
+      \@ovro#3\relax \@ovri0.55228474983\@ovro
+      \@tempdimc\@ovri \advance\@tempdimc-\@ovro
+      \ifnum#1>\z@ \@tempswatrue \else \@tempswafalse \fi
+      \ifcase#2\relax
+        \pIIe@@qcircle\@ovro\z@\z@\@ovri\@tempdimc\@ovro{-\@ovro}\@ovro
+      \or
+        \pIIe@@qcircle\z@\@ovro{-\@ovri}\z@{-\@ovro}\@tempdimc{-\@ovro}{-\@ovro}%
+      \or
+        \pIIe@@qcircle{-\@ovro}\z@\z@{-\@ovri}{-\@tempdimc}{-\@ovro}\@ovro{-\@ovro}%
+      \or
+        \pIIe@@qcircle\z@{-\@ovro}\@ovri\z@\@ovro{-\@tempdimc}\@ovro\@ovro
+      \fi
+    \endgroup}
+  \newcommand*\pIIe@@qcircle[8]{%
+    \if@tempswa\pIIe@moveto{#1}{#2}\fi \pIIe@rcurveto{#3}{#4}{#5}{#6}{#7}{#8}}
+  \newcommand*\pIIe@badcircarg{%
+    \PackageError{pict2e}%
+      {Illegal argument in \protect\circle(*), \protect\oval, \protect\arc(*) or
+      \protect\circlearc.}%
+      {The radius of a circle, dot, arc or oval corner must be greater than zero.}}%
+  \newcommand*\maxovalrad{20pt}
+  \newcommand*\pIIe@defaultUL[2]{%
+    \@defaultunitsset\pIIe@tempdima{#2}\unitlength
+    \edef#1{\the\pIIe@tempdima}}
+  \newcommand*\pIIe@def@UL{}
+  \def\pIIe@def@UL#1\relax#2#3{%
+    \edef#2{\the\dimen@}}
+  \newcommand*\pIIe@maxovalrad{}
+  \newcommand*\pIIe@oval{}
+  \def\pIIe@oval#1(#2,#3){\@ifnextchar[{\@oval(#2,#3)}{\@oval(#2,#3)[]}}
+  \renewcommand*\oval[1][\maxovalrad]{%
+    \begingroup \pIIe@defaultUL\pIIe@maxovalrad{#1}%
+      \ifdim\pIIe@maxovalrad<\z@ \pIIe@badcircarg \fi
+      \pIIe@oval}
+  \def\@oval(#1,#2)[#3]{%
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength \divide\pIIe@tempdima\tw@
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength \divide\pIIe@tempdimb\tw@
+    \pIIe@tempdimc \ifdim\pIIe@tempdimb>\pIIe@tempdima \pIIe@tempdima \else \pIIe@tempdimb \fi
+    \ifdim\pIIe@maxovalrad<\pIIe@tempdimc \pIIe@tempdimc\pIIe@maxovalrad\relax \fi
+    \pIIe@tempdimd\pIIe@tempdima \advance\pIIe@tempdimd-\pIIe@tempdimc
+    \pIIe@tempdime\pIIe@tempdimb \advance\pIIe@tempdime-\pIIe@tempdimc
+    \pIIe@get@quadrants{#3}%
+    \ifnum15=\@tempcnta \pIIe@buttcap \fi
+    \@tempswatrue
+    \ifnum9=\@tempcnta
+      \pIIe@qoval\z@{-\pIIe@tempdimb}{\pIIe@tempdimd}{-\pIIe@tempdimb}%
+        \thr@@\pIIe@tempdimc\pIIe@tempdima\z@
+      \@tempcnta\@ne
+    \fi
+    \pIIe@qoval\pIIe@tempdima\z@\pIIe@tempdima\pIIe@tempdime%
+      \z@\pIIe@tempdimc\z@\pIIe@tempdimb
+    \pIIe@qoval\z@\pIIe@tempdimb{-\pIIe@tempdimd}\pIIe@tempdimb%
+      \@ne\pIIe@tempdimc{-\pIIe@tempdima}\z@
+    \pIIe@qoval{-\pIIe@tempdima}\z@{-\pIIe@tempdima}{-\pIIe@tempdime}%
+      \tw@\pIIe@tempdimc\z@{-\pIIe@tempdimb}%
+    \pIIe@qoval\z@{-\pIIe@tempdimb}{\pIIe@tempdimd}{-\pIIe@tempdimb}%
+      \thr@@\pIIe@tempdimc\pIIe@tempdima\z@
+    \pIIe@strokeGraph
+    \endgroup}
+  \newcommand*\pIIe@qoval[8]{%
+    \ifodd\@tempcnta
+      \if@tempswa\pIIe@moveto{#1}{#2}\fi
+      \pIIe@lineto{#3}{#4}\pIIe@qcircle{#5}{#6}\pIIe@lineto{#7}{#8}%
+      \@tempswafalse
+    \else
+      \@tempswatrue
+    \fi
+    \divide\@tempcnta\tw@}
+  \newcommand*\pIIe@get@quadrants[1]{%
+    \@ovttrue \@ovbtrue \@ovltrue \@ovrtrue \@tempcnta\z@
+    \@tfor\reserved@a:=#1\do{\csname @ov\reserved@a false\endcsname}%
+    \if@ovr \if@ovb\@setfpsbit2\fi \if@ovt\@setfpsbit4\fi \fi
+    \if@ovl \if@ovb\@setfpsbit1\fi \if@ovt\@setfpsbit8\fi \fi}
+  \def\@bezier#1(#2,#3)(#4,#5)(#6,#7){%
+    \ifnum #1=\z@
+      \@killglue
+      \begingroup
+      \@defaultunitsset\pIIe@tempdima{#2}\unitlength
+      \@defaultunitsset\pIIe@tempdimb{#3}\unitlength
+      \@defaultunitsset\pIIe@tempdimc{#4}\unitlength
+      \@defaultunitsset\pIIe@tempdimd{#5}\unitlength
+      \@defaultunitsset\pIIe@tempdime{#6}\unitlength
+      \@defaultunitsset\pIIe@tempdimf{#7}\unitlength
+      \pIIe@bezier@QtoC\pIIe@tempdima\pIIe@tempdimc\@ovro
+      \pIIe@bezier@QtoC\pIIe@tempdimb\pIIe@tempdimd\@ovri
+      \pIIe@bezier@QtoC\pIIe@tempdime\pIIe@tempdimc\@clnwd
+      \pIIe@bezier@QtoC\pIIe@tempdimf\pIIe@tempdimd\@clnht
+      \pIIe@moveto\pIIe@tempdima\pIIe@tempdimb
+      \pIIe@curveto\@ovro\@ovri\@clnwd\@clnht\pIIe@tempdime\pIIe@tempdimf
+      \pIIe@strokeGraph
+      \endgroup
+      \ignorespaces
+    \else
+      \pIIe@old@bezier{#1}(#2,#3)(#4,#5)(#6,#7)
+    \fi}
+  \newcommand*\pIIe@bezier@QtoC[3]{%
+    \@tempdimc#1\relax      \advance\@tempdimc-#2\relax
+    \divide\@tempdimc\thr@@ \advance\@tempdimc #2\relax
+    #3\@tempdimc}
+  \ifx\undefined\@arclen \newdimen\@arclen \fi
+  \ifx\undefined\@arcrad \newdimen\@arcrad \fi
+  \ifx\undefined\pIIe@tempdima \newdimen\pIIe@tempdima \fi
+  \ifx\undefined\pIIe@tempdimb \newdimen\pIIe@tempdimb \fi
+  \ifx\undefined\pIIe@tempdimc \newdimen\pIIe@tempdimc \fi
+  \ifx\undefined\pIIe@tempdimd \newdimen\pIIe@tempdimd \fi
+  \ifx\undefined\pIIe@tempdime \newdimen\pIIe@tempdime \fi
+  \ifx\undefined\pIIe@tempdimf \newdimen\pIIe@tempdimf \fi
+  \newcommand*\pIIe@arc[6][0]{%
+    \@arcrad #4\relax
+    \ifdim \@arcrad<\z@ \pIIe@badcircarg \else
+      \@arclen #6\p@ \advance\@arclen -#5\p@
+      \ifdim \@arclen<\z@ \def\sign{-}\else\def\sign{}\fi
+      \ifdim \sign\@arclen>720\p@
+        \PackageWarning {pict2e}{The arc angle is reduced to -720..720}%
+        \@whiledim \sign\@arclen>720\p@ \do {\advance\@arclen-\sign360\p@}%
+        \@tempdima #5\p@ \advance\@tempdima \@arclen
+        \edef\@angleend{\strip@pt\@tempdima}%
+        \pIIe@@arc{#1}{#2}{#3}{#4}{#5}{\@angleend}%
+      \else
+        \pIIe@@arc{#1}{#2}{#3}{#4}{#5}{#6}%
+      \fi
+    \fi}
+  \newcommand*\pIIe@@arc[6]{%
+    \begingroup
+    \ifdim \sign\@arclen>90\p@
+      \divide\@arclen 2
+      \@tempdima #5\p@ \advance\@tempdima \@arclen
+      \edef\@anglemid{\strip@pt\@tempdima}%
+      \def\@temp{\pIIe@@arc{#1}{#2}{#3}{#4}{#5}}%
+      \expandafter\@temp\expandafter{\@anglemid}%
+      \def\@temp{\pIIe@@arc{2}{#2}{#3}{#4}}%
+      \expandafter\@temp\expandafter{\@anglemid}{#6}%
+    \else
+      \CalculateSin{#5}\CalculateCos{#5}%
+      \@tempdima\UseCos{#5}\@arcrad \advance\@tempdima #2\relax
+      \@tempdimb\UseSin{#5}\@arcrad \advance\@tempdimb #3\relax
+      \ifcase #1\relax
+          \pIIe@lineto\@tempdima\@tempdimb
+      \or \pIIe@moveto\@tempdima\@tempdimb
+      \or
+      \else \PackageWarning {pict2e}%
+            {Illegal obligatory argument in \protect\circlearc.}%
+      \fi
+      \@tempdimc\@arclen \divide\@tempdimc\@iv
+      \edef\@angle{\strip@pt\@tempdimc}\CalculateTan{\@angle}%
+      \@linelen\UseTan{\@angle}\@arcrad \@linelen4\@linelen \divide\@linelen\thr@@
+      \advance\@tempdima-\UseSin{#5}\@linelen
+      \advance\@tempdimb \UseCos{#5}\@linelen
+      \pIIe@add@nums\@tempdima\@tempdimb
+      \CalculateSin{#6}\CalculateCos{#6}%
+      \@tempdima \UseCos{#6}\@arcrad \advance\@tempdima #2\relax
+      \@tempdimb \UseSin{#6}\@arcrad \advance\@tempdimb #3\relax
+      \@tempdimc \UseSin{#6}\@linelen \advance\@tempdimc \@tempdima
+      \pIIe@tempdimd-\UseCos{#6}\@linelen \advance\pIIe@tempdimd \@tempdimb
+      \pIIe@add@nums\@tempdimc\pIIe@tempdimd
+      \pIIe@add@CP\@tempdima\@tempdimb
+      \pIIe@addtoGraph\pIIe@curveto@op
+    \fi
+    \endgroup}
+  \newcommand*\pIIearc
+    {\@ifstar{\@tempswatrue\pIIe@arc@}{\@tempswafalse\pIIe@arc@}}
+  \newcommand*\pIIe@arc@[2][0,360]{\pIIe@arc@@(#1){#2}}
+  \def\pIIe@arc@@(#1,#2)#3{%
+    \@defaultunitsset\pIIe@tempdima{#3}\unitlength
+    \if@tempswa
+      \pIIe@moveto\z@\z@
+      \pIIe@arc{\z@}{\z@}{\pIIe@tempdima}{#1}{#2}%
+      \pIIe@closepath\pIIe@fillGraph
+    \else
+      \pIIe@arc[1]{\z@}{\z@}{\pIIe@tempdima}{#1}{#2}%
+      \pIIe@strokeGraph
+    \fi}
+  \ifx\undefined\arc
+  \else
+      \PackageWarning{pict2e}{\protect\arc\space redefined}%
+  \fi
+  \let\arc\pIIearc
+  \let\lp@r( \let\rp@r)
+  \def\Line(#1,#2)(#3,#4){\polyline(#1,#2)(#3,#4)}
+  \def\polyline(#1,#2){%
+    \@killglue
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \pIIe@moveto{\pIIe@tempdima}{\pIIe@tempdimb}%
+    \@ifnextchar\lp@r{\@polyline}{\PackageWarning{pict2e}%
+      {Polygonal lines require at least two vertices!}%
+    \ignorespaces}}
+  \def\@polyline(#1,#2){%
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \pIIe@lineto{\pIIe@tempdima}{\pIIe@tempdimb}%
+    \@ifnextchar\lp@r{\@polyline}{\pIIe@strokeGraph\ignorespaces}}
+  \def\Vector(#1,#2)(#3,#4){\polyvector(#1,#2)(#3,#4)}
+  \def\polyvector(#1,#2){%
+    \@killglue
+    \@ifnextchar\lp@r{\begingroup\@polyvector(#1,#2)}{%
+      \PackageWarning{pict2e}%
+      {Polygonal vectors require at least two vertices!}\ignorespaces}}
+  \def\@polyvector(#1,#2)(#3,#4){%
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \@defaultunitsset\pIIe@tempdimc{#3}\unitlength
+    \@defaultunitsset\pIIe@tempdimd{#4}\unitlength
+    \advance\pIIe@tempdimc-\pIIe@tempdima \advance\pIIe@tempdimd-\pIIe@tempdimb
+    \ifdim\pIIe@tempdimc=\z@ \@linelen\pIIe@tempdimd \else
+      \ifdim\pIIe@tempdimd=\z@ \@linelen\pIIe@tempdimc \else
+        \pIIe@pyth\pIIe@tempdimc\pIIe@tempdimd\@linelen
+      \fi
+    \fi
+    \ifdim\@linelen<\z@ \@linelen-\@linelen\fi
+    \pIIe@divide{\pIIe@tempdimc}\@linelen\pIIe@tempdime
+    \pIIe@divide{\pIIe@tempdimd}\@linelen\pIIe@tempdimf
+    \pIIe@concat\pIIe@tempdime\pIIe@tempdimf{-\pIIe@tempdimf}\pIIe@tempdime\pIIe@tempdima\pIIe@tempdimb
+    \pIIe@vector \pIIe@fillGraph
+    \@ifnextchar\lp@r{\@polyvector(#3,#4)}{\endgroup\ignorespaces}}
+  \def\polygon{%
+    \@killglue
+    \@ifstar{\begingroup\@tempswatrue\@polygon}%
+      {\begingroup\@tempswafalse\@polygon}}
+  \def\@polygon(#1,#2){%
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \pIIe@moveto{\pIIe@tempdima}{\pIIe@tempdimb}%
+    \@ifnextchar\lp@r{\@@polygon}{\PackageWarning{pict2e}%
+      {Polygons require at least two vertices!}%
+    \ignorespaces}}
+  \def\@@polygon(#1,#2){%
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \pIIe@lineto{\pIIe@tempdima}{\pIIe@tempdimb}%
+    \@ifnextchar\lp@r{\@@polygon}{\pIIe@closepath
+      \if@tempswa\pIIe@fillGraph\else\pIIe@strokeGraph\fi
+      \endgroup
+      \ignorespaces}}
+  \def\moveto(#1,#2){%
+    \@killglue
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \pIIe@moveto{\pIIe@tempdima}{\pIIe@tempdimb}%
+    \ignorespaces}
+  \def\lineto(#1,#2){%
+    \@killglue
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \pIIe@lineto{\pIIe@tempdima}{\pIIe@tempdimb}%
+    \ignorespaces}
+  \def\curveto(#1,#2)(#3,#4)(#5,#6){%
+    \@killglue
+    \@defaultunitsset\pIIe@tempdima{#1}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#2}\unitlength
+    \@defaultunitsset\pIIe@tempdimc{#3}\unitlength
+    \@defaultunitsset\pIIe@tempdimd{#4}\unitlength
+    \@defaultunitsset\pIIe@tempdime{#5}\unitlength
+    \@defaultunitsset\pIIe@tempdimf{#6}\unitlength
+    \pIIe@curveto{\pIIe@tempdima}{\pIIe@tempdimb}{\pIIe@tempdimc}%
+      {\pIIe@tempdimd}{\pIIe@tempdime}{\pIIe@tempdimf}%
+    \ignorespaces}
+  \newcommand*\circlearc[6][0]{%
+    \@killglue
+    \@defaultunitsset\pIIe@tempdima{#2}\unitlength
+    \@defaultunitsset\pIIe@tempdimb{#3}\unitlength
+    \@defaultunitsset\pIIe@tempdimc{#4}\unitlength
+    \pIIe@arc[#1]{\pIIe@tempdima}{\pIIe@tempdimb}{\pIIe@tempdimc}{#5}{#6}%
+    \ignorespaces}
+  \def\closepath{\pIIe@closepath}
+  \def\strokepath{\pIIe@strokeGraph}
+  \def\fillpath{\pIIe@fillGraph}
+  \ifcase\pIIe@mode\relax
+  \or
+    \newcommand*\pIIe@linecap@op{setlinecap}
+    \newcommand*\pIIe@linejoin@op{setlinejoin}
+  \or
+    \newcommand*\pIIe@linecap@op{J}
+    \newcommand*\pIIe@linejoin@op{j}
+  \fi
+  \def\pIIe@linecap{}
+  \def\pIIe@linejoin{}
+  \def\buttcap{\edef\pIIe@linecap{ 0 \pIIe@linecap@op}}
+  \def\roundcap{\edef\pIIe@linecap{ 1 \pIIe@linecap@op}}
+  \def\squarecap{\edef\pIIe@linecap{ 2 \pIIe@linecap@op}}
+  \def\miterjoin{\edef\pIIe@linejoin{ 0 \pIIe@linejoin@op}}
+  \def\roundjoin{\edef\pIIe@linejoin{ 1 \pIIe@linejoin@op}}
+  \def\beveljoin{\edef\pIIe@linejoin{ 2 \pIIe@linejoin@op}}
+  \AtBeginDocument{\@ifundefined{cbezier}{\newcommand}{\renewcommand}*%
+      \cbezier[2][0]{\pIIe@@cbezier[#1]#2}%
+    \@ifdefinable\pIIe@@cbezier{}%
+    \def\pIIe@@cbezier#1)#2(#3)#4(#5)#6({\@cbezier#1)(#3)(#5)(}%
+    \def\@cbezier[#1](#2,#3)(#4,#5)(#6,#7)(#8,#9){%
+      \@killglue
+      \@defaultunitsset\pIIe@tempdima{#2}\unitlength
+      \@defaultunitsset\pIIe@tempdimb{#3}\unitlength
+      \pIIe@moveto{\pIIe@tempdima}{\pIIe@tempdimb}%
+      \@defaultunitsset\pIIe@tempdima{#4}\unitlength
+      \@defaultunitsset\pIIe@tempdimb{#5}\unitlength
+      \@defaultunitsset\pIIe@tempdimc{#6}\unitlength
+      \@defaultunitsset\pIIe@tempdimd{#7}\unitlength
+      \@defaultunitsset\pIIe@tempdime{#8}\unitlength
+      \@defaultunitsset\pIIe@tempdimf{#9}\unitlength
+      \pIIe@curveto{\pIIe@tempdima}{\pIIe@tempdimb}{\pIIe@tempdimc}%
+        {\pIIe@tempdimd}{\pIIe@tempdime}{\pIIe@tempdimf}%
+      \pIIe@strokeGraph
+      \ignorespaces}%
+  }
+\else
+  \renewcommand*\oval[1][]{\pIIe@oldoval}
+  \newcommand*\maxovalrad{20pt}
+  \newcommand*\OriginalPictureCmds{}
+\fi
+\Gin@codes
+\let\Gin@codes\relax
+\endinput
+%%
+%% End of file `pict2e.sty'.

examples/PACKAGES/cgdna/util/generate.py

@@ -22,22 +22,26 @@
 """
 Import basic modules
 """
+
+# for python2/3 compatibility
+from __future__ import print_function
+
 import sys, os, timeit
 
 from timeit import default_timer as timer
 start_time = timer()
 """
-Try to import numpy; if failed, import a local version mynumpy 
+Try to import numpy; if failed, import a local version mynumpy
 which needs to be provided
 """
 try:
     import numpy as np
 except:
-    print >> sys.stderr, "numpy not found. Exiting."
+    print("numpy not found. Exiting.", file=sys.stderr)
     sys.exit(1)
 
 """
-Check that the required arguments (box offset and size in simulation units 
+Check that the required arguments (box offset and size in simulation units
 and the sequence file were provided
 """
 try:
@@ -45,8 +49,8 @@ try:
     box_length = float(sys.argv[2])
     infile = sys.argv[3]
 except:
-    print >> sys.stderr, "Usage: %s <%s> <%s> <%s>" % (sys.argv[0], \
-	"box offset", "box length", "file with sequences")
+    print( "Usage: %s <%s> <%s> <%s>" % (sys.argv[0], \
+	"box offset", "box length", "file with sequences"), file=sys.stderr)
     sys.exit(1)
 box = np.array ([box_length, box_length, box_length])
 
@@ -57,8 +61,7 @@ try:
     inp = open (infile, 'r')
     inp.close()
 except:
-    print >> sys.stderr, "Could not open file '%s' for reading. \
-					      Aborting." % infile
+    print( "Could not open file '%s' for reading. Aborting." % infile, file=sys.stderr)
     sys.exit(2)
 
 # return parts of a string
@@ -86,7 +89,7 @@ Define auxiliary variables for the construction of a helix
 # center of the double strand
 CM_CENTER_DS = POS_BASE + 0.2
 
-# ideal distance between base sites of two nucleotides 
+# ideal distance between base sites of two nucleotides
 # which are to be base paired in a duplex
 BASE_BASE = 0.3897628551303122
 
@@ -118,7 +121,7 @@ strandnum = []
 
 bonds = []
 
-""" 
+"""
 Convert local body frame to quaternion DOF
 """
 def exyz_to_quat (mya1, mya3):
@@ -135,25 +138,25 @@ def exyz_to_quat (mya1, mya3):
     # compute other components from it
 
     if q0sq >= 0.25:
-	myquat[0] = np.sqrt(q0sq)
-	myquat[1] = (mya2[2] - mya3[1]) / (4.0*myquat[0])
-	myquat[2] = (mya3[0] - mya1[2]) / (4.0*myquat[0])
-	myquat[3] = (mya1[1] - mya2[0]) / (4.0*myquat[0])
+        myquat[0] = np.sqrt(q0sq)
+        myquat[1] = (mya2[2] - mya3[1]) / (4.0*myquat[0])
+        myquat[2] = (mya3[0] - mya1[2]) / (4.0*myquat[0])
+        myquat[3] = (mya1[1] - mya2[0]) / (4.0*myquat[0])
     elif q1sq >= 0.25:
-	myquat[1] = np.sqrt(q1sq)
-	myquat[0] = (mya2[2] - mya3[1]) / (4.0*myquat[1])
-	myquat[2] = (mya2[0] + mya1[1]) / (4.0*myquat[1])
-	myquat[3] = (mya1[2] + mya3[0]) / (4.0*myquat[1])
+        myquat[1] = np.sqrt(q1sq)
+        myquat[0] = (mya2[2] - mya3[1]) / (4.0*myquat[1])
+        myquat[2] = (mya2[0] + mya1[1]) / (4.0*myquat[1])
+        myquat[3] = (mya1[2] + mya3[0]) / (4.0*myquat[1])
     elif q2sq >= 0.25:
-	myquat[2] = np.sqrt(q2sq)
-	myquat[0] = (mya3[0] - mya1[2]) / (4.0*myquat[2])
-	myquat[1] = (mya2[0] + mya1[1]) / (4.0*myquat[2])
-	myquat[3] = (mya3[1] + mya2[2]) / (4.0*myquat[2])
+        myquat[2] = np.sqrt(q2sq)
+        myquat[0] = (mya3[0] - mya1[2]) / (4.0*myquat[2])
+        myquat[1] = (mya2[0] + mya1[1]) / (4.0*myquat[2])
+        myquat[3] = (mya3[1] + mya2[2]) / (4.0*myquat[2])
     elif q3sq >= 0.25:
-	myquat[3] = np.sqrt(q3sq)
-	myquat[0] = (mya1[1] - mya2[0]) / (4.0*myquat[3])
-	myquat[1] = (mya3[0] + mya1[2]) / (4.0*myquat[3])
-	myquat[2] = (mya3[1] + mya2[2]) / (4.0*myquat[3])
+        myquat[3] = np.sqrt(q3sq)
+        myquat[0] = (mya1[1] - mya2[0]) / (4.0*myquat[3])
+        myquat[1] = (mya3[0] + mya1[2]) / (4.0*myquat[3])
+        myquat[2] = (mya3[1] + mya2[2]) / (4.0*myquat[3])
 
     norm = 1.0/np.sqrt(myquat[0]*myquat[0] + myquat[1]*myquat[1] + \
 			  myquat[2]*myquat[2] + myquat[3]*myquat[3])
@@ -169,62 +172,62 @@ Adds a strand to the system by appending it to the array of previous strands
 """
 def add_strands (mynewpositions, mynewa1s, mynewa3s):
     overlap = False
-	
-    # This is a simple check for each of the particles where for previously 
-    # placed particles i we check whether it overlaps with any of the 
+
+    # This is a simple check for each of the particles where for previously
+    # placed particles i we check whether it overlaps with any of the
     # newly created particles j
 
-    print >> sys.stdout, "## Checking for overlaps"
+    print( "## Checking for overlaps", file=sys.stdout)
 
-    for i in xrange(len(positions)):
+    for i in range(len(positions)):
 
-	p = positions[i]
-	pa1 = a1s[i]
+        p = positions[i]
+        pa1 = a1s[i]
 
-	for j in xrange (len(mynewpositions)):
+        for j in range (len(mynewpositions)):
 
-	    q = mynewpositions[j]
-	    qa1 = mynewa1s[j]
+            q = mynewpositions[j]
+            qa1 = mynewa1s[j]
 
-	    # skip particles that are anyway too far away
-	    dr = p - q
-	    dr -= box * np.rint (dr / box)
-	    if np.dot(dr, dr) > RC2:
-		continue
+            # skip particles that are anyway too far away
+            dr = p - q
+            dr -= box * np.rint(dr / box)
+            if np.dot(dr, dr) > RC2:
+                continue
 
-	    # base site and backbone site of the two particles
+            # base site and backbone site of the two particles
             p_pos_back = p + pa1 * POS_BACK
             p_pos_base = p + pa1 * POS_BASE
             q_pos_back = q + qa1 * POS_BACK
             q_pos_base = q + qa1 * POS_BASE
 
-	    # check for no overlap between the two backbone sites
+            # check for no overlap between the two backbone sites
             dr = p_pos_back - q_pos_back
-            dr -= box * np.rint (dr / box)
+            dr -= box * np.rint(dr / box)
             if np.dot(dr, dr) < RC2_BACK:
                 overlap = True
 
-	    # check for no overlap between the two base sites
+            # check for no overlap between the two base sites
             dr = p_pos_base -  q_pos_base
-            dr -= box * np.rint (dr / box)
+            dr -= box * np.rint(dr / box)
             if np.dot(dr, dr) < RC2_BASE:
                 overlap = True
 
-	    # check for no overlap between backbone site of particle p 
-	    # with base site of particle q
+            # check for no overlap between backbone site of particle p
+            # with base site of particle q
             dr = p_pos_back - q_pos_base
             dr -= box * np.rint (dr / box)
             if np.dot(dr, dr) < RC2_BACK_BASE:
                 overlap = True
 
-	    # check for no overlap between base site of particle p and 
-	    # backbone site of particle q
+            # check for no overlap between base site of particle p and
+            # backbone site of particle q
             dr = p_pos_base - q_pos_back
             dr -= box * np.rint (dr / box)
             if np.dot(dr, dr) < RC2_BACK_BASE:
                 overlap = True
 
-	    # exit if there is an overlap
+            # exit if there is an overlap
             if overlap:
                 return False
 
@@ -237,10 +240,10 @@ def add_strands (mynewpositions, mynewa1s, mynewa3s):
             a1s.append (p)
         for p in mynewa3s:
             a3s.append (p)
-	# calculate quaternion from local body frame and append
-	for ia in xrange(len(mynewpositions)):
-	    mynewquaternions = exyz_to_quat(mynewa1s[ia],mynewa3s[ia])
-	    quaternions.append(mynewquaternions)
+        # calculate quaternion from local body frame and append
+        for ia in range(len(mynewpositions)):
+            mynewquaternions = exyz_to_quat(mynewa1s[ia],mynewa3s[ia])
+            quaternions.append(mynewquaternions)
 
     return True
 
@@ -281,7 +284,7 @@ def get_rotation_matrix(axis, anglest):
                     [olc*x*z-st*y, olc*y*z+st*x, olc*z*z+ct]])
 
 """
-Generates the position and orientation vectors of a 
+Generates the position and orientation vectors of a
 (single or double) strand from a sequence string
 """
 def generate_strand(bp, sequence=None, start_pos=np.array([0, 0, 0]), \
@@ -295,76 +298,75 @@ def generate_strand(bp, sequence=None, start_pos=np.array([0, 0, 0]), \
     # overall direction of the helix
     dir = np.array(dir, dtype=float)
     if sequence == None:
-	sequence = np.random.randint(1, 5, bp)
+        sequence = np.random.randint(1, 5, bp)
 
-    # the elseif here is most likely redundant 
+    # the elseif here is most likely redundant
     elif len(sequence) != bp:
-	n = bp - len(sequence)
-	sequence += np.random.randint(1, 5, n)
-	print >> sys.stderr, "sequence is too short, adding %d random bases" % n
+        n = bp - len(sequence)
+        sequence += np.random.randint(1, 5, n)
+        print( "sequence is too short, adding %d random bases" % n, file=sys.stderr)
 
     # normalize direction
     dir_norm = np.sqrt(np.dot(dir,dir))
     if dir_norm < 1e-10:
-	print >> sys.stderr, "direction must be a valid vector, \
-			      defaulting to (0, 0, 1)"
-	dir = np.array([0, 0, 1])
+        print( "direction must be a valid vector, defaulting to (0, 0, 1)", file=sys.stderr)
+        dir = np.array([0, 0, 1])
     else: dir /= dir_norm
 
     # find a vector orthogonal to dir to act as helix direction,
     # if not provided switch off random orientation
     if perp is None or perp is False:
-	v1 = np.random.random_sample(3)
-	v1 -= dir * (np.dot(dir, v1))
-	v1 /= np.sqrt(sum(v1*v1))
+        v1 = np.random.random_sample(3)
+        v1 -= dir * (np.dot(dir, v1))
+        v1 /= np.sqrt(sum(v1*v1))
     else:
-	v1 = perp;
+        v1 = perp;
 
     # generate rotational matrix representing the overall rotation of the helix
     R0 = get_rotation_matrix(dir, rot)
-	    
+
     # rotation matrix corresponding to one step along the helix
     R = get_rotation_matrix(dir, [1, "bp"])
 
-    # set the vector a1 (backbone to base) to v1 
+    # set the vector a1 (backbone to base) to v1
     a1 = v1
-    
-    # apply the global rotation to a1 
+
+    # apply the global rotation to a1
     a1 = np.dot(R0, a1)
-    
+
     # set the position of the fist backbone site to start_pos
     rb = np.array(start_pos)
-	    
+
     # set a3 to the direction of the helix
     a3 = dir
     for i in range(bp):
     # work out the position of the centre of mass of the nucleotide
-	rcdm = rb - CM_CENTER_DS * a1
-	
-	# append to newpositions
-	mynewpositions.append(rcdm)
-	mynewa1s.append(a1)
-	mynewa3s.append(a3)
-	
-	# if we are not at the end of the helix, we work out a1 and rb for the 
-	# next nucleotide along the helix
-	if i != bp - 1:
-	    a1 = np.dot(R, a1)
-	    rb += a3 * BASE_BASE
+        rcdm = rb - CM_CENTER_DS * a1
 
-    # if we are working on a double strand, we do a cycle similar 
+        # append to newpositions
+        mynewpositions.append(rcdm)
+        mynewa1s.append(a1)
+        mynewa3s.append(a3)
+
+        # if we are not at the end of the helix, we work out a1 and rb for the
+        # next nucleotide along the helix
+        if i != bp - 1:
+            a1 = np.dot(R, a1)
+            rb += a3 * BASE_BASE
+
+    # if we are working on a double strand, we do a cycle similar
     # to the previous one but backwards
     if double == True:
-	a1 = -a1
-	a3 = -dir
-	R = R.transpose()
-	for i in range(bp):
-	    rcdm = rb - CM_CENTER_DS * a1
-	    mynewpositions.append (rcdm)
-	    mynewa1s.append (a1)
-	    mynewa3s.append (a3)
-	    a1 = np.dot(R, a1)
-	    rb += a3 * BASE_BASE
+        a1 = -a1
+        a3 = -dir
+        R = R.transpose()
+        for i in range(bp):
+            rcdm = rb - CM_CENTER_DS * a1
+            mynewpositions.append (rcdm)
+            mynewa1s.append (a1)
+            mynewa3s.append (a3)
+            a1 = np.dot(R, a1)
+            rb += a3 * BASE_BASE
 
     assert (len (mynewpositions) > 0)
 
@@ -391,10 +393,10 @@ def read_strands(filename):
     try:
         infile = open (filename)
     except:
-        print >> sys.stderr, "Could not open file '%s'. Aborting." % filename
+        print( "Could not open file '%s'. Aborting." % filename, file=sys.stderr )
         sys.exit(2)
 
-    # This block works out the number of nucleotides and strands by reading 
+    # This block works out the number of nucleotides and strands by reading
     # the number of non-empty lines in the input file and the number of letters,
     # taking the possible DOUBLE keyword into account.
     nstrands, nnucl, nbonds = 0, 0, 0
@@ -406,30 +408,29 @@ def read_strands(filename):
         if line[:6] == 'DOUBLE':
             line = line.split()[1]
             length = len(line)
-            print >> sys.stdout, "## Found duplex of %i base pairs" % length
+            print( "## Found duplex of %i base pairs" % length, file=sys.stdout)
             nnucl += 2*length
             nstrands += 2
-	    nbonds += (2*length-2)
+            nbonds += (2*length-2)
         else:
             line = line.split()[0]
             length = len(line)
-            print >> sys.stdout, \
-		    "## Found single strand of %i bases" % length
+            print( "## Found single strand of %i bases" % length, file=sys.stdout)
             nnucl += length
             nstrands += 1
-	    nbonds += length-1
+            nbonds += length-1
     # rewind the sequence input file
     infile.seek(0)
 
-    print >> sys.stdout, "## nstrands, nnucl = ", nstrands, nnucl
+    print( "## nstrands, nnucl = ", nstrands, nnucl, file=sys.stdout)
 
     # generate the data file in LAMMPS format
     try:
         out = open ("data.oxdna", "w")
     except:
-        print >> sys.stderr, "Could not open data file for writing. Aborting."
+        print( "Could not open data file for writing. Aborting.", file=sys.stderr)
         sys.exit(2)
-	
+
     lines = infile.readlines()
     nlines = len(lines)
     i = 1
@@ -440,115 +441,114 @@ def read_strands(filename):
         line = line.upper().strip()
 
         # skip empty lines
-        if len(line) == 0: 
-	    i += 1
-	    continue
+        if len(line) == 0:
+            i += 1
+            continue
 
-	# block for duplexes: last argument of the generate function 
-	# is set to 'True'
+        # block for duplexes: last argument of the generate function
+        # is set to 'True'
         if line[:6] == 'DOUBLE':
             line = line.split()[1]
             length = len(line)
             seq = [(base_to_number[x]) for x in line]
 
-	    myns += 1
-	    for b in xrange(length):
-		basetype.append(seq[b])
-		strandnum.append(myns)
+            myns += 1
+            for b in range(length):
+                basetype.append(seq[b])
+                strandnum.append(myns)
 
-	    for b in xrange(length-1):
-		bondpair = [noffset + b, noffset + b + 1]
-		bonds.append(bondpair)
-	    noffset += length
+            for b in range(length-1):
+                bondpair = [noffset + b, noffset + b + 1]
+                bonds.append(bondpair)
+            noffset += length
 
-	    # create the sequence of the second strand as made of 
-	    # complementary bases
-	    seq2 = [5-s for s in seq]
-	    seq2.reverse()
+            # create the sequence of the second strand as made of
+            # complementary bases
+            seq2 = [5-s for s in seq]
+            seq2.reverse()
 
-	    myns += 1
-	    for b in xrange(length):
-		basetype.append(seq2[b])
-		strandnum.append(myns)
+            myns += 1
+            for b in range(length):
+                basetype.append(seq2[b])
+                strandnum.append(myns)
 
-	    for b in xrange(length-1):
-		bondpair = [noffset + b, noffset + b + 1]
-		bonds.append(bondpair)
-	    noffset += length
- 
-            print >> sys.stdout, "## Created duplex of %i bases" % (2*length)
+            for b in range(length-1):
+                bondpair = [noffset + b, noffset + b + 1]
+                bonds.append(bondpair)
+            noffset += length
 
-	    # generate random position of the first nucleotide
+            print( "## Created duplex of %i bases" % (2*length), file=sys.stdout)
+
+            # generate random position of the first nucleotide
             cdm = box_offset + np.random.random_sample(3) * box
 
-            # generate the random direction of the helix 
+            # generate the random direction of the helix
             axis = np.random.random_sample(3)
             axis /= np.sqrt(np.dot(axis, axis))
 
-            # use the generate function defined above to create 
-	    # the position and orientation vector of the strand 
+            # use the generate function defined above to create
+            # the position and orientation vector of the strand
             newpositions, newa1s, newa3s = generate_strand(len(line), \
-		    sequence=seq, dir=axis, start_pos=cdm, double=True)
+                sequence=seq, dir=axis, start_pos=cdm, double=True)
 
             # generate a new position for the strand until it does not overlap
-	    # with anything already present
-	    start = timer()
+            # with anything already present
+            start = timer()
             while not add_strands(newpositions, newa1s, newa3s):
                 cdm = box_offset + np.random.random_sample(3) * box
                 axis = np.random.random_sample(3)
                 axis /= np.sqrt(np.dot(axis, axis))
                 newpositions, newa1s, newa3s = generate_strand(len(line), \
-		      sequence=seq, dir=axis, start_pos=cdm, double=True)
-                print >> sys.stdout, "## Trying %i" % i
-	    end = timer()
-            print >> sys.stdout, "## Added duplex of %i bases (line %i/%i) in %.2fs, now at %i/%i" % \
-				      (2*length, i, nlines, end-start, len(positions), nnucl)
+                    sequence=seq, dir=axis, start_pos=cdm, double=True)
+                print( "## Trying %i" % i, file=sys.stdout)
+            end = timer()
+            print( "## Added duplex of %i bases (line %i/%i) in %.2fs, now at %i/%i" % \
+				      (2*length, i, nlines, end-start, len(positions), nnucl), file=sys.stdout)
 
-	# block for single strands: last argument of the generate function 
-	# is set to 'False'
+        # block for single strands: last argument of the generate function
+        # is set to 'False'
         else:
             length = len(line)
             seq = [(base_to_number[x]) for x in line]
 
-	    myns += 1
-	    for b in xrange(length):
-		basetype.append(seq[b])
-		strandnum.append(myns)
+            myns += 1
+            for b in range(length):
+                basetype.append(seq[b])
+                strandnum.append(myns)
 
-	    for b in xrange(length-1):
-		bondpair = [noffset + b, noffset + b + 1]
-		bonds.append(bondpair)
-	    noffset += length
+            for b in range(length-1):
+                bondpair = [noffset + b, noffset + b + 1]
+                bonds.append(bondpair)
+            noffset += length
 
-	    # generate random position of the first nucleotide
+            # generate random position of the first nucleotide
             cdm = box_offset + np.random.random_sample(3) * box
 
-            # generate the random direction of the helix 
+            # generate the random direction of the helix
             axis = np.random.random_sample(3)
             axis /= np.sqrt(np.dot(axis, axis))
 
-            print >> sys.stdout, \
-		      "## Created single strand of %i bases" % length
+            print("## Created single strand of %i bases" % length, file=sys.stdout)
 
             newpositions, newa1s, newa3s = generate_strand(length, \
 		      sequence=seq, dir=axis, start_pos=cdm, double=False)
-	    start = timer()
+            start = timer()
             while not add_strands(newpositions, newa1s, newa3s):
                 cdm = box_offset + np.random.random_sample(3) * box
                 axis = np.random.random_sample(3)
-		axis /= np.sqrt(np.dot(axis, axis))
+                axis /= np.sqrt(np.dot(axis, axis))
                 newpositions, newa1s, newa3s = generate_strand(length, \
-			  sequence=seq, dir=axis, start_pos=cdm, double=False)
+                    sequence=seq, dir=axis, start_pos=cdm, double=False)
                 print >> sys.stdout, "## Trying  %i" % (i)
-	    end = timer()
-            print >> sys.stdout, "## Added single strand of %i bases (line %i/%i) in %.2fs, now at %i/%i" % \
-				      (length, i, nlines, end-start,len(positions), nnucl)
+            end = timer()
+            print( "## Added single strand of %i bases (line %i/%i) in %.2fs, now at %i/%i" % \
+				      (length, i, nlines, end-start,len(positions), nnucl), file=sys.stdout)
 
         i += 1
 
     # sanity check
     if not len(positions) == nnucl:
-        print len(positions), nnucl
+        print( len(positions), nnucl )
         raise AssertionError
 
     out.write('# LAMMPS data file\n')
@@ -580,44 +580,41 @@ def read_strands(filename):
     out.write('Atoms\n')
     out.write('\n')
 
-    for i in xrange(nnucl):
-	out.write('%d %d %22.15le %22.15le %22.15le %d 1 1\n' \
-		  % (i+1, basetype[i], \
-		     positions[i][0], positions[i][1], positions[i][2], \
-		     strandnum[i]))
+    for i in range(nnucl):
+        out.write('%d %d %22.15le %22.15le %22.15le %d 1 1\n' \
+            % (i+1, basetype[i], positions[i][0], positions[i][1], positions[i][2], strandnum[i]))
 
     out.write('\n')
     out.write('# Atom-ID, translational, rotational velocity\n')
     out.write('Velocities\n')
     out.write('\n')
 
-    for i in xrange(nnucl):
-	out.write("%d %22.15le %22.15le %22.15le %22.15le %22.15le %22.15le\n" \
-		  % (i+1,0.0,0.0,0.0,0.0,0.0,0.0))
+    for i in range(nnucl):
+        out.write("%d %22.15le %22.15le %22.15le %22.15le %22.15le %22.15le\n" \
+            % (i+1,0.0,0.0,0.0,0.0,0.0,0.0))
 
     out.write('\n')
     out.write('# Atom-ID, shape, quaternion\n')
     out.write('Ellipsoids\n')
     out.write('\n')
 
-    for i in xrange(nnucl):
-	out.write(\
-    "%d %22.15le %22.15le %22.15le %22.15le %22.15le %22.15le %22.15le\n"  \
-      % (i+1,1.1739845031423408,1.1739845031423408,1.1739845031423408, \
-	quaternions[i][0],quaternions[i][1], quaternions[i][2],quaternions[i][3]))
- 
+    for i in range(nnucl):
+        out.write("%d %22.15le %22.15le %22.15le %22.15le %22.15le %22.15le %22.15le\n"  \
+            % (i+1,1.1739845031423408,1.1739845031423408,1.1739845031423408, \
+            quaternions[i][0],quaternions[i][1], quaternions[i][2],quaternions[i][3]))
+
     out.write('\n')
     out.write('# Bond topology\n')
     out.write('Bonds\n')
     out.write('\n')
 
-    for i in xrange(nbonds):
-	out.write("%d  %d  %d  %d\n" % (i+1,1,bonds[i][0],bonds[i][1]))
+    for i in range(nbonds):
+        out.write("%d  %d  %d  %d\n" % (i+1,1,bonds[i][0],bonds[i][1]))
 
     out.close()
 
-    print >> sys.stdout, "## Wrote data to 'data.oxdna'"
-    print >> sys.stdout, "## DONE"
+    print("## Wrote data to 'data.oxdna'", file=sys.stdout)
+    print("## DONE", file=sys.stdout)
 
 # call the above main() function, which executes the program
 read_strands (infile)
@@ -627,4 +624,6 @@ runtime = end_time-start_time
 hours = runtime/3600
 minutes = (runtime-np.rint(hours)*3600)/60
 seconds = (runtime-np.rint(hours)*3600-np.rint(minutes)*60)%60
-print >> sys.stdout, "## Total runtime %ih:%im:%.2fs" % (hours,minutes,seconds)
+print( "## Total runtime %ih:%im:%.2fs" % (hours,minutes,seconds), file=sys.stdout)
+
+

examples/PACKAGES/cgdna/util/generate_simple.py

@@ -1,5 +1,8 @@
 # Setup tool for oxDNA input in LAMMPS format.
 
+# for python2/3 compatibility
+from __future__ import print_function
+
 import math,numpy as np,sys,os
 
 # system size
@@ -250,59 +253,59 @@ def duplex_array():
       qrot3=math.sin(0.5*twist)
 
       for letter in strand[2]:
-	temp1=[]
-	temp2=[]
+        temp1=[]
+        temp2=[]
 
-	temp1.append(nt2num[letter])
-	temp2.append(compnt2num[letter])
+        temp1.append(nt2num[letter])
+        temp2.append(compnt2num[letter])
 
-	temp1.append([posx1,posy1,posz1])
-	temp2.append([posx2,posy2,posz2])
+        temp1.append([posx1,posy1,posz1])
+        temp2.append([posx2,posy2,posz2])
 
-	vel=[0,0,0,0,0,0]
-	temp1.append(vel)
-	temp2.append(vel)
+        vel=[0,0,0,0,0,0]
+        temp1.append(vel)
+        temp2.append(vel)
 
-	temp1.append(shape)
-	temp2.append(shape)
+        temp1.append(shape)
+        temp2.append(shape)
 
-	temp1.append(quat1)
-	temp2.append(quat2)
+        temp1.append(quat1)
+        temp2.append(quat2)
 
-	quat1_0 = quat1[0]*qrot0 - quat1[1]*qrot1 - quat1[2]*qrot2 - quat1[3]*qrot3 
-	quat1_1 = quat1[0]*qrot1 + quat1[1]*qrot0 + quat1[2]*qrot3 - quat1[3]*qrot2 
-	quat1_2 = quat1[0]*qrot2 + quat1[2]*qrot0 + quat1[3]*qrot1 - quat1[1]*qrot3 
-	quat1_3 = quat1[0]*qrot3 + quat1[3]*qrot0 + quat1[1]*qrot2 + quat1[2]*qrot1 
+        quat1_0 = quat1[0]*qrot0 - quat1[1]*qrot1 - quat1[2]*qrot2 - quat1[3]*qrot3
+        quat1_1 = quat1[0]*qrot1 + quat1[1]*qrot0 + quat1[2]*qrot3 - quat1[3]*qrot2
+        quat1_2 = quat1[0]*qrot2 + quat1[2]*qrot0 + quat1[3]*qrot1 - quat1[1]*qrot3
+        quat1_3 = quat1[0]*qrot3 + quat1[3]*qrot0 + quat1[1]*qrot2 + quat1[2]*qrot1
 
-	quat1 = [quat1_0,quat1_1,quat1_2,quat1_3]
+        quat1 = [quat1_0,quat1_1,quat1_2,quat1_3]
 
-	posx1=axisx - dcomh*(quat1[0]**2+quat1[1]**2-quat1[2]**2-quat1[3]**2)
-	posy1=axisy - dcomh*(2*(quat1[1]*quat1[2]+quat1[0]*quat1[3]))
-	posz1=posz1+risez
+        posx1=axisx - dcomh*(quat1[0]**2+quat1[1]**2-quat1[2]**2-quat1[3]**2)
+        posy1=axisy - dcomh*(2*(quat1[1]*quat1[2]+quat1[0]*quat1[3]))
+        posz1=posz1+risez
 
-	quat2_0 = quat2[0]*qrot0 - quat2[1]*qrot1 - quat2[2]*qrot2 + quat2[3]*qrot3 
-	quat2_1 = quat2[0]*qrot1 + quat2[1]*qrot0 - quat2[2]*qrot3 - quat2[3]*qrot2 
-	quat2_2 = quat2[0]*qrot2 + quat2[2]*qrot0 + quat2[3]*qrot1 + quat2[1]*qrot3 
-	quat2_3 =-quat2[0]*qrot3 + quat2[3]*qrot0 + quat2[1]*qrot2 + quat2[2]*qrot1 
+        quat2_0 = quat2[0]*qrot0 - quat2[1]*qrot1 - quat2[2]*qrot2 + quat2[3]*qrot3
+        quat2_1 = quat2[0]*qrot1 + quat2[1]*qrot0 - quat2[2]*qrot3 - quat2[3]*qrot2
+        quat2_2 = quat2[0]*qrot2 + quat2[2]*qrot0 + quat2[3]*qrot1 + quat2[1]*qrot3
+        quat2_3 =-quat2[0]*qrot3 + quat2[3]*qrot0 + quat2[1]*qrot2 + quat2[2]*qrot1
 
-	quat2 = [quat2_0,quat2_1,quat2_2,quat2_3]
+        quat2 = [quat2_0,quat2_1,quat2_2,quat2_3]
 
-	posx2=axisx + dcomh*(quat1[0]**2+quat1[1]**2-quat1[2]**2-quat1[3]**2)
-	posy2=axisy + dcomh*(2*(quat1[1]*quat1[2]+quat1[0]*quat1[3]))
-	posz2=posz1
+        posx2=axisx + dcomh*(quat1[0]**2+quat1[1]**2-quat1[2]**2-quat1[3]**2)
+        posy2=axisy + dcomh*(2*(quat1[1]*quat1[2]+quat1[0]*quat1[3]))
+        posz2=posz1
 
-	if (len(nucleotide)+1 > strandstart):
-	  topology.append([1,len(nucleotide),len(nucleotide)+1])
-	  comptopo.append([1,len(nucleotide)+len(strand[2]),len(nucleotide)+len(strand[2])+1])
+        if (len(nucleotide)+1 > strandstart):
+          topology.append([1,len(nucleotide),len(nucleotide)+1])
+          comptopo.append([1,len(nucleotide)+len(strand[2]),len(nucleotide)+len(strand[2])+1])
 
-	nucleotide.append(temp1)
-	compstrand.append(temp2)
+        nucleotide.append(temp1)
+        compstrand.append(temp2)
 
       for ib in range(len(compstrand)):
-	nucleotide.append(compstrand[len(compstrand)-1-ib])
+        nucleotide.append(compstrand[len(compstrand)-1-ib])
 
       for ib in range(len(comptopo)):
-	topology.append(comptopo[ib])
+        topology.append(comptopo[ib])
 
   return
 

examples/PACKAGES/latboltz/toycar/in.toycar

@@ -457,9 +457,9 @@ mass *4 ${wheel_mass}
 variable dx equal 1.0
 variable density equal 0.001184
 
-neighbor 0.5 bin            
+neighbor 0.5 bin
 neigh_modify delay 0 every 1 check yes
-comm_modify cutoff 3.0     
+comm_modify cutoff 3.0
 
 pair_style lj/cut 1.2
 pair_coeff * * 0.0 0.0
@@ -474,7 +474,7 @@ reset_timestep 0
 #variable node_force equal "v_total_force / 178"
 #fix drag all addforce 0.0 0.0 0.0
 
-velocity all set 0.0 -7.5.0 0.0 units box
+velocity all set 0.0 -7.5 0.0 units box
 
 # viscosity of air is 0.0001847
 

examples/PACKAGES/pace/plugin/paceplugin.nsis

@@ -39,7 +39,7 @@ InstallDir "$LOCALAPPDATA\${PACEPLUGIN}"
 
 ShowInstDetails show
 ShowUninstDetails show
-SetCompressor lzma
+SetCompressor zlib
 
 !define MUI_ABORTWARNING
 

examples/PACKAGES/plumed/plugin/CMakeLists.txt

@@ -0,0 +1,59 @@
+# -*- CMake -*- build system for plugin examples.
+# The is meant to be used as a template for plugins that are
+# distributed independent from the LAMMPS package.
+##########################################
+
+cmake_minimum_required(VERSION 3.16)
+
+project(plumedplugin VERSION 1.0 LANGUAGES CXX)
+
+set(CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR})
+include(CheckIncludeFileCXX)
+include(LAMMPSInterfacePlugin)
+include(PLUMED)
+
+##########################
+# building the plugins
+
+add_library(plumedplugin MODULE plumedplugin.cpp ${LAMMPS_SOURCE_DIR}/PLUMED/fix_plumed.cpp)
+target_link_libraries(plumedplugin PRIVATE LAMMPS::PLUMED)
+target_link_libraries(plumedplugin PRIVATE lammps)
+target_include_directories(plumedplugin PRIVATE ${LAMMPS_SOURCE_DIR}/PLUMED)
+set_target_properties(plumedplugin PROPERTIES PREFIX "" SUFFIX ".so")
+
+# MacOS seems to need this
+if(CMAKE_SYSTEM_NAME STREQUAL Darwin)
+  set_target_properties(plumedplugin PROPERTIES LINK_FLAGS "-Wl,-undefined,dynamic_lookup")
+elseif(CMAKE_SYSTEM_NAME STREQUAL "Windows")
+# tell CMake to export all symbols to a .dll on Windows with special case for MinGW cross-compilers
+  set_target_properties(plumedplugin PROPERTIES WINDOWS_EXPORT_ALL_SYMBOLS TRUE)
+  if(CMAKE_CROSSCOMPILING)
+    set_target_properties(plumedplugin  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}/plumedplugin.nsis
+      ${CMAKE_BINARY_DIR})
+    if(BUILD_MPI)
+      if(USE_MSMPI)
+        add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION}-MSMPI plumedplugin.nsis
+          DEPENDS plumedplugin plumed_copy lammps.ico lammps-text-logo-wide.bmp plumedplugin.nsis
+          BYPRODUCTS LAMMPS-PLUMED-plugin-${LAMMPS_VERSION}-MSMPI.exe)
+      else()
+        add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION}-MPI plumedplugin.nsis
+          DEPENDS plumedplugin plumed_copy lammps.ico lammps-text-logo-wide.bmp plumedplugin.nsis
+          BYPRODUCTS LAMMPS-PLUMED-plugin-${LAMMPS_VERSION}-MPI.exe)
+      endif()
+    else()
+      add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION} plumedplugin.nsis
+        COMMAND ${CMAKE_COMMAND} -E echo ${PWD}
+        DEPENDS plumedplugin plumed_copy lammps.ico lammps-text-logo-wide.bmp plumedplugin.nsis
+        BYPRODUCTS LAMMPS-PLUMED-plugin-${LAMMPS_VERSION}.exe)
+    endif()
+  endif()
+else()
+  set_target_properties(plumedplugin PROPERTIES LINK_FLAGS "-rdynamic")
+endif()

examples/PACKAGES/plumed/plugin/LAMMPSInterfacePlugin.cmake

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

examples/PACKAGES/plumed/plugin/PLUMED.cmake

@@ -0,0 +1 @@
+../../../../cmake/Modules/Packages/PLUMED.cmake

examples/PACKAGES/plumed/plugin/README.txt

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

examples/PACKAGES/plumed/plugin/plumedplugin.cpp

@@ -0,0 +1,28 @@
+
+#include "lammpsplugin.h"
+#include "version.h"
+
+#include "fix_plumed.h"
+
+using namespace LAMMPS_NS;
+
+static Fix *fix_plumed_creator(LAMMPS *lmp, int argc, char **argv)
+{
+  return new FixPlumed(lmp, argc, argv);
+}
+
+extern "C" void lammpsplugin_init(void *lmp, void *handle, void *regfunc)
+{
+  lammpsplugin_t plugin;
+  lammpsplugin_regfunc register_plugin = (lammpsplugin_regfunc) regfunc;
+
+  // register plumed fix style
+  plugin.version = LAMMPS_VERSION;
+  plugin.style = "fix";
+  plugin.name = "plumed";
+  plugin.info = "Plumed2 plugin fix style v1.0";
+  plugin.author = "Axel Kohlmeyer (akohlmey@gmail.com)";
+  plugin.creator.v1 = (lammpsplugin_factory1 *) &fix_plumed_creator;
+  plugin.handle = handle;
+  (*register_plugin)(&plugin, lmp);
+}

examples/PACKAGES/plumed/plugin/plumedplugin.nsis

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

examples/meam/msmeam/HGa.meam

@@ -1,30 +0,0 @@
-bkgd_dyn        =       1
-emb_lin_neg = 1
-augt1=0 
-ialloy=1 
-rc	=	 5.9 
-#H
-attrac(1,1)=0.460 
-repuls(1,1)=0.460 
-Cmin(1,1,1)=1.3 # PuMS
-Cmax(1,1,1)= 2.80 
-nn2(1,1)=1
-#Ga
-rho0(2)         =       0.6
-attrac(2,2)=0.097 
-repuls(2,2)=0.097 
-nn2(2,2)=1
-#HGa
-attrac(1,2)=0.300 
-repuls(1,2)=0.300 
-lattce(1,2)=l12 
-re(1,2)=3.19 
-delta(1,2)=-0.48  
-alpha(1,2)=6.6 
-Cmin(1,1,2)=2.0 
-Cmin(2,1,2)= 2.0 
-Cmin(1,2,1)=2.0 
-Cmin(2,2,1)     =       1.4
-Cmin(1,2,2)     =       1.4
-Cmin(1,1,2)     =       1.4
-nn2(1,2)=1

examples/meam/msmeam/HGa.msmeam

@@ -0,0 +1 @@
+../../../potentials/HGa.msmeam

examples/meam/msmeam/data.msmeam.bu

@@ -1,25 +0,0 @@
-LAMMPS data file via write_data, version 16 Feb 2016, timestep = 1
-
-3 atoms
-2 atom types
-
--4.0000000000000000e+00 4.0000000000000000e+00 xlo xhi
--4.0000000000000000e+00 4.0000000000000000e+00 ylo yhi
--4.0000000000000000e+00 4.0000000000000000e+00 zlo zhi
-
-Masses
-
-1 1.0079
-2 69.723
-
-Atoms # atomic
-
-1 1 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0 0 0
-2 2 2.2000000000000002e+00 0.0000000000000000e+00 0.0000000000000000e+00 0 0 0
-3 2 2.9999999999999999e-01 2.2999999999999998e+00 0.0000000000000000e+00 0 0 0
-
-Velocities
-
-1 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00
-2 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00
-3 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00

examples/meam/msmeam/in.msmeam

@@ -1,5 +1,3 @@
-echo	both
-log	log.msmeam
 # Test of MEAM potential for HGa
 
 # ------------------------ INITIALIZATION ----------------------------
@@ -21,11 +19,11 @@ create_atoms    1 single 0 0 0  units box
 create_atoms    2 single 2.2 0 0  units box
 create_atoms    2 single 0.3 2.3 0  units box
 # ---------- Define Settings ---------------------
-variable	teng equal "c_eatoms"
+variable        teng equal "c_eatoms"
 compute pot_energy all pe/atom
 compute stress all stress/atom NULL
 dump 1 all custom 1 dump.msmeam id x y z fx fy fz c_pot_energy c_stress[1] c_stress[2] c_stress[3] c_stress[4] c_stress[5] c_stress[6]
-run	1
-write_data	data.msmeam
+run     1
+#write_data      data.msmeam
 
 print "All done!"

examples/meam/msmeam/library.msmeam

@@ -1,14 +0,0 @@
-# DATE: 2018-09-22 UNITS: metal CONTRIBUTOR: Steve Valone, smv@lanl.gov CITATION: Baskes, PRB 1992; smv, sr, mib, JNM 2010
-# ms-meam data format May 2010
-#  elt   lat    z     ielement      atwt
-#  alpha b0     b1    b2     b3     b1m     b2m   b3m      alat   esub   asub
-#    -   t0     t1    t2     t3     t1m     t2m   t3m      rozero ibar
-#  NOTE:  leading character cannot be a space
-
-'H'    'dim'  1.0   1      1.0079
-2.960  2.960  3.0   1.0    1.0    1.0    3.0  1.0       0.741  2.235  2.50
-1.0    0.44721 0.0  0.00   0.0    0.31623 0  6.70 0
-
-'Ga4'  'fcc'  12.0  31     69.723
-4.42   4.80   3.10  6.00   0.00   0.0    0.0  0.5       4.247  2.897  0.97
-1.0    1.649 1.435  0.00   0.0    0.0  2.0       0.70   0

examples/meam/msmeam/library.msmeam

@@ -0,0 +1 @@
+../../../potentials/library.msmeam

examples/meam/msmeam/log.1Mar2024.msmeam.g++.1

@@ -0,0 +1,126 @@
+LAMMPS (7 Feb 2024 - Development - patch_7Feb2024_update1-182-g93942f2013-modified)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (src/comm.cpp:98)
+  using 1 OpenMP thread(s) per MPI task
+# Test of MEAM potential for HGa
+
+# ------------------------ INITIALIZATION ----------------------------
+units           metal
+dimension       3
+boundary        p       p       p
+atom_style      atomic
+variable latparam equal 4.646
+variable ncell equal 3
+
+# ----------------------- ATOM DEFINITION ----------------------------
+region          box block -4 4 -4 4 -4 4
+create_box      2 box
+Created orthogonal box = (-4 -4 -4) to (4 4 4)
+  1 by 1 by 1 MPI processor grid
+
+#
+
+include potential.mod
+# NOTE: This script can be modified for different pair styles
+# See in.elastic for more info.
+
+variable Pu string H
+print "potential chosen ${Pu}"
+potential chosen H
+# Choose potential
+pair_style meam/ms
+print		"we just executed"
+we just executed
+
+pair_coeff      * * library.msmeam ${Pu} Ga4  HGa.msmeam ${Pu} Ga4
+pair_coeff      * * library.msmeam H Ga4  HGa.msmeam ${Pu} Ga4
+pair_coeff      * * library.msmeam H Ga4  HGa.msmeam H Ga4
+Reading MEAM library file library.msmeam with DATE: 2018-09-22
+# Setup neighbor style
+neighbor 1.0 bin
+neigh_modify once no every 1 delay 0 check yes
+
+# Setup minimization style
+variable dmax equal 1.0e-2
+min_style	     cg
+min_modify	     dmax ${dmax} line quadratic
+min_modify	     dmax 0.01 line quadratic
+compute eng all pe/atom
+compute eatoms all reduce sum c_eng
+
+# Setup output
+thermo		100
+thermo_style custom step temp etotal  press pxx pyy pzz pxy pxz pyz lx ly lz vol c_eatoms
+thermo_modify norm yes
+create_atoms    1 single 0 0 0  units box
+Created 1 atoms
+  using box units in orthogonal box = (-4 -4 -4) to (4 4 4)
+  create_atoms CPU = 0.000 seconds
+create_atoms    2 single 2.2 0 0  units box
+Created 1 atoms
+  using box units in orthogonal box = (-4 -4 -4) to (4 4 4)
+  create_atoms CPU = 0.000 seconds
+create_atoms    2 single 0.3 2.3 0  units box
+Created 1 atoms
+  using box units in orthogonal box = (-4 -4 -4) to (4 4 4)
+  create_atoms CPU = 0.000 seconds
+# ---------- Define Settings ---------------------
+variable        teng equal "c_eatoms"
+compute pot_energy all pe/atom
+compute stress all stress/atom NULL
+# dump 1 all custom 1 dump.msmeam id x y z fx fy fz c_pot_energy c_stress[1] c_stress[2] c_stress[3] c_stress[4] c_stress[5] c_stress[6]
+run     1
+WARNING: No fixes with time integration, atoms won't move (src/verlet.cpp:60)
+Neighbor list info ...
+  update: every = 1 steps, delay = 0 steps, check = yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 6.9
+  ghost atom cutoff = 6.9
+  binsize = 3.45, bins = 3 3 3
+  2 neighbor lists, perpetual/occasional/extra = 2 0 0
+  (1) pair meam/ms, perpetual
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+  (2) pair meam/ms, perpetual, half/full from (1)
+      attributes: half, newton on
+      pair build: halffull/newton
+      stencil: none
+      bin: none
+Per MPI rank memory allocation (min/avg/max) = 8.587 | 8.587 | 8.587 Mbytes
+   Step          Temp          TotEng         Press           Pxx            Pyy            Pzz            Pxy            Pxz            Pyz             Lx             Ly             Lz           Volume        c_eatoms   
+         0   0              15.433079      491354.7       838670.96      635393.15      0              80195.797      0              0              8              8              8              512            15.433079    
+         1   0              15.433079      491354.7       838670.96      635393.15      0              80195.797      0              0              8              8              8              512            15.433079    
+Loop time of 4.4446e-05 on 1 procs for 1 steps with 3 atoms
+
+Performance: 1943.932 ns/day, 0.012 hours/ns, 22499.213 timesteps/s, 67.498 katom-step/s
+31.5% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 2.9908e-05 | 2.9908e-05 | 2.9908e-05 |   0.0 | 67.29
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 1.033e-06  | 1.033e-06  | 1.033e-06  |   0.0 |  2.32
+Output  | 9.347e-06  | 9.347e-06  | 9.347e-06  |   0.0 | 21.03
+Modify  | 2.02e-07   | 2.02e-07   | 2.02e-07   |   0.0 |  0.45
+Other   |            | 3.956e-06  |            |       |  8.90
+
+Nlocal:              3 ave           3 max           3 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:             78 ave          78 max          78 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:              7 ave           7 max           7 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:           14 ave          14 max          14 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 14
+Ave neighs/atom = 4.6666667
+Neighbor list builds = 0
+Dangerous builds = 0
+#write_data      data.msmeam
+
+print "All done!"
+All done!
+Total wall time: 0:00:00

examples/meam/msmeam/log.1Mar2024.msmeam.g++.4

@@ -0,0 +1,126 @@
+LAMMPS (7 Feb 2024 - Development - patch_7Feb2024_update1-182-g93942f2013-modified)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (src/comm.cpp:98)
+  using 1 OpenMP thread(s) per MPI task
+# Test of MEAM potential for HGa
+
+# ------------------------ INITIALIZATION ----------------------------
+units           metal
+dimension       3
+boundary        p       p       p
+atom_style      atomic
+variable latparam equal 4.646
+variable ncell equal 3
+
+# ----------------------- ATOM DEFINITION ----------------------------
+region          box block -4 4 -4 4 -4 4
+create_box      2 box
+Created orthogonal box = (-4 -4 -4) to (4 4 4)
+  1 by 2 by 2 MPI processor grid
+
+#
+
+include potential.mod
+# NOTE: This script can be modified for different pair styles
+# See in.elastic for more info.
+
+variable Pu string H
+print "potential chosen ${Pu}"
+potential chosen H
+# Choose potential
+pair_style meam/ms
+print		"we just executed"
+we just executed
+
+pair_coeff      * * library.msmeam ${Pu} Ga4  HGa.msmeam ${Pu} Ga4
+pair_coeff      * * library.msmeam H Ga4  HGa.msmeam ${Pu} Ga4
+pair_coeff      * * library.msmeam H Ga4  HGa.msmeam H Ga4
+Reading MEAM library file library.msmeam with DATE: 2018-09-22
+# Setup neighbor style
+neighbor 1.0 bin
+neigh_modify once no every 1 delay 0 check yes
+
+# Setup minimization style
+variable dmax equal 1.0e-2
+min_style	     cg
+min_modify	     dmax ${dmax} line quadratic
+min_modify	     dmax 0.01 line quadratic
+compute eng all pe/atom
+compute eatoms all reduce sum c_eng
+
+# Setup output
+thermo		100
+thermo_style custom step temp etotal  press pxx pyy pzz pxy pxz pyz lx ly lz vol c_eatoms
+thermo_modify norm yes
+create_atoms    1 single 0 0 0  units box
+Created 1 atoms
+  using box units in orthogonal box = (-4 -4 -4) to (4 4 4)
+  create_atoms CPU = 0.000 seconds
+create_atoms    2 single 2.2 0 0  units box
+Created 1 atoms
+  using box units in orthogonal box = (-4 -4 -4) to (4 4 4)
+  create_atoms CPU = 0.000 seconds
+create_atoms    2 single 0.3 2.3 0  units box
+Created 1 atoms
+  using box units in orthogonal box = (-4 -4 -4) to (4 4 4)
+  create_atoms CPU = 0.000 seconds
+# ---------- Define Settings ---------------------
+variable        teng equal "c_eatoms"
+compute pot_energy all pe/atom
+compute stress all stress/atom NULL
+# dump 1 all custom 1 dump.msmeam id x y z fx fy fz c_pot_energy c_stress[1] c_stress[2] c_stress[3] c_stress[4] c_stress[5] c_stress[6]
+run     1
+WARNING: No fixes with time integration, atoms won't move (src/verlet.cpp:60)
+Neighbor list info ...
+  update: every = 1 steps, delay = 0 steps, check = yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 6.9
+  ghost atom cutoff = 6.9
+  binsize = 3.45, bins = 3 3 3
+  2 neighbor lists, perpetual/occasional/extra = 2 0 0
+  (1) pair meam/ms, perpetual
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+  (2) pair meam/ms, perpetual, half/full from (1)
+      attributes: half, newton on
+      pair build: halffull/newton
+      stencil: none
+      bin: none
+Per MPI rank memory allocation (min/avg/max) = 7.965 | 8.123 | 8.594 Mbytes
+   Step          Temp          TotEng         Press           Pxx            Pyy            Pzz            Pxy            Pxz            Pyz             Lx             Ly             Lz           Volume        c_eatoms   
+         0   0              15.433079      491354.7       838670.96      635393.15      0              80195.797      0              0              8              8              8              512            15.433079    
+         1   0              15.433079      491354.7       838670.96      635393.15      0              80195.797      0              0              8              8              8              512            15.433079    
+Loop time of 8.70645e-05 on 4 procs for 1 steps with 3 atoms
+
+Performance: 992.368 ns/day, 0.024 hours/ns, 11485.738 timesteps/s, 34.457 katom-step/s
+29.0% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 4.3957e-05 | 4.67e-05   | 5.1056e-05 |   0.0 | 53.64
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 1.105e-05  | 1.3822e-05 | 1.7033e-05 |   0.0 | 15.88
+Output  | 1.5765e-05 | 1.9045e-05 | 2.5216e-05 |   0.0 | 21.87
+Modify  | 2.58e-07   | 3.465e-07  | 3.81e-07   |   0.0 |  0.40
+Other   |            | 7.151e-06  |            |       |  8.21
+
+Nlocal:           0.75 ave           3 max           0 min
+Histogram: 3 0 0 0 0 0 0 0 0 1
+Nghost:          38.25 ave          42 max          36 min
+Histogram: 2 0 0 0 0 1 0 0 0 1
+Neighs:           1.75 ave           7 max           0 min
+Histogram: 3 0 0 0 0 0 0 0 0 1
+FullNghs:          3.5 ave          14 max           0 min
+Histogram: 3 0 0 0 0 0 0 0 0 1
+
+Total # of neighbors = 14
+Ave neighs/atom = 4.6666667
+Neighbor list builds = 0
+Dangerous builds = 0
+#write_data      data.msmeam
+
+print "All done!"
+All done!
+Total wall time: 0:00:00

examples/meam/msmeam/log.msmeam.bu

@@ -1,107 +0,0 @@
-# Test of MEAM potential for HGa
-
-# ------------------------ INITIALIZATION ----------------------------
-units           metal
-dimension       3
-boundary        p       p       p
-atom_style      atomic
-variable latparam equal 4.646
-variable ncell equal 3
-
-# ----------------------- ATOM DEFINITION ----------------------------
-region          box block -4 4 -4 4 -4 4
-create_box      2 box
-Created orthogonal box = (-4 -4 -4) to (4 4 4)
-  1 by 1 by 1 MPI processor grid
-
-#
-
-include potential.mod
-# NOTE: This script can be modified for different pair styles
-# See in.elastic for more info.
-
-variable Pu string H
-print "potential chosen ${Pu}"
-potential chosen H
-# Choose potential
-pair_style      MSmeam
-print		"we just executed"
-we just executed
-
-pair_coeff      * * library.MSmeam ${Pu} Ga4  HGaMS.meam ${Pu} Ga4
-pair_coeff      * * library.MSmeam H Ga4  HGaMS.meam ${Pu} Ga4
-pair_coeff      * * library.MSmeam H Ga4  HGaMS.meam H Ga4
-Reading potential file library.MSmeam with DATE: 2018-09-22
-# Setup neighbor style
-neighbor 1.0 nsq
-neigh_modify once no every 1 delay 0 check yes
-
-# Setup minimization style
-variable dmax equal 1.0e-2
-min_style	     cg
-min_modify	     dmax ${dmax} line quadratic
-min_modify	     dmax 0.01 line quadratic
-compute eng all pe/atom
-compute eatoms all reduce sum c_eng
-
-# Setup output
-thermo		100
-thermo_style custom step temp etotal  press pxx pyy pzz pxy pxz pyz lx ly lz vol c_eatoms
-thermo_modify norm yes
-create_atoms    1 single 0 0 0  units box
-Created 1 atoms
-create_atoms    2 single 2.2 0 0  units box
-Created 1 atoms
-create_atoms    2 single 0.3 2.3 0  units box
-Created 1 atoms
-# ---------- Define Settings ---------------------
-variable	teng equal "c_eatoms"
-compute pot_energy all pe/atom
-compute stress all stress/atom NULL
-dump 1 all custom 1 dump.msmeam id x y z fx fy fz c_pot_energy c_stress[1] c_stress[2] c_stress[3] c_stress[4] c_stress[5] c_stress[6]
-run	1
-WARNING: No fixes defined, atoms won't move (../verlet.cpp:55)
-Neighbor list info ...
-  2 neighbor list requests
-  update every 1 steps, delay 0 steps, check yes
-  max neighbors/atom: 2000, page size: 100000
-  master list distance cutoff = 6.9
-  ghost atom cutoff = 6.9
-Memory usage per processor = 12.9295 Mbytes
-Step Temp TotEng Press Pxx Pyy Pzz Pxy Pxz Pyz Lx Ly Lz Volume eatoms 
-       0            0    15.433079    491354.68    838670.91    635393.13            0    80195.793            0            0            8            8            8          512    15.433079 
-       1            0    15.433079    491354.68    838670.91    635393.13            0    80195.793            0            0            8            8            8          512    15.433079 
-Loop time of 0.000172138 on 1 procs for 1 steps with 3 atoms
-
-Performance: 501.922 ns/day, 0.048 hours/ns, 5809.285 timesteps/s
-81.3% CPU use with 1 MPI tasks x no OpenMP threads
-
-MPI task timing breakdown:
-Section |  min time  |  avg time  |  max time  |%varavg| %total
----------------------------------------------------------------
-Pair    | 6.6996e-05 | 6.6996e-05 | 6.6996e-05 |   0.0 | 38.92
-Neigh   | 0          | 0          | 0          |   0.0 |  0.00
-Comm    | 1.9073e-06 | 1.9073e-06 | 1.9073e-06 |   0.0 |  1.11
-Output  | 9.7036e-05 | 9.7036e-05 | 9.7036e-05 |   0.0 | 56.37
-Modify  | 0          | 0          | 0          |   0.0 |  0.00
-Other   |            | 6.199e-06  |            |       |  3.60
-
-Nlocal:    3 ave 3 max 3 min
-Histogram: 1 0 0 0 0 0 0 0 0 0
-Nghost:    78 ave 78 max 78 min
-Histogram: 1 0 0 0 0 0 0 0 0 0
-Neighs:    7 ave 7 max 7 min
-Histogram: 1 0 0 0 0 0 0 0 0 0
-FullNghs:  14 ave 14 max 14 min
-Histogram: 1 0 0 0 0 0 0 0 0 0
-
-Total # of neighbors = 14
-Ave neighs/atom = 4.66667
-Neighbor list builds = 0
-Dangerous builds = 0
-write_data	data.msmeam
-
-print "All done!"
-All done!
-Total wall time: 0:00:00
-

examples/meam/msmeam/msmeam.dump.bu

@@ -1,24 +0,0 @@
-ITEM: TIMESTEP
-0
-ITEM: NUMBER OF ATOMS
-3
-ITEM: BOX BOUNDS pp pp pp
--4 4
--4 4
--4 4
-ITEM: ATOMS id x y z fx fy fz c_pot_energy c_stress[1] c_stress[2] c_stress[3] c_stress[4] c_stress[5] c_stress[6] 
-1 0 0 0 -131.925 -88.3005 0 22.9153 -2.147e+08 -1.62661e+08 -0 -2.05301e+07 -0 -0 
-2 2.2 0 0 120.809 -0.482171 0 14.7692 -2.12028e+08 -0 -0 403352 -0 -0 
-3 0.3 2.3 0 11.1159 88.7827 0 8.61478 -2.67145e+06 -1.62661e+08 -0 -2.09335e+07 -0 -0 
-ITEM: TIMESTEP
-1
-ITEM: NUMBER OF ATOMS
-3
-ITEM: BOX BOUNDS pp pp pp
--4 4
--4 4
--4 4
-ITEM: ATOMS id x y z fx fy fz c_pot_energy c_stress[1] c_stress[2] c_stress[3] c_stress[4] c_stress[5] c_stress[6] 
-1 0 0 0 -131.925 -88.3005 0 22.9153 -2.147e+08 -1.62661e+08 -0 -2.05301e+07 -0 -0 
-2 2.2 0 0 120.809 -0.482171 0 14.7692 -2.12028e+08 -0 -0 403352 -0 -0 
-3 0.3 2.3 0 11.1159 88.7827 0 8.61478 -2.67145e+06 -1.62661e+08 -0 -2.09335e+07 -0 -0 

examples/meam/msmeam/potential.mod

@@ -7,7 +7,7 @@ print "potential chosen ${Pu}"
 pair_style meam/ms
 print		"we just executed"
 
-pair_coeff      * * library.msmeam ${Pu} Ga4  HGa.meam ${Pu} Ga4
+pair_coeff      * * library.msmeam ${Pu} Ga4  HGa.msmeam ${Pu} Ga4
 # Setup neighbor style
 neighbor 1.0 bin
 neigh_modify once no every 1 delay 0 check yes

fortran/lammps.f90

@@ -3633,7 +3633,7 @@ CONTAINS
 
     n = LEN_TRIM(f_string)
     ptr = lammps_malloc(n+1)
-    CALL C_F_POINTER(ptr, c_string, [1])
+    CALL C_F_POINTER(ptr, c_string, [n+1])
     DO i=1, n
         c_string(i) = f_string(i:i)
     END DO

lib/gpu/geryon/nvd_device.h

@@ -138,7 +138,7 @@ class UCL_Device {
   /** \note You cannot delete the default stream **/
   inline void pop_command_queue() {
     if (_cq.size()<2) return;
-    CU_SAFE_CALL_NS(cuStreamDestroy(_cq.back()));
+    cuStreamDestroy(_cq.back());
     _cq.pop_back();
   }
 
@@ -426,8 +426,8 @@ void UCL_Device::clear() {
   if (_device>-1) {
     for (int i=1; i<num_queues(); i++) pop_command_queue();
 #if GERYON_NVD_PRIMARY_CONTEXT
-    CU_SAFE_CALL_NS(cuCtxSetCurrent(_old_context));
-    CU_SAFE_CALL_NS(cuDevicePrimaryCtxRelease(_cu_device));
+    cuCtxSetCurrent(_old_context);
+    cuDevicePrimaryCtxRelease(_cu_device);
 #else
     cuCtxDestroy(_context);
 #endif

lib/gpu/lal_base_amoeba.h

@@ -33,6 +33,7 @@
 
 //#define ASYNC_DEVICE_COPY
 
+#if 0
 #if !defined(USE_OPENCL) && !defined(USE_HIP)
 // temporary workaround for int2 also defined in cufft
 #ifdef int2
@@ -40,6 +41,7 @@
 #endif
 #include "cufft.h"
 #endif
+#endif
 
 namespace LAMMPS_AL {
 
@@ -313,10 +315,11 @@ class BaseAmoeba {
   virtual int fphi_mpole();
   virtual int polar_real(const int eflag, const int vflag) = 0;
 
-
+#if 0
   #if !defined(USE_OPENCL) && !defined(USE_HIP)
   cufftHandle plan;
   #endif
+#endif
   bool fft_plan_created;
 };
 

lib/linalg/lmp_f2c.h

@@ -19,141 +19,13 @@ typedef long int logical;
 typedef short int shortlogical;
 typedef char logical1;
 typedef char integer1;
-#ifdef INTEGER_STAR_8   /* Adjust for integer*8. */
-typedef long long longint;              /* system-dependent */
-typedef unsigned long long ulongint;    /* system-dependent */
-#define qbit_clear(a,b) ((a) & ~((ulongint)1 << (b)))
-#define qbit_set(a,b)   ((a) |  ((ulongint)1 << (b)))
-#endif
 
 #define TRUE_ (1)
 #define FALSE_ (0)
 
-/* Extern is for use with -E */
-#ifndef Extern
-#define Extern extern
-#endif
-
-/* I/O stuff */
-
-#ifdef f2c_i2
-/* for -i2 */
-typedef short flag;
-typedef short ftnlen;
-typedef short ftnint;
-#else
-typedef long int flag;
 typedef long int ftnlen;
-typedef long int ftnint;
-#endif
-
-/*external read, write*/
-typedef struct
-{       flag cierr;
-        ftnint ciunit;
-        flag ciend;
-        char *cifmt;
-        ftnint cirec;
-} cilist;
-
-/*internal read, write*/
-typedef struct
-{       flag icierr;
-        char *iciunit;
-        flag iciend;
-        char *icifmt;
-        ftnint icirlen;
-        ftnint icirnum;
-} icilist;
-
-/*open*/
-typedef struct
-{       flag oerr;
-        ftnint ounit;
-        char *ofnm;
-        ftnlen ofnmlen;
-        char *osta;
-        char *oacc;
-        char *ofm;
-        ftnint orl;
-        char *oblnk;
-} olist;
-
-/*close*/
-typedef struct
-{       flag cerr;
-        ftnint cunit;
-        char *csta;
-} cllist;
-
-/*rewind, backspace, endfile*/
-typedef struct
-{       flag aerr;
-        ftnint aunit;
-} alist;
-
-/* inquire */
-typedef struct
-{       flag inerr;
-        ftnint inunit;
-        char *infile;
-        ftnlen infilen;
-        ftnint	*inex;	/*parameters in standard's order*/
-        ftnint	*inopen;
-        ftnint	*innum;
-        ftnint	*innamed;
-        char	*inname;
-        ftnlen	innamlen;
-        char	*inacc;
-        ftnlen	inacclen;
-        char	*inseq;
-        ftnlen	inseqlen;
-        char 	*indir;
-        ftnlen	indirlen;
-        char	*infmt;
-        ftnlen	infmtlen;
-        char	*inform;
-        ftnint	informlen;
-        char	*inunf;
-        ftnlen	inunflen;
-        ftnint	*inrecl;
-        ftnint	*innrec;
-        char	*inblank;
-        ftnlen	inblanklen;
-} inlist;
-
 #define VOID void
 
-union Multitype {       /* for multiple entry points */
-        integer1 g;
-        shortint h;
-        integer i;
-        /* longint j; */
-        real r;
-        doublereal d;
-        complex c;
-        doublecomplex z;
-        };
-
-typedef union Multitype Multitype;
-
-/*typedef long int Long;*/      /* No longer used; formerly in Namelist */
-
-struct Vardesc {        /* for Namelist */
-        char *name;
-        char *addr;
-        ftnlen *dims;
-        int  type;
-        };
-typedef struct Vardesc Vardesc;
-
-struct Namelist {
-        char *name;
-        Vardesc **vars;
-        int nvars;
-        };
-typedef struct Namelist Namelist;
-
 #define abs(x) ((x) >= 0 ? (x) : -(x))
 #define dabs(x) (doublereal)abs(x)
 #define min(a,b) ((a) <= (b) ? (a) : (b))

lib/linalg/s_lmp_cat.cpp

@@ -4,7 +4,7 @@
 // concatenate two strings
 
 extern "C" {
-void s_lmp_cat(char *lp, char *rpp[], ftnint rnp[], ftnint *np, ftnlen ll)
+void s_lmp_cat(char *lp, char **rpp, integer *rnp, integer *np, ftnlen ll)
 {
     ftnlen i, nc;
     char *rp;

lib/mdi/Install.py

@@ -32,7 +32,7 @@ make lib-mdi args="-m mpi" # build MDI lib with same settings as in the mpi Make
 
 # settings
 
-version = "1.4.16"
+version = "1.4.26"
 url = "https://github.com/MolSSI-MDI/MDI_Library/archive/v%s.tar.gz" % version
 
 # known checksums for different MDI versions. used to validate the download.
@@ -41,6 +41,7 @@ checksums = { \
               '1.4.12' : '7a222353ae8e03961d5365e6cd48baee', \
               '1.4.14' : '7a059bb12535360fdcb7de8402f9a0fc', \
               '1.4.16' : '407db44e2d79447ab5c1233af1965f65', \
+              '1.4.26' : '3124bb85259471e2a53a891f04bf697a', \
               }
 
 # print error message or help

lib/plumed/Install.py

@@ -16,8 +16,10 @@ parser = ArgumentParser(prog='Install.py',
                         description="LAMMPS library build wrapper script")
 
 # settings
+# Note: must also adjust check for supported API versions in
+# fix_plumed.cpp when version changes from v2.n.x to v2.n+1.y
 
-version = "2.8.2"
+version = "2.9.1"
 mode = "static"
 
 # help message
@@ -45,6 +47,10 @@ checksums = { \
         '2.7.6' : 'fb8c0ec10f97a9353eb123a5c4c35aa6', \
         '2.8.1' : '6bfe72ebdae63dc38a9ca27d9b0e08f8', \
         '2.8.2' : '599092b6a0aa6fff992612537ad98994', \
+        '2.8.3' : '76d23cd394eba9e6530316ed1184e219', \
+        '2.8.4' : '9f59c4f9bda86fe5bef19543c295a981', \
+        '2.9.0' : '661eabeebee05cf84bbf9dc23d7d5f46', \
+        '2.9.1' : 'c3b2d31479c1e9ce211719d40e9efbd7', \
         }
 
 # parse and process arguments

potentials/HGa.msmeam

@@ -1,29 +1,29 @@
 bkgd_dyn        =       1
 emb_lin_neg = 1
-augt1=0 
-ialloy=1 
-rc	=	 5.9 
+augt1=0
+ialloy=1
+rc      =        5.9
 #H
-attrac(1,1)=0.460 
-repuls(1,1)=0.460 
+attrac(1,1)=0.460
+repuls(1,1)=0.460
 Cmin(1,1,1)=1.3 # PuMS
-Cmax(1,1,1)= 2.80 
+Cmax(1,1,1)= 2.80
 nn2(1,1)=1
 #Ga
 rho0(2)         =       0.6
-attrac(2,2)=0.097 
-repuls(2,2)=0.097 
+attrac(2,2)=0.097
+repuls(2,2)=0.097
 nn2(2,2)=1
 #HGa
-attrac(1,2)=0.300 
-repuls(1,2)=0.300 
-lattce(1,2)=l12 
-re(1,2)=3.19 
-delta(1,2)=-0.48  
-alpha(1,2)=6.6 
-Cmin(1,1,2)=2.0 
-Cmin(2,1,2)= 2.0 
-Cmin(1,2,1)=2.0 
+attrac(1,2)=0.300
+repuls(1,2)=0.300
+lattce(1,2)=l12
+re(1,2)=3.19
+delta(1,2)=-0.48
+alpha(1,2)=6.6
+Cmin(1,1,2)=2.0
+Cmin(2,1,2)= 2.0
+Cmin(1,2,1)=2.0
 Cmin(2,2,1)     =       1.4
 Cmin(1,2,2)     =       1.4
 Cmin(1,1,2)     =       1.4

python/lammps/core.py

@@ -334,8 +334,8 @@ class lammps(object):
     if self.has_mpi_support:
       try:
         from mpi4py import __version__ as mpi4py_version
-        # tested to work with mpi4py versions 2 and 3
-        self.has_mpi4py = mpi4py_version.split('.')[0] in ['2','3']
+        # tested to work with mpi4py versions 2, 3, and 4
+        self.has_mpi4py = mpi4py_version.split('.')[0] in ['2','3','4']
       except ImportError:
         # ignore failing import
         pass
@@ -361,7 +361,7 @@ class lammps(object):
         if not self.has_mpi_support:
           raise Exception('LAMMPS not compiled with real MPI library')
         if not self.has_mpi4py:
-          raise Exception('Python mpi4py version is not 2 or 3')
+          raise Exception('Python mpi4py version is not 2, 3, or 4')
         if self.MPI._sizeof(self.MPI.Comm) == sizeof(c_int):
           MPI_Comm = c_int
         else:
@@ -374,12 +374,16 @@ class lammps(object):
         narg = 0
         cargs = None
         if cmdargs is not None:
-          cmdargs.insert(0,"lammps")
-          narg = len(cmdargs)
-          for i in range(narg):
-            if type(cmdargs[i]) is str:
-              cmdargs[i] = cmdargs[i].encode()
-          cargs = (c_char_p*(narg+1))(*cmdargs)
+          myargs = ["lammps".encode()]
+          narg = len(cmdargs) + 1
+          for arg in cmdargs:
+            if type(arg) is str:
+              myargs.append(arg.encode())
+            elif type(arg) is bytes:
+              myargs.append(arg)
+            else:
+              raise TypeError('Unsupported cmdargs type ', type(arg))
+          cargs = (c_char_p*(narg+1))(*myargs)
           cargs[narg] = None
           self.lib.lammps_open.argtypes = [c_int, c_char_p*(narg+1), MPI_Comm, c_void_p]
         else:
@@ -395,12 +399,16 @@ class lammps(object):
           self.comm = self.MPI.COMM_WORLD
         self.opened = 1
         if cmdargs is not None:
-          cmdargs.insert(0,"lammps")
-          narg = len(cmdargs)
-          for i in range(narg):
-            if type(cmdargs[i]) is str:
-              cmdargs[i] = cmdargs[i].encode()
-          cargs = (c_char_p*(narg+1))(*cmdargs)
+          myargs = ["lammps".encode()]
+          narg = len(cmdargs) + 1
+          for arg in cmdargs:
+            if type(arg) is str:
+              myargs.append(arg.encode())
+            elif type(arg) is bytes:
+              myargs.append(arg)
+            else:
+              raise TypeError('Unsupported cmdargs type ', type(arg))
+          cargs = (c_char_p*(narg+1))(*myargs)
           cargs[narg] = None
           self.lib.lammps_open_no_mpi.argtypes = [c_int, c_char_p*(narg+1), c_void_p]
           self.lmp = c_void_p(self.lib.lammps_open_no_mpi(narg,cargs,None))

python/lammps/pylammps.py

@@ -192,11 +192,23 @@ class Atom(object):
   @property
   def mass(self):
     """
-    Return the atom mass
+    Return the atom mass as a per-atom property.
+    This returns either the per-type mass or the per-atom
+    mass (AKA 'rmass') depending on what is available with
+    preference for the per-atom mass.
+
+    .. versionchanged:: 2Aug2023_update_4
+
+       Support both per-type and per-atom masses. With
+       per-type return "mass[type[i]]" else return "rmass[i]".
+       Per-atom mass is preferred if available.
 
     :type: float
     """
-    return self.get("mass", self.index)
+    if self._pylmp.lmp.extract_setting('rmass_flag'):
+      return self.get("rmass", self.index)
+    else:
+      return self.get("mass", self.type)
 
   @property
   def radius(self):

src/ASPHERE/pair_ylz.cpp

@@ -300,6 +300,7 @@ double PairYLZ::init_one(int i, int j)
   zeta[j][i] = zeta[i][j];
   mu[j][i] = mu[i][j];
   beta[j][i] = beta[i][j];
+  cut[j][i] = cut[i][j];
 
   return cut[i][j];
 }
@@ -409,7 +410,7 @@ void PairYLZ::write_data_all(FILE *fp)
 {
   for (int i = 1; i <= atom->ntypes; i++)
     for (int j = i; j <= atom->ntypes; j++)
-      fprintf(fp, "%d %d %g %g %g %g %g %g\n", i, j, epsilon[i][i], sigma[i][i], cut[i][j],
+      fprintf(fp, "%d %d %g %g %g %g %g %g\n", i, j, epsilon[i][j], sigma[i][j], cut[i][j],
               zeta[i][j], mu[i][j], beta[i][j]);
 }
 

src/BODY/body_rounded_polyhedron.cpp

@@ -99,10 +99,9 @@ int BodyRoundedPolyhedron::nedges(AtomVecBody::Bonus *bonus)
 {
   int nvertices = bonus->ivalue[0];
   int nedges = bonus->ivalue[1];
-  //int nfaces = bonus->ivalue[2];
   if (nvertices == 1) return 0;
   else if (nvertices == 2) return 1;
-  return nedges; //(nvertices+nfaces-2); // Euler formula: V-E+F=2
+  return nedges;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -116,6 +115,9 @@ double *BodyRoundedPolyhedron::edges(AtomVecBody::Bonus *bonus)
 
 int BodyRoundedPolyhedron::nfaces(AtomVecBody::Bonus *bonus)
 {
+  int nvertices = bonus->ivalue[0];
+  if (nvertices < 3) return 0;
+
   return bonus->ivalue[2];
 }
 

src/COLLOID/pair_brownian.cpp

@@ -480,15 +480,18 @@ void PairBrownian::init_style()
   // are re-calculated at every step.
 
   flagdeform = flagwall = 0;
-  for (int i = 0; i < modify->nfix; i++) {
-    if (strcmp(modify->fix[i]->style, "deform") == 0)
-      flagdeform = 1;
-    else if (strstr(modify->fix[i]->style, "wall") != nullptr) {
-      if (flagwall) error->all(FLERR, "Cannot use multiple fix wall commands with pair brownian");
-      flagwall = 1;    // Walls exist
-      wallfix = dynamic_cast<FixWall *>(modify->fix[i]);
-      if (wallfix->xflag) flagwall = 2;    // Moving walls exist
-    }
+  wallfix = nullptr;
+
+  if (modify->get_fix_by_style("^deform").size() > 0) flagdeform = 1;
+  auto fixes = modify->get_fix_by_style("^wall");
+  if (fixes.size() > 1)
+    error->all(FLERR, "Cannot use multiple fix wall commands with pair brownian");
+  else if (fixes.size() == 1) {
+    wallfix = dynamic_cast<FixWall *>(fixes[0]);
+    if (!wallfix)
+      error->all(FLERR, "Fix {} is not compatible with pair brownian", fixes[0]->style);
+    flagwall = 1;
+    if (wallfix->xflag) flagwall = 2; // Moving walls exist
   }
 
   // set the isotropic constants depending on the volume fraction

src/COLLOID/pair_brownian_poly.cpp

@@ -349,17 +349,18 @@ void PairBrownianPoly::init_style()
   // are re-calculated at every step.
 
   flagdeform = flagwall = 0;
-  for (int i = 0; i < modify->nfix; i++) {
-    if (strcmp(modify->fix[i]->style,"deform") == 0)
-      flagdeform = 1;
-    else if (strstr(modify->fix[i]->style,"wall") != nullptr) {
-      if (flagwall)
-        error->all(FLERR,
-                   "Cannot use multiple fix wall commands with pair brownian");
-      flagwall = 1; // Walls exist
-      wallfix = dynamic_cast<FixWall *>(modify->fix[i]);
-      if (wallfix->xflag) flagwall = 2; // Moving walls exist
-    }
+  wallfix = nullptr;
+
+  if (modify->get_fix_by_style("^deform").size() > 0) flagdeform = 1;
+  auto fixes = modify->get_fix_by_style("^wall");
+  if (fixes.size() > 1)
+    error->all(FLERR, "Cannot use multiple fix wall commands with pair brownian/poly");
+  else if (fixes.size() == 1) {
+    wallfix = dynamic_cast<FixWall *>(fixes[0]);
+    if (!wallfix)
+      error->all(FLERR, "Fix {} is not compatible with pair brownian/poly", fixes[0]->style);
+    flagwall = 1;
+    if (wallfix->xflag) flagwall = 2; // Moving walls exist
   }
 
   // set the isotropic constants that depend on the volume fraction

src/COLLOID/pair_lubricate.cpp

@@ -560,21 +560,23 @@ void PairLubricate::init_style()
   // are re-calculated at every step.
 
   shearing = flagdeform = flagwall = 0;
-  for (int i = 0; i < modify->nfix; i++) {
-    if (strcmp(modify->fix[i]->style,"deform") == 0) {
-      shearing = flagdeform = 1;
-      if ((dynamic_cast<FixDeform *>(modify->fix[i]))->remapflag != Domain::V_REMAP)
-        error->all(FLERR,"Using pair lubricate with inconsistent "
-                   "fix deform remap option");
-    }
-    if (strstr(modify->fix[i]->style,"wall") != nullptr) {
-      if (flagwall)
-        error->all(FLERR,
-                   "Cannot use multiple fix wall commands with pair lubricate");
-      flagwall = 1; // Walls exist
-      wallfix = dynamic_cast<FixWall *>(modify->fix[i]);
-      if (wallfix->xflag) flagwall = 2; // Moving walls exist
-    }
+
+  auto fixes = modify->get_fix_by_style("^deform");
+  if (fixes.size() > 0) {
+    shearing = flagdeform = 1;
+    auto *myfix = dynamic_cast<FixDeform *>(fixes[0]);
+    if (myfix && (myfix->remapflag != Domain::V_REMAP))
+      error->all(FLERR,"Using pair lubricate with inconsistent fix deform remap option");
+  }
+  fixes = modify->get_fix_by_style("^wall");
+  if (fixes.size() > 1)
+    error->all(FLERR, "Cannot use multiple fix wall commands with pair lubricate");
+  else if (fixes.size() == 1) {
+    wallfix = dynamic_cast<FixWall *>(fixes[0]);
+    if (!wallfix)
+      error->all(FLERR, "Fix {} is not compatible with pair lubricate", fixes[0]->style);
+    flagwall = 1;
+    if (wallfix->xflag) flagwall = 2; // Moving walls exist
   }
 
   // set the isotropic constants that depend on the volume fraction

src/COLLOID/pair_lubricateU.cpp

@@ -1792,18 +1792,18 @@ void PairLubricateU::init_style()
   // are re-calculated at every step.
 
   flagdeform = flagwall = 0;
-  for (int i = 0; i < modify->nfix; i++) {
-    if (strcmp(modify->fix[i]->style,"deform") == 0)
-      flagdeform = 1;
-    else if (strstr(modify->fix[i]->style,"wall") != nullptr) {
-      if (flagwall)
-        error->all(FLERR,
-                   "Cannot use multiple fix wall commands with "
-                   "pair lubricateU");
-      flagwall = 1; // Walls exist
-      wallfix = dynamic_cast<FixWall *>(modify->fix[i]);
-      if (wallfix->xflag) flagwall = 2; // Moving walls exist
-    }
+  wallfix = nullptr;
+
+  if (modify->get_fix_by_style("^deform").size() > 0) flagdeform = 1;
+  auto fixes = modify->get_fix_by_style("^wall");
+  if (fixes.size() > 1)
+    error->all(FLERR, "Cannot use multiple fix wall commands with pair lubricateU");
+  else if (fixes.size() == 1) {
+    wallfix = dynamic_cast<FixWall *>(fixes[0]);
+    if (!wallfix)
+      error->all(FLERR, "Fix {} is not compatible with pair lubricateU", fixes[0]->style);
+    flagwall = 1;
+    if (wallfix->xflag) flagwall = 2; // Moving walls exist
   }
 
   // set the isotropic constants depending on the volume fraction

src/COLLOID/pair_lubricateU_poly.cpp

@@ -1156,18 +1156,18 @@ void PairLubricateUPoly::init_style()
   // are re-calculated at every step.
 
   flagdeform = flagwall = 0;
-  for (int i = 0; i < modify->nfix; i++) {
-    if (strcmp(modify->fix[i]->style,"deform") == 0)
-      flagdeform = 1;
-    else if (strstr(modify->fix[i]->style,"wall") != nullptr) {
-      if (flagwall)
-        error->all(FLERR,
-                   "Cannot use multiple fix wall commands with "
-                   "pair lubricateU");
-      flagwall = 1; // Walls exist
-      wallfix = dynamic_cast<FixWall *>(modify->fix[i]);
-      if (wallfix->xflag) flagwall = 2; // Moving walls exist
-    }
+  wallfix = nullptr;
+
+  if (modify->get_fix_by_style("^deform").size() > 0) flagdeform = 1;
+  auto fixes = modify->get_fix_by_style("^wall");
+  if (fixes.size() > 1)
+    error->all(FLERR, "Cannot use multiple fix wall commands with pair lubricateU/poly");
+  else if (fixes.size() == 1) {
+    wallfix = dynamic_cast<FixWall *>(fixes[0]);
+    if (!wallfix)
+      error->all(FLERR, "Fix {} is not compatible with pair lubricateU/poly", fixes[0]->style);
+    flagwall = 1;
+    if (wallfix->xflag) flagwall = 2; // Moving walls exist
   }
 
   // set the isotropic constants depending on the volume fraction

src/COLLOID/pair_lubricate_poly.cpp

@@ -460,30 +460,23 @@ void PairLubricatePoly::init_style()
   // are re-calculated at every step.
 
   shearing = flagdeform = flagwall = 0;
-  for (int i = 0; i < modify->nfix; i++) {
-    if (strcmp(modify->fix[i]->style,"deform") == 0) {
-      shearing = flagdeform = 1;
-      if ((dynamic_cast<FixDeform *>(modify->fix[i]))->remapflag != Domain::V_REMAP)
-        error->all(FLERR,"Using pair lubricate with inconsistent "
-                   "fix deform remap option");
-    }
-    if (strstr(modify->fix[i]->style,"wall") != nullptr) {
-      if (flagwall)
-        error->all(FLERR,
-                   "Cannot use multiple fix wall commands with "
-                   "pair lubricate/poly");
-      flagwall = 1; // Walls exist
-      wallfix = dynamic_cast<FixWall *>(modify->fix[i]);
-      if (wallfix->xflag) flagwall = 2; // Moving walls exist
-    }
 
-    if (strstr(modify->fix[i]->style,"wall") != nullptr) {
-      flagwall = 1; // Walls exist
-      if ((dynamic_cast<FixWall *>(modify->fix[i]))->xflag) {
-        flagwall = 2; // Moving walls exist
-        wallfix = dynamic_cast<FixWall *>(modify->fix[i]);
-      }
-    }
+  auto fixes = modify->get_fix_by_style("^deform");
+  if (fixes.size() > 0) {
+    shearing = flagdeform = 1;
+    auto *myfix = dynamic_cast<FixDeform *>(fixes[0]);
+    if (myfix && (myfix->remapflag != Domain::V_REMAP))
+      error->all(FLERR,"Using pair lubricate/poly with inconsistent fix deform remap option");
+  }
+  fixes = modify->get_fix_by_style("^wall");
+  if (fixes.size() > 1)
+    error->all(FLERR, "Cannot use multiple fix wall commands with pair lubricate/poly");
+  else if (fixes.size() == 1) {
+    wallfix = dynamic_cast<FixWall *>(fixes[0]);
+    if (!wallfix)
+      error->all(FLERR, "Fix {} is not compatible with pair lubricate/poly", fixes[0]->style);
+    flagwall = 1;
+    if (wallfix->xflag) flagwall = 2; // Moving walls exist
   }
 
   double vol_T;

src/DIFFRACTION/compute_xrd.cpp

@@ -261,7 +261,7 @@ void ComputeXRD::init()
   double ang = 0.0;
 
   double convf = 360 / MY_PI;
-  if (radflag ==1) convf = 1;
+  if (radflag == 1) convf = 2;
 
   int n = 0;
   for (int m = 0; m < mmax; m++) {

src/DPD-MESO/fix_mvv_dpd.cpp

@@ -24,6 +24,7 @@
 #include "fix_mvv_dpd.h"
 
 #include "atom.h"
+#include "domain.h"
 #include "error.h"
 #include "force.h"
 #include "update.h"
@@ -65,6 +66,11 @@ void FixMvvDPD::init()
   if (!atom->vest_flag)
     error->all(FLERR,"Fix mvv/dpd requires atom attribute vest e.g. from atom style mdpd");
 
+  // Cannot use vremap since its effects aren't propagated to vest
+  //   see RHEO or SPH packages for examples of patches
+  if (domain->deform_vremap)
+    error->all(FLERR, "Fix mvv/dpd cannot be used with velocity remapping");
+
   if (!force->pair_match("^mdpd",0) && !force->pair_match("^dpd",0)) {
     if (force->pair_match("^hybrid",0)) {
       if (!(force->pair_match("^mdpd",0,1) || force->pair_match("^dpd",0),1)) {

src/DPD-MESO/fix_mvv_edpd.cpp

@@ -33,6 +33,7 @@
 #include "fix_mvv_edpd.h"
 
 #include "atom.h"
+#include "domain.h"
 #include "error.h"
 #include "force.h"
 #include "update.h"
@@ -73,6 +74,11 @@ void FixMvvEDPD::init()
 {
   if (!atom->edpd_flag) error->all(FLERR,"Fix mvv/edpd requires atom style edpd");
 
+  // Cannot use vremap since its effects aren't propagated to vest
+  //   see RHEO or SPH packages for examples of patches
+  if (domain->deform_vremap)
+    error->all(FLERR, "Fix mvv/edpd cannot be used with velocity remapping");
+
   if (!force->pair_match("^edpd",0)) {
     if (force->pair_match("^hybrid",0)) {
       if (!force->pair_match("^edpd",0,1)) {

src/DPD-MESO/fix_mvv_tdpd.cpp

@@ -29,6 +29,7 @@
 #include "fix_mvv_tdpd.h"
 
 #include "atom.h"
+#include "domain.h"
 #include "error.h"
 #include "force.h"
 #include "update.h"
@@ -71,6 +72,11 @@ void FixMvvTDPD::init()
 {
   if (!atom->tdpd_flag) error->all(FLERR,"Fix mvv/tdpd requires atom style tdpd");
 
+  // Cannot use vremap since its effects aren't propagated to vest
+  //   see RHEO or SPH packages for examples of patches
+  if (domain->deform_vremap)
+    error->all(FLERR, "Fix mvv/tdpd cannot be used with velocity remapping");
+
   if (!force->pair_match("^tdpd",0)) {
     if (force->pair_match("^hybrid",0)) {
       if (!force->pair_match("^tdpd",0,1)) {

src/DPD-SMOOTH/fix_meso_move.cpp

@@ -350,7 +350,14 @@ void FixMesoMove::init () {
 }
 
 void FixMesoMove::setup_pre_force (int /*vflag*/) {
+
+  // Cannot use vremap since its effects aren't propagated to vest
+  //   see RHEO or SPH packages for examples of patches
+  if (domain->deform_vremap)
+    error->all(FLERR, "Fix meso/move cannot be used with velocity remapping");
+
   // set vest equal to v
+
   double **v = atom->v;
   double **vest = atom->vest;
   int *mask = atom->mask;

src/DPD-SMOOTH/fix_rigid_meso.cpp

@@ -29,11 +29,12 @@
 ------------------------------------------------------------------------- */
 
 #include "fix_rigid_meso.h"
+
 #include "math_extra.h"
 #include "atom.h"
 #include "domain.h"
-#include "memory.h"
 #include "error.h"
+#include "memory.h"
 
 using namespace LAMMPS_NS;
 using namespace FixConst;
@@ -92,6 +93,11 @@ void FixRigidMeso::setup (int vflag) {
     conjqm[ibody][2] *= 2.0;
     conjqm[ibody][3] *= 2.0;
   }
+
+  // Cannot use vremap since its effects aren't propagated to vest
+  //   see RHEO or SPH packages for examples of patches
+  if (domain->deform_vremap)
+    error->all(FLERR, "Fix rigid/meso cannot be used with velocity remapping");
 }
 
 /* ----------------------------------------------------------------------

src/EFF/fix_langevin_eff.cpp

@@ -137,7 +137,7 @@ void FixLangevinEff::post_force_no_tally()
   dof = domain->dimension * particles;
   fix_dof = 0;
   for (int i = 0; i < modify->nfix; i++)
-    fix_dof += modify->fix[i]->dof(igroup);
+    fix_dof += (int)modify->fix[i]->dof(igroup);
 
   // extra_dof = domain->dimension
   dof -= domain->dimension + fix_dof;
@@ -306,7 +306,7 @@ void FixLangevinEff::post_force_tally()
   dof = domain->dimension * particles;
   fix_dof = 0;
   for (int i = 0; i < modify->nfix; i++)
-    fix_dof += modify->fix[i]->dof(igroup);
+    fix_dof += (int)modify->fix[i]->dof(igroup);
 
   // extra_dof = domain->dimension
   dof -= domain->dimension + fix_dof;

src/EFF/fix_nh_eff.cpp

@@ -16,7 +16,6 @@
    Contributing author: Andres Jaramillo-Botero (Caltech)
 ------------------------------------------------------------------------- */
 
-
 #include "fix_nh_eff.h"
 
 #include "atom.h"
@@ -62,7 +61,7 @@ void FixNHEff::nve_v()
     if (mask[i] & groupbit) {
       if (abs(spin[i])==1) {
         dtfm = dtf / mass[type[i]];
-        ervel[i] = dtfm * erforce[i] / mefactor;
+        ervel[i] += dtfm * erforce[i] / mefactor;
       }
     }
   }
@@ -79,15 +78,26 @@ void FixNHEff::nve_x()
   FixNH::nve_x();
 
   double *eradius = atom->eradius;
+  double *erforce = atom->erforce;
   double *ervel = atom->ervel;
+  double *mass = atom->mass;
+  int *type = atom->type;
   int *spin = atom->spin;
   int *mask = atom->mask;
   int nlocal = atom->nlocal;
   if (igroup == atom->firstgroup) nlocal = atom->nfirst;
 
+  double mefactor = domain->dimension/4.0;
+  double dtfm;
+
   for (int i = 0; i < nlocal; i++)
-    if (mask[i] & groupbit)
-      if (abs(spin[i])==1) eradius[i] += dtv * ervel[i];
+    if (mask[i] & groupbit) {
+      dtfm = dtf / mass[type[i]];
+      if (abs(spin[i])==1) {
+        ervel[i] += dtfm * erforce[i] / mefactor;
+        eradius[i] += dtv * ervel[i];
+      }
+    }
 }
 
 /* ----------------------------------------------------------------------

src/ELECTRODE/fix_electrode_conp.cpp

@@ -494,7 +494,7 @@ void FixElectrodeConp::setup_post_neighbor()
     if (read_mat)
       read_from_file(input_file_mat, elastance, "elastance");
     else if (!read_inv) {
-      if (etypes_neighlists) neighbor->build_one(mat_neighlist, 0);
+      if (etypes_neighlists) neighbor->build_one(mat_neighlist);
       auto array_compute = std::unique_ptr<ElectrodeMatrix>(new ElectrodeMatrix(lmp, igroup, eta));
       array_compute->setup(tag_to_iele, pair, mat_neighlist);
       if (tfflag) { array_compute->setup_tf(tf_types); }

src/ELECTRODE/pppm_electrode.cpp

@@ -633,7 +633,9 @@ void PPPMElectrode::project_psi(double *vec, int sensor_grpbit)
   // project u_brick with weight matrix
   double **x = atom->x;
   int *mask = atom->mask;
-  double const scaleinv = 1.0 / (nx_pppm * ny_pppm * nz_pppm);
+  const bigint ngridtotal = (bigint) nx_pppm * ny_pppm * nz_pppm;
+  const double scaleinv = 1.0 / ngridtotal;
+
   for (int i = 0; i < atom->nlocal; i++) {
     if (!(mask[i] & sensor_grpbit)) continue;
     double v = 0.;
@@ -1362,7 +1364,7 @@ double PPPMElectrode::compute_qopt()
   // each proc calculates contributions from every Pth grid point
 
   bigint ngridtotal = (bigint) nx_pppm * ny_pppm * nz_pppm;
-  int nxy_pppm = nx_pppm * ny_pppm;
+  bigint nxy_pppm = (bigint) nx_pppm * ny_pppm;
 
   double qopt = 0.0;
 

src/EXTRA-COMPUTE/compute_hexorder_atom.cpp

@@ -147,10 +147,8 @@ void ComputeHexOrderAtom::compute_peratom()
   }
 
   // invoke full neighbor list (will copy or build if necessary)
-  // on the first step of a run, set preflag to one in neighbor->build_one(...)
 
-  if (update->firststep == update->ntimestep) neighbor->build_one(list,1);
-  else neighbor->build_one(list);
+  neighbor->build_one(list);
 
   inum = list->inum;
   ilist = list->ilist;

src/EXTRA-COMPUTE/compute_hma.cpp

@@ -291,8 +291,7 @@ void ComputeHMA::compute_vector()
       double *special_coul = force->special_coul;
       int newton_pair = force->newton_pair;
 
-      if (update->firststep == update->ntimestep) neighbor->build_one(list,1);
-      else neighbor->build_one(list);
+      neighbor->build_one(list);
       int inum = list->inum;
       int *ilist = list->ilist;
       int *numneigh = list->numneigh;

src/EXTRA-COMPUTE/compute_stress_cartesian.cpp

@@ -386,7 +386,7 @@ void ComputeStressCartesian::compute_array()
         if (tag[i] > tag[j]) {
           if ((tag[i] + tag[j]) % 2 == 0) continue;
         } else if (tag[i] < tag[j]) {
-          if ((tag[i] < tag[j]) % 2 == 1) continue;
+          if ((tag[i] + tag[j]) % 2 == 1) continue;
         }
       }
 

src/EXTRA-DUMP/dump_xtc.cpp

@@ -14,7 +14,7 @@
 
 /* ----------------------------------------------------------------------
    Contributing authors: Naveen Michaud-Agrawal (Johns Hopkins U)
-                         open-source XDR routines from
+                         Open Source XDR based I/O routines from
                            Frans van Hoesel (https://www.rug.nl/staff/f.h.j.van.hoesel/)
                            are included in this file
                          Axel Kohlmeyer (Temple U)
@@ -35,27 +35,29 @@
 #include "output.h"
 #include "update.h"
 
+#include "xdr_compat.h"
+
 #include <climits>
 #include <cmath>
 #include <cstring>
 
 using namespace LAMMPS_NS;
 
-#define EPS 1e-5
-#define XTC_MAGIC 1995
+static constexpr double EPS = 1.0e-5;
+static constexpr int XTC_MAGIC = 1995;
 
 #define MYMIN(a,b) ((a) < (b) ? (a) : (b))
 #define MYMAX(a,b) ((a) > (b) ? (a) : (b))
 
-int xdropen(XDR *, const char *, const char *);
-int xdrclose(XDR *);
-void xdrfreebuf();
-int xdr3dfcoord(XDR *, float *, int *, float *);
+static int xdropen(XDR *, const char *, const char *);
+static int xdrclose(XDR *);
+static void xdrfreebuf();
+static int xdr3dfcoord(XDR *, float *, int *, float *);
 
 /* ---------------------------------------------------------------------- */
 
-DumpXTC::DumpXTC(LAMMPS *lmp, int narg, char **arg) : Dump(lmp, narg, arg),
-  coords(nullptr)
+DumpXTC::DumpXTC(LAMMPS *lmp, int narg, char **arg)
+  : Dump(lmp, narg, arg), coords(nullptr), xd(nullptr)
 {
   if (narg != 5) error->all(FLERR,"Illegal dump xtc command");
   if (binary || compressed || multifile || multiproc)
@@ -68,6 +70,7 @@ DumpXTC::DumpXTC(LAMMPS *lmp, int narg, char **arg) : Dump(lmp, narg, arg),
   flush_flag = 0;
   unwrap_flag = 0;
   precision = 1000.0;
+  xd = new XDR;
 
   // allocate global array for atom coords
 
@@ -105,9 +108,10 @@ DumpXTC::~DumpXTC()
   memory->destroy(coords);
 
   if (me == 0) {
-    xdrclose(&xd);
+    xdrclose(xd);
     xdrfreebuf();
   }
+  delete xd;
 }
 
 /* ---------------------------------------------------------------------- */
@@ -150,7 +154,8 @@ void DumpXTC::openfile()
 
   fp = nullptr;
   if (me == 0)
-    if (xdropen(&xd,filename,"w") == 0) error->one(FLERR,"Cannot open dump file");
+    if (xdropen(xd,filename,"w") == 0)
+      error->one(FLERR,"Cannot open XTC format dump file {}: {}", filename, utils::getsyserror());
 }
 
 /* ---------------------------------------------------------------------- */
@@ -176,11 +181,11 @@ void DumpXTC::write_header(bigint nbig)
   if (me != 0) return;
 
   int tmp = XTC_MAGIC;
-  xdr_int(&xd,&tmp);
-  xdr_int(&xd,&n);
-  xdr_int(&xd,&ntimestep);
+  xdr_int(xd,&tmp);
+  xdr_int(xd,&n);
+  xdr_int(xd,&ntimestep);
   float time_value = ntimestep * tfactor * update->dt;
-  xdr_float(&xd,&time_value);
+  xdr_float(xd,&time_value);
 
   // cell basis vectors
   if (domain->triclinic) {
@@ -192,18 +197,18 @@ void DumpXTC::write_header(bigint nbig)
     float xz = sfactor * domain->xz;
     float yz = sfactor * domain->yz;
 
-    xdr_float(&xd,&xdim); xdr_float(&xd,&zero); xdr_float(&xd,&zero);
-    xdr_float(&xd,&xy  ); xdr_float(&xd,&ydim); xdr_float(&xd,&zero);
-    xdr_float(&xd,&xz  ); xdr_float(&xd,&yz  ); xdr_float(&xd,&zdim);
+    xdr_float(xd,&xdim); xdr_float(xd,&zero); xdr_float(xd,&zero);
+    xdr_float(xd,&xy  ); xdr_float(xd,&ydim); xdr_float(xd,&zero);
+    xdr_float(xd,&xz  ); xdr_float(xd,&yz  ); xdr_float(xd,&zdim);
   } else {
     float zero = 0.0;
     float xdim = sfactor * (domain->boxhi[0] - domain->boxlo[0]);
     float ydim = sfactor * (domain->boxhi[1] - domain->boxlo[1]);
     float zdim = sfactor * (domain->boxhi[2] - domain->boxlo[2]);
 
-    xdr_float(&xd,&xdim); xdr_float(&xd,&zero); xdr_float(&xd,&zero);
-    xdr_float(&xd,&zero); xdr_float(&xd,&ydim); xdr_float(&xd,&zero);
-    xdr_float(&xd,&zero); xdr_float(&xd,&zero); xdr_float(&xd,&zdim);
+    xdr_float(xd,&xdim); xdr_float(xd,&zero); xdr_float(xd,&zero);
+    xdr_float(xd,&zero); xdr_float(xd,&ydim); xdr_float(xd,&zero);
+    xdr_float(xd,&zero); xdr_float(xd,&zero); xdr_float(xd,&zdim);
   }
 }
 
@@ -328,7 +333,7 @@ double DumpXTC::memory_usage()
 
 void DumpXTC::write_frame()
 {
-  xdr3dfcoord(&xd,coords,&natoms,&precision);
+  xdr3dfcoord(xd,coords,&natoms,&precision);
 }
 
 // ----------------------------------------------------------------------
@@ -406,7 +411,7 @@ static int magicints[] = {
  |
  | xdropen - open xdr file
  |
- | This versions differs from xdrstdio_create, because I need to know
+ | This version differs from xdrstdio_create, because I need to know
  | the state of the file (read or write) so I can use xdr3dfcoord
  | in eigther read or write mode, and the file descriptor
  | so I can close the file (something xdr_destroy doesn't do).
@@ -1048,7 +1053,7 @@ int xdr3dfcoord(XDR *xdrs, float *fp, int *size, float *precision)
     }
     if (buf[1] != 0) buf[0]++;
     xdr_int(xdrs, &(buf[0])); /* buf[0] holds the length in bytes */
-    return errval * (xdr_opaque(xdrs, (caddr_t)&(buf[3]), (u_int)buf[0]));
+    return errval * (xdr_opaque(xdrs, (char *)&(buf[3]), (unsigned int)buf[0]));
   } else {
 
     /* xdrs is open for reading */
@@ -1129,7 +1134,7 @@ int xdr3dfcoord(XDR *xdrs, float *fp, int *size, float *precision)
 
     if (xdr_int(xdrs, &(buf[0])) == 0)
       return 0;
-    if (xdr_opaque(xdrs, (caddr_t)&(buf[3]), (u_int)buf[0]) == 0)
+    if (xdr_opaque(xdrs, (char *)&(buf[3]), (unsigned int)buf[0]) == 0)
       return 0;
     buf[0] = buf[1] = buf[2] = 0;
 

src/EXTRA-DUMP/dump_xtc.h

@@ -21,8 +21,8 @@ DumpStyle(xtc,DumpXTC);
 #define LMP_DUMP_XTC_H
 
 #include "dump.h"
-#include "xdr_compat.h"
 
+struct XDR;
 namespace LAMMPS_NS {
 
 class DumpXTC : public Dump {
@@ -37,7 +37,7 @@ class DumpXTC : public Dump {
   float precision;    // user-adjustable precision setting
   float *coords;
   double sfactor, tfactor;    // scaling factors for positions and time unit
-  XDR xd;
+  XDR *xd;
 
   void init_style() override;
   int modify_param(int, char **) override;

src/EXTRA-DUMP/xdr_compat.cpp

@@ -1,54 +1,63 @@
-// clang-format off
 #include "xdr_compat.h"
+
 #include <cstdlib>
 #include <cstring>
 
-/* This file is needed for systems, that do not provide XDR support
- * in their system libraries. It was written for windows, but will
- * most probably work on other platforms too. better make sure you
- * test that the xtc files produced are ok before using it.
+/*
+ * This file contains an implementation of the Sun External Data Representation (XDR)
+ * routines.  They have been adapted specifically for the use with the LAMMPS xtc dump
+ * style to produce compressed trajectory files in the Gromacs XTC format.
  *
- * It is also needed on BG/L and Cray XT3/XT4 as we don't have
- * XDR support in the lightweight kernel runtimes either.
+ * The XDR sources are available under the BSD 3-clause license for example in
+ * the MIT Kerberos 5 distribution with the following copyright notice and license.
  *
- * This file contains the definitions for Sun External Data
- * Representation (XDR) headers and routines.
+ * @(#)xdr.h    2.2 88/07/29 4.0 RPCSRC
  *
- * Although the rest of LAMPPS is GPL, you can copy and use the XDR
- * routines in any way you want as long as you obey Sun's license:
+ * Copyright (c) 2010, Oracle America, Inc.
  *
- * Sun RPC is a product of Sun Microsystems, Inc. and is provided for
- * unrestricted use provided that this legend is included on all tape
- * media and as a part of the software program in whole or part.  Users
- * may copy or modify Sun RPC without charge, but are not authorized
- * to license or distribute it to anyone else except as part of a product or
- * program developed by the user.
+ * All rights reserved.
  *
- * SUN RPC IS PROVIDED AS IS WITH NO WARRANTIES OF ANY KIND INCLUDING THE
- * WARRANTIES OF DESIGN, MERCHANTIBILITY AND FITNESS FOR A PARTICULAR
- * PURPOSE, OR ARISING FROM A COURSE OF DEALING, USAGE OR TRADE PRACTICE.
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions are met:
  *
- * Sun RPC is provided with no support and without any obligation on the
- * part of Sun Microsystems, Inc. to assist in its use, correction,
- * modification or enhancement.
+ *     * Redistributions of source code must retain the above copyright
+ *       notice, this list of conditions and the following disclaimer.
  *
- * SUN MICROSYSTEMS, INC. SHALL HAVE NO LIABILITY WITH RESPECT TO THE
- * INFRINGEMENT OF COPYRIGHTS, TRADE SECRETS OR ANY PATENTS BY SUN RPC
- * OR ANY PART THEREOF.
+ *     * Redistributions in binary form must reproduce the above copyright
+ *       notice, this list of conditions and the following disclaimer in
+ *       the documentation and/or other materials provided with the
+ *       distribution.
  *
- * In no event will Sun Microsystems, Inc. be liable for any lost revenue
- * or profits or other special, indirect and consequential damages, even if
- * Sun has been advised of the possibility of such damages.
+ *     * Neither the name of the "Oracle America, Inc." nor the names of
+ *       its contributors may be used to endorse or promote products
+ *       derived from this software without specific prior written permission.
  *
- * Sun Microsystems, Inc.
- * 2550 Garcia Avenue
- * Mountain View, California  94043
+ * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+ * IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+ * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+ * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
+ * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+ * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
+ * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+ * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+ * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+ * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  */
 
 #ifdef __cplusplus
 extern "C" {
 #endif
 
+#ifndef FALSE
+#define FALSE (0)
+#endif
+#ifndef TRUE
+#define TRUE (1)
+#endif
+
+#define BYTES_PER_XDR_UNIT (4)
+
 /*
  * for unit alignment
  */
@@ -58,19 +67,18 @@ static xdr_uint32_t xdr_swapbytes(xdr_uint32_t x)
 {
   xdr_uint32_t y;
   int i;
-  char *px=(char *)&x;
-  char *py=(char *)&y;
+  char *px = (char *) &x;
+  char *py = (char *) &y;
 
-  for (i=0;i<4;i++)
-    py[i]=px[3-i];
+  for (i = 0; i < 4; i++) py[i] = px[3 - i];
 
   return y;
 }
 
 static xdr_uint32_t xdr_htonl(xdr_uint32_t x)
 {
-  short s=0x0F00;
-  if (*((char *)&s)==(char)0x0F) {
+  short s = 0x0F00;
+  if (*((char *) &s) == (char) 0x0F) {
     /* bigendian, do nothing */
     return x;
   } else {
@@ -81,8 +89,8 @@ static xdr_uint32_t xdr_htonl(xdr_uint32_t x)
 
 static xdr_uint32_t xdr_ntohl(xdr_uint32_t x)
 {
-  short s=0x0F00;
-  if (*((char *)&s)==(char)0x0F) {
+  short s = 0x0F00;
+  if (*((char *) &s) == (char) 0x0F) {
     /* bigendian, do nothing */
     return x;
   } else {
@@ -91,221 +99,39 @@ static xdr_uint32_t xdr_ntohl(xdr_uint32_t x)
   }
 }
 
-
-/*
- * Free a data structure using XDR
- * Not a filter, but a convenient utility nonetheless
- */
-void
-xdr_free (xdrproc_t proc, char *objp)
-{
-  XDR x;
-
-  x.x_op = XDR_FREE;
-  (*proc) (&x, objp);
-}
-
-/*
- * XDR nothing
- */
-bool_t
-xdr_void (void)
-{
-  return TRUE;
-}
-
 /*
  * XDR integers
  */
-bool_t
-xdr_int (XDR *xdrs, int *ip)
+bool_t xdr_int(XDR *xdrs, int *ip)
 {
   xdr_int32_t l;
 
-  switch (xdrs->x_op)
-  {
+  switch (xdrs->x_op) {
+
     case XDR_ENCODE:
       l = (xdr_int32_t) (*ip);
-      return xdr_putint32 (xdrs, &l);
+      return xdr_putint32(xdrs, &l);
+      break;
 
     case XDR_DECODE:
-      if (!xdr_getint32 (xdrs, &l))
-          {
-            return FALSE;
-          }
+      if (!xdr_getint32(xdrs, &l)) return FALSE;
       *ip = (int) l;
+      return TRUE;
+      break;
 
     case XDR_FREE:
       return TRUE;
+      break;
   }
   return FALSE;
 }
 
-
-/*
- * XDR unsigned integers
- */
-bool_t
-xdr_u_int (XDR *xdrs, unsigned int *up)
-{
-  xdr_uint32_t l;
-
-  switch (xdrs->x_op)
-  {
-    case XDR_ENCODE:
-      l = (xdr_uint32_t) (*up);
-      return xdr_putuint32 (xdrs, &l);
-
-    case XDR_DECODE:
-      if (!xdr_getuint32 (xdrs, &l))
-          {
-            return FALSE;
-          }
-      *up = (unsigned int) l;
-
-    case XDR_FREE:
-      return TRUE;
-  }
-  return FALSE;
-}
-
-
-
-
-/*
- * XDR short integers
- */
-bool_t
-xdr_short (XDR *xdrs, short *sp)
-{
-  xdr_int32_t l;
-
-  switch (xdrs->x_op)
-    {
-    case XDR_ENCODE:
-      l = (xdr_int32_t) *sp;
-      return xdr_putint32 (xdrs, &l);
-
-    case XDR_DECODE:
-      if (!xdr_getint32 (xdrs, &l))
-        {
-          return FALSE;
-        }
-      *sp = (short) l;
-      return TRUE;
-
-    case XDR_FREE:
-      return TRUE;
-    }
-  return FALSE;
-}
-
-
-/*
- * XDR unsigned short integers
- */
-bool_t
-xdr_u_short (XDR *xdrs, unsigned short *usp)
-{
-  xdr_uint32_t l;
-
-  switch (xdrs->x_op)
-    {
-    case XDR_ENCODE:
-      l = (xdr_uint32_t) *usp;
-      return xdr_putuint32 (xdrs, &l);
-
-    case XDR_DECODE:
-      if (!xdr_getuint32 (xdrs, &l))
-        {
-          return FALSE;
-        }
-          *usp = (unsigned short) l;
-      return TRUE;
-
-    case XDR_FREE:
-      return TRUE;
-    }
-  return FALSE;
-}
-
-
-/*
- * XDR a char
- */
-bool_t
-xdr_char (XDR *xdrs, char *cp)
-{
-  int i;
-
-  i = (*cp);
-  if (!xdr_int (xdrs, &i))
-    {
-      return FALSE;
-    }
-  *cp = i;
-  return TRUE;
-}
-
-/*
- * XDR an unsigned char
- */
-bool_t
-xdr_u_char (XDR *xdrs, unsigned char *cp)
-{
-  unsigned int u;
-
-  u = (*cp);
-  if (!xdr_u_int (xdrs, &u))
-    {
-      return FALSE;
-    }
-  *cp = u;
-  return TRUE;
-}
-
-/*
- * XDR booleans
- */
-bool_t
-xdr_bool (XDR *xdrs, int *bp)
-{
-#define XDR_FALSE        ((xdr_int32_t) 0)
-#define XDR_TRUE        ((xdr_int32_t) 1)
-
-  xdr_int32_t lb;
-
-  switch (xdrs->x_op)
-    {
-    case XDR_ENCODE:
-      lb = *bp ? XDR_TRUE : XDR_FALSE;
-      return xdr_putint32 (xdrs, &lb);
-
-    case XDR_DECODE:
-      if (!xdr_getint32 (xdrs, &lb))
-        {
-          return FALSE;
-        }
-      *bp = (lb == XDR_FALSE) ? FALSE : TRUE;
-      return TRUE;
-
-    case XDR_FREE:
-      return TRUE;
-    }
-  return FALSE;
-#undef XDR_FALSE
-#undef XDR_TRUE
-}
-
-
-
 /*
  * XDR opaque data
  * Allows the specification of a fixed size sequence of opaque bytes.
  * cp points to the opaque object and cnt gives the byte length.
  */
-bool_t
-xdr_opaque (XDR *xdrs, char *cp, unsigned int cnt)
+bool_t xdr_opaque(XDR *xdrs, char *cp, unsigned int cnt)
 {
   unsigned int rndup;
   static char crud[BYTES_PER_XDR_UNIT];
@@ -313,219 +139,61 @@ xdr_opaque (XDR *xdrs, char *cp, unsigned int cnt)
   /*
    * if no data we are done
    */
-  if (cnt == 0)
-    return TRUE;
+  if (cnt == 0) return TRUE;
 
   /*
    * round byte count to full xdr units
    */
   rndup = cnt % BYTES_PER_XDR_UNIT;
-  if (rndup > 0)
-    rndup = BYTES_PER_XDR_UNIT - rndup;
-
-  switch (xdrs->x_op)
-    {
-    case XDR_DECODE:
-      if (!xdr_getbytes (xdrs, cp, cnt))
-        {
-          return FALSE;
-        }
-      if (rndup == 0)
-        return TRUE;
-      return xdr_getbytes (xdrs, (char *)crud, rndup);
-
-    case XDR_ENCODE:
-      if (!xdr_putbytes (xdrs, cp, cnt))
-        {
-          return FALSE;
-        }
-      if (rndup == 0)
-        return TRUE;
-      return xdr_putbytes (xdrs, xdr_zero, rndup);
-
-    case XDR_FREE:
-      return TRUE;
-    }
-  return FALSE;
-}
-
-
-/*
- * XDR null terminated ASCII strings
- * xdr_string deals with "C strings" - arrays of bytes that are
- * terminated by a nullptr character.  The parameter cpp references a
- * pointer to storage; If the pointer is null, then the necessary
- * storage is allocated.  The last parameter is the max allowed length
- * of the string as specified by a protocol.
- */
-bool_t
-xdr_string (XDR *xdrs, char **cpp, unsigned int maxsize)
-{
-  char *sp = *cpp;        /* sp is the actual string pointer */
-  unsigned int size = 0;
-  unsigned int nodesize = 0;
-
-  /*
-   * first deal with the length since xdr strings are counted-strings
-   */
-  switch (xdrs->x_op)
-    {
-    case XDR_FREE:
-      if (sp == nullptr)
-        {
-          return TRUE;                /* already free */
-        }
-      /* fall through... */
-    case XDR_ENCODE:
-      if (sp == nullptr)
-            return FALSE;
-      size = strlen (sp);
-      break;
-    case XDR_DECODE:
-      break;
-    }
-
-  if (!xdr_u_int (xdrs, &size))
-    {
-      return FALSE;
-    }
-  if (size > maxsize)
-    {
-      return FALSE;
-    }
-  nodesize = size + 1;
-
-  /*
-   * now deal with the actual bytes
-   */
-  switch (xdrs->x_op)
-    {
-    case XDR_DECODE:
-      if (nodesize == 0)
-        {
-          return TRUE;
-        }
-      if (sp == nullptr)
-        *cpp = sp = (char *) malloc (nodesize);
-      if (sp == nullptr)
-        {
-          (void) fputs ("xdr_string: out of memory\n", stderr);
-          return FALSE;
-        }
-      sp[size] = 0;
-      /* fall into ... */
-
-    case XDR_ENCODE:
-      return xdr_opaque (xdrs, sp, size);
-
-    case XDR_FREE:
-      free (sp);
-      *cpp = nullptr;
-      return TRUE;
-    }
-  return FALSE;
-}
-
-
-
-/* Floating-point stuff */
-
-bool_t
-xdr_float(XDR *xdrs, float *fp)
-{
-        xdr_int32_t tmp;
-
-        switch (xdrs->x_op) {
-
-        case XDR_ENCODE:
-                tmp = *(xdr_int32_t *)fp;
-               return (xdr_putint32(xdrs, &tmp));
-
-                break;
-
-        case XDR_DECODE:
-                        if (xdr_getint32(xdrs, &tmp)) {
-                                *(xdr_int32_t *)fp = tmp;
-                                return (TRUE);
-                        }
-
-                break;
-
-        case XDR_FREE:
-                return (TRUE);
-        }
-        return (FALSE);
-}
-
-
-bool_t
-xdr_double(XDR *xdrs, double *dp)
-{
-
-  /* Windows and some other systems dont define double-precision
-   * word order in the header files, so unfortunately we have
-   * to calculate it!
-   */
-  static int LSW=-1; /* Least significant fp word */
-  int *ip;
-  xdr_int32_t tmp[2];
-
-  if (LSW<0) {
-    double x=0.987654321; /* Just a number */
-
-    /* Possible representations in IEEE double precision:
-     * (S=small endian, B=big endian)
-     *
-     * Byte order, Word order, Hex
-     *     S           S       b8 56 0e 3c dd 9a ef 3f
-     *     B           S       3c 0e 56 b8 3f ef 9a dd
-     *     S           B       dd 9a ef 3f b8 56 0e 3c
-     *     B           B       3f ef 9a dd 3c 0e 56 b8
-     */
-
-    unsigned char ix = *((char *)&x);
-
-    if (ix==0xdd || ix==0x3f)
-      LSW=1;  /* Big endian word order */
-    else if (ix==0xb8 || ix==0x3c)
-      LSW=0;  /* Small endian word order */
-    else { /* Catch strange errors */
-      printf("Error when detecting floating-point word order.\n"
-             "Do you have a non-IEEE system?\n"
-             "If possible, use the XDR libraries provided with your system,\n"
-             "instead of the Gromacs fallback XDR source.\n");
-      exit(0);
-    }
-  }
+  if (rndup > 0) rndup = BYTES_PER_XDR_UNIT - rndup;
 
   switch (xdrs->x_op) {
 
-  case XDR_ENCODE:
-    ip = (int *)dp;
-    tmp[0] = ip[!LSW];
-    tmp[1] = ip[LSW];
-    return (xdr_putint32(xdrs, tmp) &&
-               xdr_putint32(xdrs, tmp+1));
+    case XDR_DECODE:
+      if (!xdr_getbytes(xdrs, cp, cnt)) { return FALSE; }
+      if (rndup == 0) return TRUE;
+      return xdr_getbytes(xdrs, (char *) crud, rndup);
+      break;
 
-    break;
+    case XDR_ENCODE:
+      if (!xdr_putbytes(xdrs, cp, cnt)) { return FALSE; }
+      if (rndup == 0) return TRUE;
+      return xdr_putbytes(xdrs, xdr_zero, rndup);
+      break;
 
-  case XDR_DECODE:
-    ip = (int *)dp;
-    if (xdr_getint32(xdrs, tmp+!LSW) &&
-           xdr_getint32(xdrs, tmp+LSW)) {
-        ip[0] = tmp[0];
-        ip[1] = tmp[1];
-        return (TRUE);
-    }
-
-    break;
-
-  case XDR_FREE:
-    return (TRUE);
+    case XDR_FREE:
+      return TRUE;
+      break;
   }
-  return (FALSE);
+  return FALSE;
 }
 
+/* Floating-point stuff */
+
+bool_t xdr_float(XDR *xdrs, float *fp)
+{
+  xdr_int32_t tmp;
+
+  switch (xdrs->x_op) {
+
+    case XDR_ENCODE:
+      tmp = *(xdr_int32_t *) fp;
+      return (xdr_putint32(xdrs, &tmp));
+      break;
+
+    case XDR_DECODE:
+      if (xdr_getint32(xdrs, &tmp)) {
+        *(xdr_int32_t *) fp = tmp;
+        return TRUE;
+      }
+      break;
+
+    case XDR_FREE:
+      return TRUE;
+      break;
+  }
+  return FALSE;
+}
 
 /* Array routines */
 
@@ -539,55 +207,37 @@ xdr_double(XDR *xdrs, double *dp)
  * > elemsize: size of each element
  * > xdr_elem: routine to XDR each element
  */
-bool_t
-xdr_vector (XDR *xdrs, char *basep, unsigned int nelem,
-            unsigned int elemsize, xdrproc_t xdr_elem)
+bool_t xdr_vector(XDR *xdrs, char *basep, unsigned int nelem, unsigned int elemsize,
+                  xdrproc_t xdr_elem)
 {
-#define LASTUNSIGNED        ((unsigned int)0-1)
+#define LASTUNSIGNED ((unsigned int) 0 - 1)
   unsigned int i;
   char *elptr;
 
   elptr = basep;
-  for (i = 0; i < nelem; i++)
-    {
-      if (!(*xdr_elem) (xdrs, elptr, LASTUNSIGNED))
-        {
-          return FALSE;
-        }
-      elptr += elemsize;
-    }
+  for (i = 0; i < nelem; i++) {
+    if (!(*xdr_elem)(xdrs, elptr, LASTUNSIGNED)) { return FALSE; }
+    elptr += elemsize;
+  }
   return TRUE;
 #undef LASTUNSIGNED
 }
 
-
-
-static bool_t xdrstdio_getbytes (XDR *, char *, unsigned int);
-static bool_t xdrstdio_putbytes (XDR *, char *, unsigned int);
-static unsigned int xdrstdio_getpos (XDR *);
-static bool_t xdrstdio_setpos (XDR *, unsigned int);
-static xdr_int32_t *xdrstdio_inline (XDR *, int);
-static void xdrstdio_destroy (XDR *);
-static bool_t xdrstdio_getint32 (XDR *, xdr_int32_t *);
-static bool_t xdrstdio_putint32 (XDR *, xdr_int32_t *);
-static bool_t xdrstdio_getuint32 (XDR *, xdr_uint32_t *);
-static bool_t xdrstdio_putuint32 (XDR *, xdr_uint32_t *);
+static bool_t xdrstdio_getbytes(XDR *, char *, unsigned int);
+static bool_t xdrstdio_putbytes(XDR *, char *, unsigned int);
+static void xdrstdio_destroy(XDR *);
+static bool_t xdrstdio_getint32(XDR *, xdr_int32_t *);
+static bool_t xdrstdio_putint32(XDR *, xdr_int32_t *);
 
 /*
  * Ops vector for stdio type XDR
  */
-static const struct xdr_ops xdrstdio_ops =
-{
-  xdrstdio_getbytes,               /* deserialize counted bytes */
-  xdrstdio_putbytes,             /* serialize counted bytes */
-  xdrstdio_getpos,                /* get offset in the stream */
-  xdrstdio_setpos,                /* set offset in the stream */
-  xdrstdio_inline,                /* prime stream for inline macros */
-  xdrstdio_destroy,                /* destroy stream */
-  xdrstdio_getint32,        /* deserialize a int */
-  xdrstdio_putint32,        /* serialize a int */
-  xdrstdio_getuint32,        /* deserialize a int */
-  xdrstdio_putuint32                /* serialize a int */
+static const struct xdr_ops xdrstdio_ops = {
+    xdrstdio_getbytes, /* deserialize counted bytes */
+    xdrstdio_putbytes, /* serialize counted bytes */
+    xdrstdio_destroy,  /* destroy stream */
+    xdrstdio_getint32, /* deserialize a int */
+    xdrstdio_putint32, /* serialize a int */
 };
 
 /*
@@ -595,8 +245,7 @@ static const struct xdr_ops xdrstdio_ops =
  * Sets the xdr stream handle xdrs for use on the stream file.
  * Operation flag is set to op.
  */
-void
-xdrstdio_create (XDR *xdrs, FILE *file, enum xdr_op op)
+void xdrstdio_create(XDR *xdrs, FILE *file, enum xdr_op op)
 {
   xdrs->x_op = op;
   /* We have to add the const since the `struct xdr_ops' in `struct XDR'
@@ -611,104 +260,42 @@ xdrstdio_create (XDR *xdrs, FILE *file, enum xdr_op op)
  * Destroy a stdio xdr stream.
  * Cleans up the xdr stream handle xdrs previously set up by xdrstdio_create.
  */
-static void
-xdrstdio_destroy (XDR *xdrs)
+static void xdrstdio_destroy(XDR *xdrs)
 {
-  (void) fflush ((FILE *) xdrs->x_private);
+  (void) fflush((FILE *) xdrs->x_private);
   /* xx should we close the file ?? */
 }
 
-
-static bool_t
-xdrstdio_getbytes (XDR *xdrs, char *addr, unsigned int len)
+static bool_t xdrstdio_getbytes(XDR *xdrs, char *addr, unsigned int len)
 {
-  if ((len != 0) && (fread (addr, (int) len, 1,
-                            (FILE *) xdrs->x_private) != 1))
-    return FALSE;
+  if ((len != 0) && (fread(addr, (int) len, 1, (FILE *) xdrs->x_private) != 1)) return FALSE;
   return TRUE;
 }
 
-static bool_t
-xdrstdio_putbytes (XDR *xdrs, char *addr, unsigned int len)
+static bool_t xdrstdio_putbytes(XDR *xdrs, char *addr, unsigned int len)
 {
-  if ((len != 0) && (fwrite (addr, (int) len, 1,
-                             (FILE *) xdrs->x_private) != 1))
-    return FALSE;
+  if ((len != 0) && (fwrite(addr, (int) len, 1, (FILE *) xdrs->x_private) != 1)) return FALSE;
   return TRUE;
 }
 
-static unsigned int
-xdrstdio_getpos (XDR *xdrs)
-{
-  return (unsigned int) ftell ((FILE *) xdrs->x_private);
-}
-
-static bool_t
-xdrstdio_setpos (XDR *xdrs, unsigned int pos)
-{
-  return fseek ((FILE *) xdrs->x_private, (xdr_int32_t) pos, 0) < 0 ? FALSE : TRUE;
-}
-
-static xdr_int32_t *
-xdrstdio_inline (XDR * /*xdrs*/, int /*len*/)
-{
-  /*
-   * Must do some work to implement this: must ensure
-   * enough data in the underlying stdio buffer,
-   * that the buffer is aligned so that we can indirect through a
-   * long *, and stuff this pointer in xdrs->x_buf.  Doing
-   * a fread or fwrite to a scratch buffer would defeat
-   * most of the gains to be had here and require storage
-   * management on this buffer, so we don't do this.
-   */
-  return nullptr;
-}
-
-static bool_t
-xdrstdio_getint32 (XDR *xdrs, xdr_int32_t *ip)
+static bool_t xdrstdio_getint32(XDR *xdrs, xdr_int32_t *ip)
 {
   xdr_int32_t mycopy;
 
-  if (fread ((char *) &mycopy, 4, 1, (FILE *) xdrs->x_private) != 1)
-    return FALSE;
-  *ip = xdr_ntohl (mycopy);
+  if (fread((char *) &mycopy, 4, 1, (FILE *) xdrs->x_private) != 1) return FALSE;
+  *ip = xdr_ntohl(mycopy);
   return TRUE;
 }
 
-static bool_t
-xdrstdio_putint32 (XDR *xdrs, xdr_int32_t *ip)
+static bool_t xdrstdio_putint32(XDR *xdrs, xdr_int32_t *ip)
 {
-  xdr_int32_t mycopy = xdr_htonl (*ip);
+  xdr_int32_t mycopy = xdr_htonl(*ip);
 
   ip = &mycopy;
-  if (fwrite ((char *) ip, 4, 1, (FILE *) xdrs->x_private) != 1)
-    return FALSE;
+  if (fwrite((char *) ip, 4, 1, (FILE *) xdrs->x_private) != 1) return FALSE;
   return TRUE;
 }
 
-static bool_t
-xdrstdio_getuint32 (XDR *xdrs, xdr_uint32_t *ip)
-{
-        xdr_uint32_t mycopy;
-
-        if (fread ((char *) &mycopy, 4, 1, (FILE *) xdrs->x_private) != 1)
-                return FALSE;
-        *ip = xdr_ntohl (mycopy);
-        return TRUE;
-}
-
-static bool_t
-xdrstdio_putuint32 (XDR *xdrs, xdr_uint32_t *ip)
-{
-        xdr_uint32_t mycopy = xdr_htonl (*ip);
-
-        ip = &mycopy;
-        if (fwrite ((char *) ip, 4, 1, (FILE *) xdrs->x_private) != 1)
-                return FALSE;
-        return TRUE;
-}
-
 #ifdef __cplusplus
 }
 #endif
-

src/EXTRA-DUMP/xdr_compat.h

@@ -1,7 +1,7 @@
 #ifndef LMP_XDR_COMPAT_H
 #define LMP_XDR_COMPAT_H
 
-#include <climits>
+#include <cstdint>
 #include <cstdio>
 
 #ifdef __cplusplus
@@ -9,47 +9,53 @@ extern "C" {
 #endif
 
 /*
- * This file is needed for systems, that do not provide XDR support
- * in their system libraries. It was written for windows, but will
- * most probably work on other platforms too. better make sure you
- * test that the xtc files produced are ok before using it.
+ * This file contains the definitions for Sun External Data Representation (XDR).
+ * They have been adapted specifically for the use with the LAMMPS xtc dump style
+ * to produce compressed trajectory files in the Gromacs XTC format.
  *
- * It is also needed on BG/L, BG/P and Cray XT3/XT4/XT5 as we don't
- * have XDR support in the lightweight kernel runtimes either.
+ * The XDR sources are available under the BSD 3-clause license for example in
+ * the MIT Kerberos 5 distribution with the following copyright notice and license.
  *
- * This file contains the definitions for Sun External Data
- * Representation (XDR) headers and routines.
+ * @(#)xdr.h    2.2 88/07/29 4.0 RPCSRC
  *
- * Although the rest of LAMPPS is GPL, you can copy and use the XDR
- * routines in any way you want as long as you obey Sun's license:
- * Sun RPC is a product of Sun Microsystems, Inc. and is provided for
- * unrestricted use provided that this legend is included on all tape
- * media and as a part of the software program in whole or part.  Users
- * may copy or modify Sun RPC without charge, but are not authorized
- * to license or distribute it to anyone else except as part of a product or
- * program developed by the user.
+ * Copyright (c) 2010, Oracle America, Inc.
  *
- * SUN RPC IS PROVIDED AS IS WITH NO WARRANTIES OF ANY KIND INCLUDING THE
- * WARRANTIES OF DESIGN, MERCHANTIBILITY AND FITNESS FOR A PARTICULAR
- * PURPOSE, OR ARISING FROM A COURSE OF DEALING, USAGE OR TRADE PRACTICE.
+ * All rights reserved.
  *
- * Sun RPC is provided with no support and without any obligation on the
- * part of Sun Microsystems, Inc. to assist in its use, correction,
- * modification or enhancement.
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions are met:
  *
- * SUN MICROSYSTEMS, INC. SHALL HAVE NO LIABILITY WITH RESPECT TO THE
- * INFRINGEMENT OF COPYRIGHTS, TRADE SECRETS OR ANY PATENTS BY SUN RPC
- * OR ANY PART THEREOF.
+ *     * Redistributions of source code must retain the above copyright
+ *       notice, this list of conditions and the following disclaimer.
  *
- * In no event will Sun Microsystems, Inc. be liable for any lost revenue
- * or profits or other special, indirect and consequential damages, even if
- * Sun has been advised of the possibility of such damages.
+ *     * Redistributions in binary form must reproduce the above copyright
+ *       notice, this list of conditions and the following disclaimer in
+ *       the documentation and/or other materials provided with the
+ *       distribution.
  *
- * Sun Microsystems, Inc.
- * 2550 Garcia Avenue
- * Mountain View, California  94043
+ *     * Neither the name of the "Oracle America, Inc." nor the names of
+ *       its contributors may be used to endorse or promote products
+ *       derived from this software without specific prior written permission.
+ *
+ * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+ * IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+ * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+ * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
+ * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+ * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
+ * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+ * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+ * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+ * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  */
 
+/* compatibility typedefs */
+typedef int bool_t;
+
+typedef int32_t xdr_int32_t;
+typedef uint32_t xdr_uint32_t;
+
 /*
  * Xdr operations.  XDR_ENCODE causes the type to be encoded into the
  * stream.  XDR_DECODE causes the type to be extracted from the stream.
@@ -57,44 +63,8 @@ extern "C" {
  * XDR_DECODE request.
  */
 
-typedef int bool_t;
-
-#if defined(_WIN32) || defined(__APPLE__) || defined(__FreeBSD__) || defined(__DragonFly__) || \
-    defined(__OpenBSD__) || defined(__NetBSD__) || (defined(__linux__) && !defined(__GLIBC_MINOR__))
-typedef char *caddr_t;
-typedef unsigned int u_int;
-#endif
-
-/*
- * Aninteger type that is 32 bits wide. Check if int,
- * long or short is 32 bits and die if none of them is :-)
- */
-#if (INT_MAX == 2147483647)
-typedef int xdr_int32_t;
-typedef unsigned int xdr_uint32_t;
-#elif (LONG_MAX == 2147483647L)
-typedef long xdr_int32_t;
-typedef unsigned long xdr_uint32_t;
-#elif (SHRT_MAX == 2147483647)
-typedef short xdr_int32_t;
-typedef unsigned short xdr_uint32_t;
-#else
-#error ERROR: No 32 bit wide integer type found!
-#endif
-
 enum xdr_op { XDR_ENCODE = 0, XDR_DECODE = 1, XDR_FREE = 2 };
 
-#ifndef FALSE
-#define FALSE (0)
-#endif
-#ifndef TRUE
-#define TRUE (1)
-#endif
-
-#define BYTES_PER_XDR_UNIT (4)
-/* Macro to round up to units of 4. */
-#define XDR_RNDUP(x) (((x) + BYTES_PER_XDR_UNIT - 1) & ~(BYTES_PER_XDR_UNIT - 1))
-
 /*
  * The XDR handle.
  * Contains operation which is being applied to the stream,
@@ -113,26 +83,16 @@ struct XDR {
 };
 
 struct xdr_ops {
+  /* get some bytes from XDR stream */
   bool_t (*x_getbytes)(XDR *__xdrs, char *__addr, unsigned int __len);
-  /* get some bytes from " */
+  /* put some bytes to XDR stream */
   bool_t (*x_putbytes)(XDR *__xdrs, char *__addr, unsigned int __len);
-  /* put some bytes to " */
-  unsigned int (*x_getpostn)(XDR *__xdrs);
-  /* returns bytes off from beginning */
-  bool_t (*x_setpostn)(XDR *__xdrs, unsigned int __pos);
-  /* lets you reposition the stream */
-  xdr_int32_t *(*x_inline)(XDR *__xdrs, int __len);
-  /* buf quick ptr to buffered data */
-  void (*x_destroy)(XDR *__xdrs);
   /* free privates of this xdr_stream */
+  void (*x_destroy)(XDR *__xdrs);
+  /* get a int from XDR stream */
   bool_t (*x_getint32)(XDR *__xdrs, xdr_int32_t *__ip);
-  /* get a int from underlying stream */
+  /* put a int to XDR stream */
   bool_t (*x_putint32)(XDR *__xdrs, xdr_int32_t *__ip);
-  /* put a int to " */
-  bool_t (*x_getuint32)(XDR *__xdrs, xdr_uint32_t *__ip);
-  /* get a unsigned int from underlying stream */
-  bool_t (*x_putuint32)(XDR *__xdrs, xdr_uint32_t *__ip);
-  /* put a int to " */
 };
 
 /*
@@ -151,53 +111,25 @@ typedef bool_t (*xdrproc_t)(XDR *, void *, ...);
  *
  * XDR          *xdrs;
  * xdr_int32_t  *int32p;
- * long         *longp;
- * char         *addr;
  * unsigned int  len;
- * unsigned int  pos;
  */
 
 #define xdr_getint32(xdrs, int32p) (*(xdrs)->x_ops->x_getint32)(xdrs, int32p)
-
 #define xdr_putint32(xdrs, int32p) (*(xdrs)->x_ops->x_putint32)(xdrs, int32p)
-
-#define xdr_getuint32(xdrs, uint32p) (*(xdrs)->x_ops->x_getuint32)(xdrs, uint32p)
-
-#define xdr_putuint32(xdrs, uint32p) (*(xdrs)->x_ops->x_putuint32)(xdrs, uint32p)
-
 #define xdr_getbytes(xdrs, addr, len) (*(xdrs)->x_ops->x_getbytes)(xdrs, addr, len)
-
 #define xdr_putbytes(xdrs, addr, len) (*(xdrs)->x_ops->x_putbytes)(xdrs, addr, len)
-
-#define xdr_getpos(xdrs) (*(xdrs)->x_ops->x_getpostn)(xdrs)
-
-#define xdr_setpos(xdrs, pos) (*(xdrs)->x_ops->x_setpostn)(xdrs, pos)
-
-#define xdr_inline(xdrs, len) (*(xdrs)->x_ops->x_inline)(xdrs, len)
-
 #define xdr_destroy(xdrs)                                            \
   do {                                                               \
     if ((xdrs)->x_ops->x_destroy) (*(xdrs)->x_ops->x_destroy)(xdrs); \
   } while (0)
 
 extern bool_t xdr_int(XDR *__xdrs, int *__ip);
-extern bool_t xdr_u_int(XDR *__xdrs, unsigned int *__ip);
-extern bool_t xdr_short(XDR *__xdrs, short *__ip);
-extern bool_t xdr_u_short(XDR *__xdrs, unsigned short *__ip);
-extern bool_t xdr_bool(XDR *__xdrs, int *__bp);
 extern bool_t xdr_opaque(XDR *__xdrs, char *__cp, unsigned int __cnt);
-extern bool_t xdr_string(XDR *__xdrs, char **__cpp, unsigned int __maxsize);
-extern bool_t xdr_char(XDR *__xdrs, char *__cp);
-extern bool_t xdr_u_char(XDR *__xdrs, unsigned char *__cp);
 extern bool_t xdr_vector(XDR *__xdrs, char *__basep, unsigned int __nelem, unsigned int __elemsize,
                          xdrproc_t __xdr_elem);
 extern bool_t xdr_float(XDR *__xdrs, float *__fp);
-extern bool_t xdr_double(XDR *__xdrs, double *__dp);
 extern void xdrstdio_create(XDR *__xdrs, FILE *__file, enum xdr_op __xop);
 
-/* free memory buffers for xdr */
-extern void xdr_free(xdrproc_t __proc, char *__objp);
-
 #ifdef __cplusplus
 }
 #endif

src/EXTRA-FIX/fix_viscosity.cpp

@@ -120,14 +120,14 @@ int FixViscosity::setmask()
 
 void FixViscosity::init()
 {
-  // warn if any fix ave/spatial comes after this fix
+  // warn if any fix ave/chunk comes after this fix
   // can cause glitch in averaging since ave will happen after swap
 
   int foundme = 0;
-  for (int i = 0; i < modify->nfix; i++) {
-    if (modify->fix[i] == this) foundme = 1;
-    if (foundme && strcmp(modify->fix[i]->style,"ave/spatial") == 0 && me == 0)
-      error->warning(FLERR,"Fix viscosity comes before fix ave/spatial");
+  for (const auto &ifix : modify->get_fix_list()) {
+    if (ifix == this) foundme = 1;
+    if (foundme && utils::strmatch(ifix->style,"^ave/chunk") && (me == 0))
+      error->warning(FLERR,"Fix viscosity comes before fix ave/chunk");
   }
 
   // set bounds of 2 slabs in pdim

src/EXTRA-MOLECULE/angle_cosine_periodic.cpp

@@ -117,7 +117,7 @@ void AngleCosinePeriodic::compute(int eflag, int vflag)
     tn = 1.0;
     tn_1 = 1.0;
     tn_2 = 0.0;
-    un = 1.0;
+    un = (m==1) ? 2.0 : 1.0;
     un_1 = 2.0;
     un_2 = 0.0;
 

src/GPU/pppm_gpu.cpp

@@ -405,7 +405,8 @@ void PPPMGPU::poisson_ik()
 
   // if requested, compute energy and virial contribution
 
-  double scaleinv = 1.0/(nx_pppm*ny_pppm*nz_pppm);
+  bigint ngridtotal = (bigint) nx_pppm * ny_pppm * nz_pppm;
+  double scaleinv = 1.0 / ngridtotal;
   double s2 = scaleinv*scaleinv;
 
   if (eflag_global || vflag_global) {

src/GRANULAR/fix_pour.cpp

@@ -701,7 +701,7 @@ void FixPour::pre_exchange()
   // rebuild atom map
 
   if (atom->map_style != Atom::MAP_NONE) {
-    if (success) atom->map_init();
+    atom->map_init();
     atom->map_set();
   }
 

src/INTEL/fix_intel.cpp

@@ -551,6 +551,9 @@ void FixIntel::kspace_init_check()
 
   if (intel_pair == 0)
     error->all(FLERR,"Intel styles for kspace require intel pair style.");
+
+  if (utils::strmatch(update->integrate_style, "^verlet/split"))
+    error->all(FLERR,"Intel styles for kspace are not compatible with run_style verlet/split");
 }
 
 /* ---------------------------------------------------------------------- */

src/INTEL/fix_nh_intel.cpp

@@ -26,6 +26,7 @@
 #include "memory.h"
 #include "modify.h"
 #include "neighbor.h"
+#include "respa.h"
 #include "update.h"
 
 #include <cmath>
@@ -339,8 +340,13 @@ void FixNHIntel::reset_dt()
 
   // If using respa, then remap is performed in innermost level
 
-  if (utils::strmatch(update->integrate_style,"^respa"))
+  if (utils::strmatch(update->integrate_style,"^respa")) {
+    auto respa_ptr = dynamic_cast<Respa *>(update->integrate);
+    if (!respa_ptr) error->all(FLERR, "Failure to access Respa style {}", update->integrate_style);
+    nlevels_respa = respa_ptr->nlevels;
+    step_respa = respa_ptr->step;
     dto = 0.5*step_respa[0];
+  }
 
   if (pstat_flag)
     pdrag_factor = 1.0 - (update->dt * p_freq_max * drag / nc_pchain);