diff --git a/.github/CODEOWNERS b/.github/CODEOWNERS index 4e8803fc5a..b856b3a8dd 100644 --- a/.github/CODEOWNERS +++ b/.github/CODEOWNERS @@ -13,40 +13,44 @@ lib/kim/* @ellio167 lib/mesont/* @iafoss # whole packages +src/ADIOS/* @pnorbert src/AMOEBA/* @sjplimp -src/COMPRESS/* @rbberger -src/GPU/* @ndtrung81 -src/KOKKOS/* @stanmoore1 -src/KIM/* @ellio167 -src/LATTE/* @cnegre -src/MLIAP/* @athomps -src/SNAP/* @athomps -src/SPIN/* @julient31 +src/BPM/* @jtclemm src/BROWNIAN/* @samueljmcameron src/CG-DNA/* @ohenrich src/CG-SPICA/* @yskmiyazaki src/COLVARS/* @giacomofiorin +src/COMPRESS/* @rbberger src/DIELECTRIC/* @ndtrung81 src/ELECTRODE/* @ludwig-ahrens src/FEP/* @agiliopadua -src/ML-HDNNP/* @singraber +src/GPU/* @ndtrung81 +src/GRANULAR/* @jtclemm @dsbolin src/INTEL/* @wmbrownintel +src/KIM/* @ellio167 +src/KOKKOS/* @stanmoore1 +src/LATTE/* @cnegre src/MANIFOLD/* @Pakketeretet2 -src/MDI/* @taylor-a-barnes +src/MDI/* @taylor-a-barnes @sjplimp src/MEAM/* @martok src/MESONT/* @iafoss +src/ML-HDNNP/* @singraber +src/ML-IAP/* @athomps +src/ML-PACE/* @yury-lysogorskiy src/MOFFF/* @hheenen src/MOLFILE/* @akohlmey src/NETCDF/* @pastewka -src/ML-PACE/* @yury-lysogorskiy -src/PLUMED/* @gtribello -src/PHONON/* @lingtikong -src/PTM/* @pmla src/OPENMP/* @akohlmey +src/PHONON/* @lingtikong +src/PLUGIN/* @akohlmey +src/PLUMED/* @gtribello +src/PTM/* @pmla src/QMMM/* @akohlmey -src/REAXFF/* @hasanmetin @stanmoore1 src/REACTION/* @jrgissing +src/REAXFF/* @hasanmetin @stanmoore1 src/SCAFACOS/* @rhalver +src/SNAP/* @athomps +src/SPIN/* @julient31 src/TALLY/* @akohlmey src/UEF/* @danicholson src/VTK/* @rbberger @@ -120,27 +124,32 @@ src/dump_movie.* @akohlmey src/exceptions.h @rbberger src/fix_nh.* @athomps src/info.* @akohlmey @rbberger -src/timer.* @akohlmey src/min* @sjplimp @stanmoore1 +src/platform.* @akohlmey +src/timer.* @akohlmey src/utils.* @akohlmey @rbberger src/verlet.* @sjplimp @stanmoore1 src/math_eigen_impl.h @jewettaij # tools -tools/msi2lmp/* @akohlmey +tools/coding_standard/* @akohlmey @rbberger tools/emacs/* @HaoZeke -tools/singularity/* @akohlmey @rbberger -tools/coding_standard/* @rbberger -tools/valgrind/* @akohlmey -tools/swig/* @akohlmey +tools/lammps-shell/* @akohlmey +tools/msi2lmp/* @akohlmey tools/offline/* @rbberger +tools/singularity/* @akohlmey @rbberger +tools/swig/* @akohlmey +tools/valgrind/* @akohlmey tools/vim/* @hammondkd # tests -unittest/* @akohlmey @rbberger +unittest/* @akohlmey # cmake cmake/* @junghans @rbberger +cmake/Modules/LAMMPSInterfacePlugin.cmake @akohlmey +cmake/Modules/MPI4WIN.cmake @akohlmey +cmake/Modules/OpenCLLoader.cmake @akohlmey cmake/Modules/Packages/COLVARS.cmake @junghans @rbberger @giacomofiorin cmake/Modules/Packages/KIM.cmake @junghans @rbberger @ellio167 cmake/presets/*.cmake @akohlmey @@ -149,13 +158,12 @@ cmake/presets/*.cmake @akohlmey python/* @rbberger # fortran -fortran/* @akohlmey +fortran/* @akohlmey @hammondkd # docs -doc/utils/*/* @rbberger -doc/Makefile @rbberger -doc/README @rbberger +doc/* @akohlmey examples/plugin/* @akohlmey +examples/PACKAGES/pace/plugin/* @akohlmey # for releases src/version.h @sjplimp diff --git a/.github/workflows/compile-msvc.yml b/.github/workflows/compile-msvc.yml index 3ec5fc4cb0..ee4044b2cd 100644 --- a/.github/workflows/compile-msvc.yml +++ b/.github/workflows/compile-msvc.yml @@ -33,6 +33,8 @@ jobs: run: | python3 -m pip install numpy python3 -m pip install pyyaml + nuget install MSMPIsdk + nuget install MSMPIDIST cmake -C cmake/presets/windows.cmake \ -D PKG_PYTHON=on \ -S cmake -B build \ diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt index c831b46120..c96f5576a5 100644 --- a/cmake/CMakeLists.txt +++ b/cmake/CMakeLists.txt @@ -317,6 +317,15 @@ if(PKG_ADIOS) # script that defines the MPI::MPI_C target enable_language(C) find_package(ADIOS2 REQUIRED) + if(BUILD_MPI) + if(NOT ADIOS2_HAVE_MPI) + message(FATAL_ERROR "ADIOS2 must be built with MPI support when LAMMPS has MPI enabled") + endif() + else() + if(ADIOS2_HAVE_MPI) + message(FATAL_ERROR "ADIOS2 must be built without MPI support when LAMMPS has MPI disabled") + endif() + endif() target_link_libraries(lammps PRIVATE adios2::adios2) endif() diff --git a/cmake/Modules/FindZMQ.cmake b/cmake/Modules/FindZMQ.cmake deleted file mode 100644 index 7d612c2eb3..0000000000 --- a/cmake/Modules/FindZMQ.cmake +++ /dev/null @@ -1,19 +0,0 @@ -find_path(ZMQ_INCLUDE_DIR zmq.h) -find_library(ZMQ_LIBRARY NAMES zmq) - -include(FindPackageHandleStandardArgs) -find_package_handle_standard_args(ZMQ DEFAULT_MSG ZMQ_LIBRARY ZMQ_INCLUDE_DIR) - -# Copy the results to the output variables and target. -if(ZMQ_FOUND) - set(ZMQ_LIBRARIES ${ZMQ_LIBRARY}) - set(ZMQ_INCLUDE_DIRS ${ZMQ_INCLUDE_DIR}) - - if(NOT TARGET ZMQ::ZMQ) - add_library(ZMQ::ZMQ UNKNOWN IMPORTED) - set_target_properties(ZMQ::ZMQ PROPERTIES - IMPORTED_LINK_INTERFACE_LANGUAGES "C" - IMPORTED_LOCATION "${ZMQ_LIBRARY}" - INTERFACE_INCLUDE_DIRECTORIES "${ZMQ_INCLUDE_DIR}") - endif() -endif() diff --git a/cmake/Modules/LAMMPSInterfacePlugin.cmake b/cmake/Modules/LAMMPSInterfacePlugin.cmake index 95bb707e86..e20ae5f5cb 100644 --- a/cmake/Modules/LAMMPSInterfacePlugin.cmake +++ b/cmake/Modules/LAMMPSInterfacePlugin.cmake @@ -112,45 +112,76 @@ if(BUILD_MPI) set(MPI_CXX_SKIP_MPICXX TRUE) # We use a non-standard procedure to cross-compile with MPI on Windows if((CMAKE_SYSTEM_NAME STREQUAL "Windows") AND CMAKE_CROSSCOMPILING) - # Download and configure custom MPICH files for Windows - message(STATUS "Downloading and configuring MPICH-1.4.1 for Windows") - set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win64-devel.tar.gz" CACHE STRING "URL for MPICH2 (win64) tarball") - set(MPICH2_WIN32_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win32-devel.tar.gz" CACHE STRING "URL for MPICH2 (win32) tarball") - set(MPICH2_WIN64_DEVEL_MD5 "4939fdb59d13182fd5dd65211e469f14" CACHE STRING "MD5 checksum of MPICH2 (win64) tarball") - set(MPICH2_WIN32_DEVEL_MD5 "a61d153500dce44e21b755ee7257e031" CACHE STRING "MD5 checksum of MPICH2 (win32) tarball") - mark_as_advanced(MPICH2_WIN64_DEVEL_URL) - mark_as_advanced(MPICH2_WIN32_DEVEL_URL) - mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) - mark_as_advanced(MPICH2_WIN32_DEVEL_MD5) + # Download and configure MinGW compatible MPICH development files for Windows + option(USE_MSMPI "Use Microsoft's MS-MPI SDK instead of MPICH2-1.4.1" OFF) + if(USE_MSMPI) + message(STATUS "Downloading and configuring MS-MPI 10.1 for Windows cross-compilation") + set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/msmpi-win64-devel.tar.gz" CACHE STRING "URL for MS-MPI (win64) tarball") + set(MPICH2_WIN64_DEVEL_MD5 "86314daf1bffb809f1fcbefb8a547490" CACHE STRING "MD5 checksum of MS-MPI (win64) tarball") + mark_as_advanced(MPICH2_WIN64_DEVEL_URL) + mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) - include(ExternalProject) - if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") - ExternalProject_Add(mpi4win_build - URL ${MPICH2_WIN64_DEVEL_URL} - URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} - CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" - BUILD_BYPRODUCTS /lib/libmpi.a) + include(ExternalProject) + if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN64_DEVEL_URL} + URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmsmpi.a) + else() + message(FATAL_ERROR "Only x86 64-bit builds are supported with MS-MPI") + endif() + + ExternalProject_get_property(mpi4win_build SOURCE_DIR) + file(MAKE_DIRECTORY "${SOURCE_DIR}/include") + add_library(MPI::MPI_CXX UNKNOWN IMPORTED) + set_target_properties(MPI::MPI_CXX PROPERTIES + IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmsmpi.a" + INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" + INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + add_dependencies(MPI::MPI_CXX mpi4win_build) + + # set variables for status reporting at the end of CMake run + set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") + set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmsmpi.a") else() - ExternalProject_Add(mpi4win_build - URL ${MPICH2_WIN32_DEVEL_URL} - URL_MD5 ${MPICH2_WIN32_DEVEL_MD5} - CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" - BUILD_BYPRODUCTS /lib/libmpi.a) + # Download and configure custom MPICH files for Windows + message(STATUS "Downloading and configuring MPICH-1.4.1 for Windows") + set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win64-devel.tar.gz" CACHE STRING "URL for MPICH2 (win64) tarball") + set(MPICH2_WIN64_DEVEL_MD5 "4939fdb59d13182fd5dd65211e469f14" CACHE STRING "MD5 checksum of MPICH2 (win64) tarball") + mark_as_advanced(MPICH2_WIN64_DEVEL_URL) + mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) + + include(ExternalProject) + if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN64_DEVEL_URL} + URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmpi.a) + else() + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN32_DEVEL_URL} + URL_MD5 ${MPICH2_WIN32_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmpi.a) + endif() + + ExternalProject_get_property(mpi4win_build SOURCE_DIR) + file(MAKE_DIRECTORY "${SOURCE_DIR}/include") + add_library(MPI::MPI_CXX UNKNOWN IMPORTED) + set_target_properties(MPI::MPI_CXX PROPERTIES + IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmpi.a" + INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" + INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + add_dependencies(MPI::MPI_CXX mpi4win_build) + + # set variables for status reporting at the end of CMake run + set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") + set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmpi.a") endif() - - ExternalProject_get_property(mpi4win_build SOURCE_DIR) - file(MAKE_DIRECTORY "${SOURCE_DIR}/include") - add_library(MPI::MPI_CXX UNKNOWN IMPORTED) - set_target_properties(MPI::MPI_CXX PROPERTIES - IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmpi.a" - INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" - INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") - add_dependencies(MPI::MPI_CXX mpi4win_build) - - # set variables for status reporting at the end of CMake run - set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") - set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") - set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmpi.a") else() find_package(MPI REQUIRED) option(LAMMPS_LONGLONG_TO_LONG "Workaround if your system or MPI version does not recognize 'long long' data types" OFF) diff --git a/cmake/Modules/MPI4WIN.cmake b/cmake/Modules/MPI4WIN.cmake index aa0c9e1833..02db6d4744 100644 --- a/cmake/Modules/MPI4WIN.cmake +++ b/cmake/Modules/MPI4WIN.cmake @@ -1,39 +1,74 @@ -# Download and configure custom MPICH files for Windows -message(STATUS "Downloading and configuring MPICH-1.4.1 for Windows") -set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win64-devel.tar.gz" CACHE STRING "URL for MPICH2 (win64) tarball") -set(MPICH2_WIN32_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win32-devel.tar.gz" CACHE STRING "URL for MPICH2 (win32) tarball") -set(MPICH2_WIN64_DEVEL_MD5 "4939fdb59d13182fd5dd65211e469f14" CACHE STRING "MD5 checksum of MPICH2 (win64) tarball") -set(MPICH2_WIN32_DEVEL_MD5 "a61d153500dce44e21b755ee7257e031" CACHE STRING "MD5 checksum of MPICH2 (win32) tarball") -mark_as_advanced(MPICH2_WIN64_DEVEL_URL) -mark_as_advanced(MPICH2_WIN32_DEVEL_URL) -mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) -mark_as_advanced(MPICH2_WIN32_DEVEL_MD5) +# Download and configure MinGW compatible MPICH development files for Windows +option(USE_MSMPI "Use Microsoft's MS-MPI SDK instead of MPICH2-1.4.1" OFF) -include(ExternalProject) -if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") - ExternalProject_Add(mpi4win_build - URL ${MPICH2_WIN64_DEVEL_URL} - URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} - CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" - BUILD_BYPRODUCTS /lib/libmpi.a) +if(USE_MSMPI) + message(STATUS "Downloading and configuring MS-MPI 10.1 for Windows cross-compilation") + set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/msmpi-win64-devel.tar.gz" CACHE STRING "URL for MS-MPI (win64) tarball") + set(MPICH2_WIN64_DEVEL_MD5 "86314daf1bffb809f1fcbefb8a547490" CACHE STRING "MD5 checksum of MS-MPI (win64) tarball") + mark_as_advanced(MPICH2_WIN64_DEVEL_URL) + mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) + + include(ExternalProject) + if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN64_DEVEL_URL} + URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmsmpi.a) + else() + message(FATAL_ERROR "Only x86 64-bit builds are supported with MS-MPI") + endif() + + ExternalProject_get_property(mpi4win_build SOURCE_DIR) + file(MAKE_DIRECTORY "${SOURCE_DIR}/include") + add_library(MPI::MPI_CXX UNKNOWN IMPORTED) + set_target_properties(MPI::MPI_CXX PROPERTIES + IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmsmpi.a" + INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" + INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + add_dependencies(MPI::MPI_CXX mpi4win_build) + + # set variables for status reporting at the end of CMake run + set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") + set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmsmpi.a") else() - ExternalProject_Add(mpi4win_build - URL ${MPICH2_WIN32_DEVEL_URL} - URL_MD5 ${MPICH2_WIN32_DEVEL_MD5} - CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" - BUILD_BYPRODUCTS /lib/libmpi.a) + message(STATUS "Downloading and configuring MPICH2-1.4.1 for Windows cross-compilation") + set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win64-devel.tar.gz" CACHE STRING "URL for MPICH2 (win64) tarball") + set(MPICH2_WIN32_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win32-devel.tar.gz" CACHE STRING "URL for MPICH2 (win32) tarball") + set(MPICH2_WIN64_DEVEL_MD5 "4939fdb59d13182fd5dd65211e469f14" CACHE STRING "MD5 checksum of MPICH2 (win64) tarball") + set(MPICH2_WIN32_DEVEL_MD5 "a61d153500dce44e21b755ee7257e031" CACHE STRING "MD5 checksum of MPICH2 (win32) tarball") + mark_as_advanced(MPICH2_WIN64_DEVEL_URL) + mark_as_advanced(MPICH2_WIN32_DEVEL_URL) + mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) + mark_as_advanced(MPICH2_WIN32_DEVEL_MD5) + + include(ExternalProject) + if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN64_DEVEL_URL} + URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmpi.a) + else() + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN32_DEVEL_URL} + URL_MD5 ${MPICH2_WIN32_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmpi.a) + endif() + + ExternalProject_get_property(mpi4win_build SOURCE_DIR) + file(MAKE_DIRECTORY "${SOURCE_DIR}/include") + add_library(MPI::MPI_CXX UNKNOWN IMPORTED) + set_target_properties(MPI::MPI_CXX PROPERTIES + IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmpi.a" + INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" + INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + add_dependencies(MPI::MPI_CXX mpi4win_build) + + # set variables for status reporting at the end of CMake run + set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") + set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmpi.a") endif() - -ExternalProject_get_property(mpi4win_build SOURCE_DIR) -file(MAKE_DIRECTORY "${SOURCE_DIR}/include") -add_library(MPI::MPI_CXX UNKNOWN IMPORTED) -set_target_properties(MPI::MPI_CXX PROPERTIES - IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmpi.a" - INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" - INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") -add_dependencies(MPI::MPI_CXX mpi4win_build) - -# set variables for status reporting at the end of CMake run -set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") -set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") -set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmpi.a") diff --git a/doc/src/bond_bpm_rotational.rst b/doc/src/bond_bpm_rotational.rst index 5c43071274..fb1db04b27 100644 --- a/doc/src/bond_bpm_rotational.rst +++ b/doc/src/bond_bpm_rotational.rst @@ -45,6 +45,8 @@ Examples Description """"""""""" +.. versionadded:: 4May2022 + The *bpm/rotational* bond style computes forces and torques based on deviations from the initial reference state of the two atoms. The reference state is stored by each bond when it is first computed in @@ -211,9 +213,9 @@ command, as *b1*, *b2*, ..., *b7*\ . Restrictions """""""""""" -This bond style can only be used if LAMMPS was built with the BPM -package. See the :doc:`Build package ` doc page for -more info. +This bond style is part of the BPM package. It is only enabled if +LAMMPS was built with that package. See the :doc:`Build package +` page for more info. By default if pair interactions are to be disabled, this bond style requires setting diff --git a/doc/src/bond_bpm_spring.rst b/doc/src/bond_bpm_spring.rst index de31357afe..46c300d364 100644 --- a/doc/src/bond_bpm_spring.rst +++ b/doc/src/bond_bpm_spring.rst @@ -45,6 +45,8 @@ Examples Description """"""""""" +.. versionadded:: 4May2022 + The *bpm/spring* bond style computes forces and torques based on deviations from the initial reference state of the two atoms. The reference state is stored by each bond when it is first computed in @@ -167,9 +169,9 @@ extra quantity can be accessed by the Restrictions """""""""""" -This bond style can only be used if LAMMPS was built with the BPM -package. See the :doc:`Build package ` doc page for -more info. +This bond style is part of the BPM package. It is only enabled if +LAMMPS was built with that package. See the :doc:`Build package +` page for more info. By default if pair interactions are to be disabled, this bond style requires setting diff --git a/doc/src/compute_nbond_atom.rst b/doc/src/compute_nbond_atom.rst index 0e53de3e49..f438836534 100644 --- a/doc/src/compute_nbond_atom.rst +++ b/doc/src/compute_nbond_atom.rst @@ -23,6 +23,8 @@ Examples Description """"""""""" +.. versionadded:: 4May2022 + Define a computation that computes the number of bonds each atom is part of. Bonds which are broken are not counted in the tally. See the :doc:`Howto broken bonds ` page for more information. @@ -40,8 +42,9 @@ LAMMPS output options. Restrictions """""""""""" -This fix can only be used if LAMMPS was built with the BPM package. -See the :doc:`Build package ` doc page for more info. +This compute is part of the BPM package. It is only enabled if LAMMPS was +built with that package. See the :doc:`Build package ` +page for more info. Related commands """""""""""""""" diff --git a/doc/src/fix_bond_react.rst b/doc/src/fix_bond_react.rst index c557d2617d..15098761cc 100644 --- a/doc/src/fix_bond_react.rst +++ b/doc/src/fix_bond_react.rst @@ -42,13 +42,16 @@ Syntax * template-ID(post-reacted) = ID of a molecule template containing post-reaction topology * map_file = name of file specifying corresponding atom-IDs in the pre- and post-reacted templates * zero or more individual keyword/value pairs may be appended to each react argument -* individual_keyword = *prob* or *max_rxn* or *stabilize_steps* or *custom_charges* or *molecule* or *modify_create* +* individual_keyword = *prob* or *rate_limit* or *max_rxn* or *stabilize_steps* or *custom_charges* or *rescale_charges* or *molecule* or *modify_create* .. parsed-literal:: *prob* values = fraction seed fraction = initiate reaction with this probability if otherwise eligible seed = random number seed (positive integer) + *rate_limit* = Nlimit Nsteps + Nlimit = maximum number of reactions allowed to occur within interval + Nsteps = the interval (number of timesteps) over which to count reactions *max_rxn* value = N N = maximum number of reactions allowed to occur *stabilize_steps* value = timesteps @@ -56,6 +59,9 @@ Syntax *custom_charges* value = *no* or fragment-ID *no* = update all atomic charges (default) fragment-ID = ID of molecule fragment whose charges are updated + *rescale_charges* value = *no* or *yes* + *no* = do not rescale atomic charges (default) + *yes* = rescale charges such that total charge does not change during reaction *molecule* value = *off* or *inter* or *intra* *off* = allow both inter- and intramolecular reactions (default) *inter* = search for reactions between molecules with different IDs @@ -514,28 +520,40 @@ example, the molecule fragment could consist of only the backbone atoms of a polymer chain. This constraint can be used to enforce a specific relative position and orientation between reacting molecules. +.. versionchanged:: TBD + The constraint of type "custom" has the following syntax: .. parsed-literal:: custom *varstring* -where "custom" is the required keyword, and *varstring* is a -variable expression. The expression must be a valid equal-style -variable formula that can be read by the :doc:`variable ` command, +where 'custom' is the required keyword, and *varstring* is a variable +expression. The expression must be a valid equal-style variable +formula that can be read by the :doc:`variable ` command, after any special reaction functions are evaluated. If the resulting expression is zero, the reaction is prevented from occurring; -otherwise, it is permitted to occur. There are two special reaction -functions available, "rxnsum" and "rxnave". These functions operate -over the atoms in a given reaction site, and have one mandatory -argument and one optional argument. The mandatory argument is the -identifier for an atom-style variable. The second, optional argument -is the name of a molecule fragment in the pre-reaction template, and -can be used to operate over a subset of atoms in the reaction site. -The "rxnsum" function sums the atom-style variable over the reaction -site, while the "rxnave" returns the average value. For example, a -constraint on the total potential energy of atoms involved in the -reaction can be imposed as follows: +otherwise, it is permitted to occur. There are three special reaction +functions available, 'rxnbond', 'rxnsum', and 'rxnave'. The 'rxnbond' +function allows per-bond values to be included in the variable strings +of the custom constraint. The 'rxnbond' function has two mandatory +arguments. The first argument is the ID of a previously defined +'compute bond/local' command. This 'compute bond/local' must compute +only one value, e.g. bond force. This value is returned by the +'rxnbond' function. The second argument is the name of a molecule +fragment in the pre-reaction template. The fragment must contain +exactly two atoms, corresponding to the atoms involved in the bond +whose value should be calculated. An example of a constraint that uses +the force experienced by a bond is provided below. The 'rxnsum' and +'rxnave' functions operate over the atoms in a given reaction site, +and have one mandatory argument and one optional argument. The +mandatory argument is the identifier for an atom-style variable. The +second, optional argument is the name of a molecule fragment in the +pre-reaction template, and can be used to operate over a subset of +atoms in the reaction site. The 'rxnsum' function sums the atom-style +variable over the reaction site, while the 'rxnave' returns the +average value. For example, a constraint on the total potential energy +of atoms involved in the reaction can be imposed as follows: .. code-block:: LAMMPS @@ -547,11 +565,32 @@ reaction can be imposed as follows: custom "rxnsum(v_my_pe) > 100" # in Constraints section of map file The above example prevents the reaction from occurring unless the -total potential energy of the reaction site is above 100. The variable -expression can be interpreted as the probability of the reaction -occurring by using an inequality and the :doc:`random(x,y,z) ` -function available for equal-style variables, similar to the 'arrhenius' -constraint above. +total potential energy of the reaction site is above 100. As a second +example, this time using the 'rxnbond' function, consider a modified +Arrhenius constraint that depends on the bond force of a specific bond: + +.. code-block:: LAMMPS + + # in LAMMPS input script + + compute bondforce all bond/local force + + compute ke_atom all ke/atom + variable ke atom c_ke_atom + + variable E_a equal 100.0 # activation energy + variable l0 equal 1.0 # characteristic length + + +.. code-block:: LAMMPS + + # in Constraints section of map file + + custom "exp(-(v_E_a-rxnbond(c_bondforce,bond1frag)*v_l0)/(2/3*rxnave(v_ke))) < random(0,1,12345)" + +By using an inequality and the 'random(x,y,z)' function, the left-hand +side can be interpreted as the probability of the reaction occurring, +similar to the 'arrhenius' constraint above. By default, all constraints must be satisfied for the reaction to occur. In other words, constraints are evaluated as a series of @@ -598,6 +637,15 @@ eligible reaction only occurs if the random number is less than the fraction. Up to :math:`N` reactions are permitted to occur, as optionally specified by the *max_rxn* keyword. +.. versionadded:: TBD + +The *rate_limit* keyword can enforce an upper limit on the overall +rate of the reaction. The number of reaction occurences is limited to +Nlimit within an interval of Nsteps timesteps. No reactions are +permitted to occur within the first Nsteps timesteps of the first run +after reading a data file. Nlimit can be specified with an equal-style +:doc:`variable `. + The *stabilize_steps* keyword allows for the specification of how many time steps a reaction site is stabilized before being returned to the overall system thermostat. In order to produce the most physical @@ -616,6 +664,19 @@ charges are updated to those specified by the post-reaction template fragment defined in the pre-reaction molecule template. In this case, only the atomic charges of atoms in the molecule fragment are updated. +.. versionadded:: TBD + +The *rescale_charges* keyword can be used to ensure the total charge +of the system does not change as reactions occur. When the argument is +set to *yes*\ , a fixed value is added to the charges of post-reaction +atoms such that their total charge equals that of the pre-reaction +site. If only a subset of atomic charges are updated via the +*custom_charges* keyword, this rescaling is applied to the subset. +This keyword could be useful for systems that contain different +molecules with the same reactive site, if the partial charges on the +reaction site vary from molecule to molecule, or when removing +reaction by-products. + The *molecule* keyword can be used to force the reaction to be intermolecular, intramolecular or either. When the value is set to *off*\ , molecule IDs are not considered when searching for reactions diff --git a/doc/src/fix_nve_bpm_sphere.rst b/doc/src/fix_nve_bpm_sphere.rst index 861586ab2a..ef170605a4 100644 --- a/doc/src/fix_nve_bpm_sphere.rst +++ b/doc/src/fix_nve_bpm_sphere.rst @@ -30,6 +30,8 @@ Examples Description """"""""""" +.. versionadded:: 4May2022 + Perform constant NVE integration to update position, velocity, angular velocity, and quaternion orientation for finite-size spherical particles in the group each timestep. V is volume; E is energy. This @@ -65,6 +67,10 @@ the :doc:`run ` command. This fix is not invoked during Restrictions """""""""""" +This fix is part of the BPM package. It is only enabled if LAMMPS was +built with that package. See the :doc:`Build package ` +page for more info. + This fix requires that atoms store torque, angular velocity (omega), a radius, and a quaternion as defined by the :doc:`atom_style bpm/sphere ` command. diff --git a/doc/src/pair_bpm_spring.rst b/doc/src/pair_bpm_spring.rst index 72e0083bd8..7dfa9bc12c 100644 --- a/doc/src/pair_bpm_spring.rst +++ b/doc/src/pair_bpm_spring.rst @@ -22,6 +22,8 @@ Examples Description """"""""""" +.. versionadded:: 4May2022 + Style *bpm/spring* computes pairwise forces with the formula .. math:: @@ -101,7 +103,11 @@ This pair style can only be used via the *pair* keyword of the Restrictions """""""""""" - none + +This pair style is part of the BPM package. It is only enabled if +LAMMPS was built with that package. See the :doc:`Build package +` page for more info. + Related commands """""""""""""""" diff --git a/doc/src/python.rst b/doc/src/python.rst index aad2f636d3..e5854faf46 100644 --- a/doc/src/python.rst +++ b/doc/src/python.rst @@ -8,14 +8,25 @@ Syntax .. parsed-literal:: - python func keyword args ... + python mode keyword args ... -* func = name of Python function -* one or more keyword/args pairs must be appended +* mode = *source* or name of Python function + + if mode is *source*: .. parsed-literal:: - keyword = *invoke* or *input* or *return* or *format* or *length* or *file* or *here* or *exists* or *source* + keyword = *here* or name of a *Python file* + *here* arg = inline + inline = one or more lines of Python code which defines func + must be a single argument, typically enclosed between triple quotes + *Python file* = name of a file with Python code which will be executed immediately + +* if *mode* is the name of a Python function, one or more keywords with/without arguments must be appended + + .. parsed-literal:: + + keyword = *invoke* or *input* or *return* or *format* or *length* or *file* or *here* or *exists* *invoke* arg = none = invoke the previously defined Python function *input* args = N i1 i2 ... iN N = # of inputs to function @@ -24,7 +35,7 @@ Syntax SELF = reference to LAMMPS itself which can be accessed by Python function variable = v_name, where name = name of LAMMPS variable, e.g. v_abc *return* arg = varReturn - varReturn = v_name = LAMMPS variable name which return value of function will be assigned to + varReturn = v_name = LAMMPS variable name which the return value of the Python function will be assigned to *format* arg = fstring with M characters M = N if no return value, where N = # of inputs M = N+1 if there is a return value @@ -38,10 +49,6 @@ Syntax inline = one or more lines of Python code which defines func must be a single argument, typically enclosed between triple quotes *exists* arg = none = Python code has been loaded by previous python command - *source* arg = *filename* or *inline* - filename = file of Python code which will be executed immediately - inline = one or more lines of Python code which will be executed immediately - must be a single argument, typically enclosed between triple quotes Examples """""""" @@ -70,80 +77,103 @@ Examples lmp.command("pair_style lj/cut ${cut}") # LAMMPS commands lmp.command("pair_coeff * * 1.0 1.0") lmp.command("run 100") - """ + """ + + python source funcdef.py + + python source here "from lammps import lammps" + Description """"""""""" -Define a Python function or execute a previously defined function or -execute some arbitrary python code. +The *python* command allows interfacing LAMMPS with an embedded Python +interpreter and enables either executing arbitrary python code in that +interpreter, registering a Python function for future execution (as a +python style variable, from a fix interfaced with python, or for direct +invocation), or invoking such a previously registered function. + Arguments, including LAMMPS variables, can be passed to the function -from the LAMMPS input script and a value returned by the Python -function to a LAMMPS variable. The Python code for the function can -be included directly in the input script or in a separate Python file. -The function can be standard Python code or it can make "callbacks" to -LAMMPS through its library interface to query or set internal values -within LAMMPS. This is a powerful mechanism for performing complex -operations in a LAMMPS input script that are not possible with the -simple input script and variable syntax which LAMMPS defines. Thus -your input script can operate more like a true programming language. +from the LAMMPS input script and a value returned by the Python function +assigned to a LAMMPS variable. The Python code for the function can be included +directly in the input script or in a separate Python file. The function +can be standard Python code or it can make "callbacks" to LAMMPS through +its library interface to query or set internal values within LAMMPS. +This is a powerful mechanism for performing complex operations in a +LAMMPS input script that are not possible with the simple input script +and variable syntax which LAMMPS defines. Thus your input script can +operate more like a true programming language. Use of this command requires building LAMMPS with the PYTHON package which links to the Python library so that the Python interpreter is embedded in LAMMPS. More details about this process are given below. There are two ways to invoke a Python function once it has been -defined. One is using the *invoke* keyword. The other is to assign +registered. One is using the *invoke* keyword. The other is to assign the function to a :doc:`python-style variable ` defined in -your input script. Whenever the variable is evaluated, it will -execute the Python function to assign a value to the variable. Note -that variables can be evaluated in many different ways within LAMMPS. -They can be substituted for directly in an input script. Or they can -be passed to various commands as arguments, so that the variable is -evaluated during a simulation run. +your input script. Whenever the variable is evaluated, it will execute +the Python function to assign a value to the variable. Note that +variables can be evaluated in many different ways within LAMMPS. They +can be substituted with their result directly in an input script, or +they can be passed to various commands as arguments, so that the +variable is evaluated during a simulation run. -A broader overview of how Python can be used with LAMMPS is given on -the :doc:`Python ` doc page. There is an examples/python -directory which illustrates use of the python command. +A broader overview of how Python can be used with LAMMPS is given in the +:doc:`Use Python with LAMMPS ` section of the +documentation. There also is an ``examples/python`` directory which +illustrates use of the python command. ---------- -The *func* setting specifies the name of the Python function. The -code for the function is defined using the *file* or *here* keywords -as explained below. In case of the *source* keyword, the name of -the function is ignored. +The first argument of the *python* command is either the *source* +keyword or the name of a Python function. This defines the mode +of the python command. + +If the *source* keyword is used, it is followed by either a file name or +the *here* keyword. No other keywords can be used. The *here* keyword +is followed by a string with python commands, either on a single line +enclosed in quotes, or as multiple lines enclosed in triple quotes. +These Python commands will be passed to the python interpreter and +executed immediately without registering a Python function for future +execution. The code will be loaded into and run in the "main" module of +the Python interpreter. This allows running arbitrary Python code at +any time while processing the LAMMPS input file. This can be used to +pre-load Python modules, initialize global variables, define functions +or classes, or perform operations using the python programming language. +The Python code will be executed in parallel on all MPI processes. No +arguments can be passed. + +In all other cases, the first argument is the name of a Python function +that will be registered with LAMMPS for future execution. The function +may already be defined (see *exists* keyword) or must be defined using +the *file* or *here* keywords as explained below. If the *invoke* keyword is used, no other keywords can be used, and a -previous python command must have defined the Python function +previous *python* command must have registered the Python function referenced by this command. This invokes the Python function with the -previously defined arguments and return value processed as explained -below. You can invoke the function as many times as you wish in your -input script. - -If the *source* keyword is used, no other keywords can be used. -The argument can be a filename or a string with python commands, -either on a single line enclosed in quotes, or as multiple lines -enclosed in triple quotes. These python commands will be passed -to the python interpreter and executed immediately without registering -a python function for future execution. +previously defined arguments and the return value is processed as +explained below. You can invoke the function as many times as you wish +in your input script. The *input* keyword defines how many arguments *N* the Python function -expects. If it takes no arguments, then the *input* keyword should -not be used. Each argument can be specified directly as a value, -e.g. 6 or 3.14159 or abc (a string of characters). The type of each +expects. If it takes no arguments, then the *input* keyword should not +be used. Each argument can be specified directly as a value, e.g. '6' +or '3.14159' or 'abc' (a string of characters). The type of each argument is specified by the *format* keyword as explained below, so that Python will know how to interpret the value. If the word SELF is used for an argument it has a special meaning. A pointer is passed to -the Python function which it converts into a reference to LAMMPS -itself. This enables the function to call back to LAMMPS through its -library interface as explained below. This allows the Python function -to query or set values internal to LAMMPS which can affect the -subsequent execution of the input script. A LAMMPS variable can also -be used as an argument, specified as v_name, where "name" is the name -of the variable. Any style of LAMMPS variable can be used, as defined -by the :doc:`variable ` command. Each time the Python -function is invoked, the LAMMPS variable is evaluated and its value is -passed to the Python function. +the Python function which it can convert into a reference to LAMMPS +itself using the :doc:`LAMMPS Python module `. This +enables the function to call back to LAMMPS through its library +interface as explained below. This allows the Python function to query +or set values internal to LAMMPS which can affect the subsequent +execution of the input script. A LAMMPS variable can also be used as an +argument, specified as v_name, where "name" is the name of the variable. +Any style of LAMMPS variable returning a scalar or a string can be used, +as defined by the :doc:`variable ` command. The *format* +keyword must be used to set the type of data that is passed to Python. +Each time the Python function is invoked, the LAMMPS variable is +evaluated and its value is passed to the Python function. The *return* keyword is only needed if the Python function returns a value. The specified *varReturn* must be of the form v_name, where @@ -153,8 +183,9 @@ numeric or string value, as specified by the *format* keyword. As explained on the :doc:`variable ` doc page, the definition of a python-style variable associates a Python function name with the -variable. This must match the *func* setting for this command. For -example these two commands would be self-consistent: +variable. This must match the *Python function name* first argument of +the *python* command. For example these two commands would be +consistent: .. code-block:: LAMMPS @@ -163,21 +194,22 @@ example these two commands would be self-consistent: The two commands can appear in either order in the input script so long as both are specified before the Python function is invoked for -the first time. +the first time. Afterwards, the variable 'foo' is associated with +the Python function 'myMultiply'. -The *format* keyword must be used if the *input* or *return* keyword -is used. It defines an *fstring* with M characters, where M = sum of +The *format* keyword must be used if the *input* or *return* keywords +are used. It defines an *fstring* with M characters, where M = sum of number of inputs and outputs. The order of characters corresponds to the N inputs, followed by the return value (if it exists). Each character must be one of the following: "i" for integer, "f" for -floating point, "s" for string, or "p" for SELF. Each character -defines the type of the corresponding input or output value of the -Python function and affects the type conversion that is performed -internally as data is passed back and forth between LAMMPS and Python. -Note that it is permissible to use a :doc:`python-style variable ` in a LAMMPS command that allows for an -equal-style variable as an argument, but only if the output of the -Python function is flagged as a numeric value ("i" or "f") via the -*format* keyword. +floating point, "s" for string, or "p" for SELF. Each character defines +the type of the corresponding input or output value of the Python +function and affects the type conversion that is performed internally as +data is passed back and forth between LAMMPS and Python. Note that it +is permissible to use a :doc:`python-style variable ` in a +LAMMPS command that allows for an equal-style variable as an argument, +but only if the output of the Python function is flagged as a numeric +value ("i" or "f") via the *format* keyword. If the *return* keyword is used and the *format* keyword specifies the output as a string, then the default maximum length of that string is @@ -192,12 +224,13 @@ truncated. Either the *file*, *here*, or *exists* keyword must be used, but only one of them. These keywords specify what Python code to load into the -Python interpreter. The *file* keyword gives the name of a file, -which should end with a ".py" suffix, which contains Python code. The -code will be immediately loaded into and run in the "main" module of -the Python interpreter. Note that Python code which contains a -function definition does not "execute" the function when it is run; it -simply defines the function so that it can be invoked later. +Python interpreter. The *file* keyword gives the name of a file +containing Python code, which should end with a ".py" suffix. The code +will be immediately loaded into and run in the "main" module of the +Python interpreter. The Python code will be executed in parallel on all +MPI processes. Note that Python code which contains a function +definition does not "execute" the function when it is run; it simply +defines the function so that it can be invoked later. The *here* keyword does the same thing, except that the Python code follows as a single argument to the *here* keyword. This can be done @@ -208,14 +241,15 @@ proper indentation, blank lines, and comments, as desired. See the how triple quotes can be used as part of input script syntax. The *exists* keyword takes no argument. It means that Python code -containing the required Python function defined by the *func* setting, -is assumed to have been previously loaded by another python command. +containing the required Python function with the given name has already +been executed, for example by a *python source* command or in the same +file that was used previously with the *file* keyword. -Note that the Python code that is loaded and run must contain a -function with the specified *func* name. To operate properly when -later invoked, the function code must match the *input* and -*return* and *format* keywords specified by the python command. -Otherwise Python will generate an error. +Note that the Python code that is loaded and run must contain a function +with the specified function name. To operate properly when later +invoked, the function code must match the *input* and *return* and +*format* keywords specified by the python command. Otherwise Python +will generate an error. ---------- @@ -225,19 +259,19 @@ LAMMPS. Whether you load Python code from a file or directly from your input script, via the *file* and *here* keywords, the code can be identical. It must be indented properly as Python requires. It can contain -comments or blank lines. If the code is in your input script, it -cannot however contain triple-quoted Python strings, since that will -conflict with the triple-quote parsing that the LAMMPS input script -performs. +comments or blank lines. If the code is in your input script, it cannot +however contain triple-quoted Python strings, since that will conflict +with the triple-quote parsing that the LAMMPS input script performs. All the Python code you specify via one or more python commands is -loaded into the Python "main" module, i.e. __main__. The code can -define global variables or statements that are outside of function -definitions. It can contain multiple functions, only one of which -matches the *func* setting in the python command. This means you can -use the *file* keyword once to load several functions, and the -*exists* keyword thereafter in subsequent python commands to access -the other functions previously loaded. +loaded into the Python "main" module, i.e. ``__name__ == '__main__'``. +The code can define global variables, define global functions, define +classes or execute statements that are outside of function definitions. +It can contain multiple functions, only one of which matches the *func* +setting in the python command. This means you can use the *file* +keyword once to load several functions, and the *exists* keyword +thereafter in subsequent python commands to register the other functions +that were previously loaded with LAMMPS. A Python function you define (or more generally, the code you load) can import other Python modules or classes, it can make calls to other @@ -264,12 +298,13 @@ outside the function: nvaluelast = nvalue return nvalue -Nsteplast stores the previous timestep the function was invoked -(passed as an argument to the function). Nvaluelast stores the return -value computed on the last function invocation. If the function is -invoked again on the same timestep, the previous value is simply -returned, without re-computing it. The "global" statement inside the -Python function allows it to overwrite the global variables. +The variable 'nsteplast' stores the previous timestep the function was +invoked (passed as an argument to the function). The variable +'nvaluelast' stores the return value computed on the last function +invocation. If the function is invoked again on the same timestep, the +previous value is simply returned, without re-computing it. The +"global" statement inside the Python function allows it to overwrite the +global variables from within the local context of the function. Note that if you load Python code multiple times (via multiple python commands), you can overwrite previously loaded variables and functions @@ -285,19 +320,39 @@ copy of the Python function(s) you define. There is no connection between the Python interpreters running on different processors. This implies three important things. -First, if you put a print statement in your Python function, you will -see P copies of the output, when running on P processors. If the -prints occur at (nearly) the same time, the P copies of the output may -be mixed together. Welcome to the world of parallel programming and -debugging. +First, if you put a print or other statement creating output to the +screen in your Python function, you will see P copies of the output, +when running on P processors. If the prints occur at (nearly) the same +time, the P copies of the output may be mixed together. When loading +the LAMMPS Python module into the embedded Python interpreter, it is +possible to pass the pointer to the current LAMMPS class instance and +via the Python interface to the LAMMPS library interface, it is possible +to determine the MPI rank of the current process and thus adapt the +Python code so that output will only appear on MPI rank 0. The +following LAMMPS input demonstrates how this could be done. The text +'Hello, LAMMPS!' should be printed only once, even when running LAMMPS +in parallel. -Second, if your Python code loads modules that are not pre-loaded by -the Python library, then it will load the module from disk. This may -be a bottleneck if 1000s of processors try to load a module at the -same time. On some large supercomputers, loading of modules from disk -by Python may be disabled. In this case you would need to pre-build a -Python library that has the required modules pre-loaded and link -LAMMPS with that library. +.. code-block:: LAMMPS + + python python_hello input 1 SELF format p here """ + def python_hello(handle): + from lammps import lammps + lmp = lammps(ptr=handle) + me = lmp.extract_setting('world_rank') + if me == 0: + print('Hello, LAMMPS!') + """ + + python python_hello invoke + +If your Python code loads Python modules that are not pre-loaded by the +Python library, then it will load the module from disk. This may be a +bottleneck if 1000s of processors try to load a module at the same time. +On some large supercomputers, loading of modules from disk by Python may +be disabled. In this case you would need to pre-build a Python library +that has the required modules pre-loaded and link LAMMPS with that +library. Third, if your Python code calls back to LAMMPS (discussed in the next section) and causes LAMMPS to perform an MPI operation requires @@ -315,22 +370,21 @@ Python function is as follows: .. code-block:: python - def foo(lmpptr,...): + def foo(handle,...): from lammps import lammps - lmp = lammps(ptr=lmpptr) + lmp = lammps(ptr=handle) lmp.command('print "Hello from inside Python"') ... -The function definition must include a variable (lmpptr in this case) -which corresponds to SELF in the python command. The first line of the -function imports the :doc:`"lammps" Python module `. -The second line creates a Python object ``lmp`` which -wraps the instance of LAMMPS that called the function. The "ptr=lmpptr" -argument is what makes that happen. The third line invokes the -command() function in the LAMMPS library interface. It takes a single -string argument which is a LAMMPS input script command for LAMMPS to -execute, the same as if it appeared in your input script. In this case, -LAMMPS should output +The function definition must include a variable ('handle' in this case) +which corresponds to SELF in the *python* command. The first line of +the function imports the :doc:`"lammps" Python module `. +The second line creates a Python object ``lmp`` which wraps the instance +of LAMMPS that called the function. The 'ptr=handle' argument is what +makes that happen. The third line invokes the command() function in the +LAMMPS library interface. It takes a single string argument which is a +LAMMPS input script command for LAMMPS to execute, the same as if it +appeared in your input script. In this case, LAMMPS should output .. parsed-literal:: @@ -344,8 +398,8 @@ The :doc:`Python_head` page describes the syntax for how Python wraps the various functions included in the LAMMPS library interface. -A more interesting example is in the examples/python/in.python script -which loads and runs the following function from examples/python/funcs.py: +A more interesting example is in the ``examples/python/in.python`` script +which loads and runs the following function from ``examples/python/funcs.py``: .. code-block:: python @@ -495,24 +549,35 @@ Restrictions """""""""""" This command is part of the PYTHON package. It is only enabled if -LAMMPS was built with that package. See the :doc:`Build package ` page for more info. +LAMMPS was built with that package. See the :doc:`Build package +` page for more info. -Building LAMMPS with the PYTHON package will link LAMMPS with the -Python library on your system. Settings to enable this are in the +Building LAMMPS with the PYTHON package will link LAMMPS with the Python +library on your system. Settings to enable this are in the lib/python/Makefile.lammps file. See the lib/python/README file for information on those settings. -If you use Python code which calls back to LAMMPS, via the SELF input argument -explained above, there is an extra step required when building LAMMPS. LAMMPS -must also be built as a shared library and your Python function must be able to -load the :doc:`"lammps" Python module ` that wraps the LAMMPS -library interface. These are the same steps required to use Python by itself -to wrap LAMMPS. Details on these steps are explained on the :doc:`Python -` doc page. Note that it is important that the stand-alone LAMMPS -executable and the LAMMPS shared library be consistent (built from the same -source code files) in order for this to work. If the two have been built at +If you use Python code which calls back to LAMMPS, via the SELF input +argument explained above, there is an extra step required when building +LAMMPS. LAMMPS must also be built as a shared library and your Python +function must be able to load the :doc:`"lammps" Python module +` that wraps the LAMMPS library interface. These are the +same steps required to use Python by itself to wrap LAMMPS. Details on +these steps are explained on the :doc:`Python ` doc page. +Note that it is important that the stand-alone LAMMPS executable and the +LAMMPS shared library be consistent (built from the same source code +files) in order for this to work. If the two have been built at different times using different source files, problems may occur. +Another limitation of calling back to Python from the LAMMPS module +using the *python* command in a LAMMPS input is that both, the Python +interpreter and LAMMPS, must be linked to the same Python runtime as a +shared library. If the Python interpreter is linked to Python +statically (which seems to happen with Conda) then loading the shared +LAMMPS library will create a second python "main" module that hides the +one from the Python interpreter and all previous defined function and +global variables will become invisible. + Related commands """""""""""""""" diff --git a/doc/utils/sphinx-config/false_positives.txt b/doc/utils/sphinx-config/false_positives.txt index 61c1d45ba7..13a3682870 100644 --- a/doc/utils/sphinx-config/false_positives.txt +++ b/doc/utils/sphinx-config/false_positives.txt @@ -2248,6 +2248,7 @@ MxN myCompute myIndex mylammps +myMultiply MyPool mysocket mySpin @@ -2473,7 +2474,7 @@ nsq Nstart nstats Nstep -Nsteplast +nsteplast Nstop nsub Nsw @@ -2503,7 +2504,7 @@ numpy Numpy Nurdin Nvalue -Nvaluelast +nvaluelast Nvalues nvc nvcc diff --git a/examples/PACKAGES/pace/plugin/CMakeLists.txt b/examples/PACKAGES/pace/plugin/CMakeLists.txt index 6ad9c791ba..25fa877ffc 100644 --- a/examples/PACKAGES/pace/plugin/CMakeLists.txt +++ b/examples/PACKAGES/pace/plugin/CMakeLists.txt @@ -38,9 +38,15 @@ elseif(CMAKE_SYSTEM_NAME STREQUAL "Windows") execute_process(COMMAND ${CMAKE_COMMAND} -E copy_if_different ${CMAKE_SOURCE_DIR}/lammps.ico ${CMAKE_SOURCE_DIR}/lammps-text-logo-wide.bmp ${CMAKE_SOURCE_DIR}/paceplugin.nsis ${CMAKE_BINARY_DIR}) if(BUILD_MPI) - add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION}-MPI paceplugin.nsis - DEPENDS paceplugin - BYPRODUCTS LAMMPS-ML-PACE-plugin-${LAMMPS_VERSION}-MPI.exe) + if(USE_MSMPI) + add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION}-MSMPI paceplugin.nsis + DEPENDS paceplugin + BYPRODUCTS LAMMPS-ML-PACE-plugin-${LAMMPS_VERSION}-MSMPI.exe) + else() + add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION}-MPI paceplugin.nsis + DEPENDS paceplugin + BYPRODUCTS LAMMPS-ML-PACE-plugin-${LAMMPS_VERSION}-MPI.exe) + endif() else() add_custom_target(package ${MAKENSIS_PATH} -V1 -DVERSION=${LAMMPS_VERSION} paceplugin.nsis COMMAND ${CMAKE_COMMAND} -E echo ${PWD} diff --git a/examples/plugins/LAMMPSInterfaceCXX.cmake b/examples/plugins/LAMMPSInterfaceCXX.cmake index a3313215d6..7eef5bd6e4 100644 --- a/examples/plugins/LAMMPSInterfaceCXX.cmake +++ b/examples/plugins/LAMMPSInterfaceCXX.cmake @@ -45,45 +45,80 @@ if(BUILD_MPI) set(MPI_CXX_SKIP_MPICXX TRUE) # We use a non-standard procedure to cross-compile with MPI on Windows if((CMAKE_SYSTEM_NAME STREQUAL "Windows") AND CMAKE_CROSSCOMPILING) - # Download and configure custom MPICH files for Windows - message(STATUS "Downloading and configuring MPICH-1.4.1 for Windows") - set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win64-devel.tar.gz" CACHE STRING "URL for MPICH2 (win64) tarball") - set(MPICH2_WIN32_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win32-devel.tar.gz" CACHE STRING "URL for MPICH2 (win32) tarball") - set(MPICH2_WIN64_DEVEL_MD5 "4939fdb59d13182fd5dd65211e469f14" CACHE STRING "MD5 checksum of MPICH2 (win64) tarball") - set(MPICH2_WIN32_DEVEL_MD5 "a61d153500dce44e21b755ee7257e031" CACHE STRING "MD5 checksum of MPICH2 (win32) tarball") - mark_as_advanced(MPICH2_WIN64_DEVEL_URL) - mark_as_advanced(MPICH2_WIN32_DEVEL_URL) - mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) - mark_as_advanced(MPICH2_WIN32_DEVEL_MD5) + # Download and configure MinGW compatible MPICH development files for Windows + option(USE_MSMPI "Use Microsoft's MS-MPI SDK instead of MPICH2-1.4.1" OFF) + if(USE_MSMPI) + message(STATUS "Downloading and configuring MS-MPI 10.1 for Windows cross-compilation") + set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/msmpi-win64-devel.tar.gz" CACHE STRING "URL for MS-MPI (win64) tarball") + set(MPICH2_WIN64_DEVEL_MD5 "86314daf1bffb809f1fcbefb8a547490" CACHE STRING "MD5 checksum of MS-MPI (win64) tarball") + mark_as_advanced(MPICH2_WIN64_DEVEL_URL) + mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) - include(ExternalProject) - if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") - ExternalProject_Add(mpi4win_build - URL ${MPICH2_WIN64_DEVEL_URL} - URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} - CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" - BUILD_BYPRODUCTS /lib/libmpi.a) + include(ExternalProject) + if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN64_DEVEL_URL} + URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmsmpi.a) + else() + message(FATAL_ERROR "Only x86 64-bit builds are supported with MS-MPI") + endif() + + ExternalProject_get_property(mpi4win_build SOURCE_DIR) + file(MAKE_DIRECTORY "${SOURCE_DIR}/include") + add_library(MPI::MPI_CXX UNKNOWN IMPORTED) + set_target_properties(MPI::MPI_CXX PROPERTIES + IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmsmpi.a" + INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" + INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + add_dependencies(MPI::MPI_CXX mpi4win_build) + + # set variables for status reporting at the end of CMake run + set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") + set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmsmpi.a") else() - ExternalProject_Add(mpi4win_build - URL ${MPICH2_WIN32_DEVEL_URL} - URL_MD5 ${MPICH2_WIN32_DEVEL_MD5} - CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" - BUILD_BYPRODUCTS /lib/libmpi.a) + # Download and configure custom MPICH files for Windows + message(STATUS "Downloading and configuring MPICH-1.4.1 for Windows") + set(MPICH2_WIN64_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win64-devel.tar.gz" CACHE STRING "URL for MPICH2 (win64) tarball") + set(MPICH2_WIN32_DEVEL_URL "${LAMMPS_THIRDPARTY_URL}/mpich2-win32-devel.tar.gz" CACHE STRING "URL for MPICH2 (win32) tarball") + set(MPICH2_WIN64_DEVEL_MD5 "4939fdb59d13182fd5dd65211e469f14" CACHE STRING "MD5 checksum of MPICH2 (win64) tarball") + set(MPICH2_WIN32_DEVEL_MD5 "a61d153500dce44e21b755ee7257e031" CACHE STRING "MD5 checksum of MPICH2 (win32) tarball") + mark_as_advanced(MPICH2_WIN64_DEVEL_URL) + mark_as_advanced(MPICH2_WIN32_DEVEL_URL) + mark_as_advanced(MPICH2_WIN64_DEVEL_MD5) + mark_as_advanced(MPICH2_WIN32_DEVEL_MD5) + + include(ExternalProject) + if(CMAKE_SYSTEM_PROCESSOR STREQUAL "x86_64") + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN64_DEVEL_URL} + URL_MD5 ${MPICH2_WIN64_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmpi.a) + else() + ExternalProject_Add(mpi4win_build + URL ${MPICH2_WIN32_DEVEL_URL} + URL_MD5 ${MPICH2_WIN32_DEVEL_MD5} + CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND "" + BUILD_BYPRODUCTS /lib/libmpi.a) + endif() + + ExternalProject_get_property(mpi4win_build SOURCE_DIR) + file(MAKE_DIRECTORY "${SOURCE_DIR}/include") + add_library(MPI::MPI_CXX UNKNOWN IMPORTED) + set_target_properties(MPI::MPI_CXX PROPERTIES + IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmpi.a" + INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" + INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + add_dependencies(MPI::MPI_CXX mpi4win_build) + + # set variables for status reporting at the end of CMake run + set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") + set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") + set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmpi.a") endif() - - ExternalProject_get_property(mpi4win_build SOURCE_DIR) - file(MAKE_DIRECTORY "${SOURCE_DIR}/include") - add_library(MPI::MPI_CXX UNKNOWN IMPORTED) - set_target_properties(MPI::MPI_CXX PROPERTIES - IMPORTED_LOCATION "${SOURCE_DIR}/lib/libmpi.a" - INTERFACE_INCLUDE_DIRECTORIES "${SOURCE_DIR}/include" - INTERFACE_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") - add_dependencies(MPI::MPI_CXX mpi4win_build) - - # set variables for status reporting at the end of CMake run - set(MPI_CXX_INCLUDE_PATH "${SOURCE_DIR}/include") - set(MPI_CXX_COMPILE_DEFINITIONS "MPICH_SKIP_MPICXX") - set(MPI_CXX_LIBRARIES "${SOURCE_DIR}/lib/libmpi.a") else() find_package(MPI REQUIRED) option(LAMMPS_LONGLONG_TO_LONG "Workaround if your system or MPI version does not recognize 'long long' data types" OFF) diff --git a/src/ADIOS/dump_atom_adios.cpp b/src/ADIOS/dump_atom_adios.cpp index a31f3f4f25..8e92300591 100644 --- a/src/ADIOS/dump_atom_adios.cpp +++ b/src/ADIOS/dump_atom_adios.cpp @@ -62,7 +62,11 @@ DumpAtomADIOS::DumpAtomADIOS(LAMMPS *lmp, int narg, char **arg) : DumpAtom(lmp, internal = new DumpAtomADIOSInternal(); try { +#if defined(MPI_STUBS) + internal->ad = new adios2::ADIOS("adios2_config.xml", adios2::DebugON); +#else internal->ad = new adios2::ADIOS("adios2_config.xml", world, adios2::DebugON); +#endif } catch (std::ios_base::failure &e) { error->all(FLERR, "ADIOS initialization failed with error: {}", e.what()); } @@ -84,11 +88,19 @@ void DumpAtomADIOS::openfile() if (multifile) { // if one file per timestep, replace '*' with current timestep auto filecurrent = utils::star_subst(filename, update->ntimestep, padflag); +#if defined(MPI_STUBS) + internal->fh = internal->io.Open(filecurrent, adios2::Mode::Write); +#else internal->fh = internal->io.Open(filecurrent, adios2::Mode::Write, world); +#endif if (!internal->fh) error->one(FLERR, "Cannot open dump file {}", filecurrent); } else { if (!singlefile_opened) { +#if defined(MPI_STUBS) + internal->fh = internal->io.Open(filename, adios2::Mode::Write); +#else internal->fh = internal->io.Open(filename, adios2::Mode::Write, world); +#endif if (!internal->fh) error->one(FLERR, "Cannot open dump file {}", filename); singlefile_opened = 1; } @@ -256,61 +268,63 @@ void DumpAtomADIOS::init_style() else if (scale_flag == 0 && image_flag == 1) pack_choice = &DumpAtomADIOS::pack_noscale_image; - /* Define the group of variables for the atom style here since it's a fixed - * set */ - internal->io = internal->ad->DeclareIO(internal->ioName); - if (!internal->io.InConfigFile()) { - // if not defined by user, we can change the default settings - // BPFile is the default writer - internal->io.SetEngine("BPFile"); - int num_aggregators = multiproc; - if (num_aggregators == 0) num_aggregators = 1; - auto nstreams = std::to_string(num_aggregators); - internal->io.SetParameters({{"substreams", nstreams}}); - if (me == 0) - utils::logmesg(lmp, "ADIOS method for {} is n-to-m (aggregation with {} writers)\n", filename, - nstreams); - } + /* Define the group of variables for the atom style here since it's a fixed set */ - internal->io.DefineVariable("ntimestep"); - internal->io.DefineVariable("natoms"); + if (!internal->io) { + internal->io = internal->ad->DeclareIO(internal->ioName); + if (!internal->io.InConfigFile()) { + // if not defined by user, we can change the default settings + // BPFile is the default writer + internal->io.SetEngine("BPFile"); + int num_aggregators = multiproc; + if (num_aggregators == 0) num_aggregators = 1; + auto nstreams = std::to_string(num_aggregators); + internal->io.SetParameters({{"substreams", nstreams}}); + if (me == 0) + utils::logmesg(lmp, "ADIOS method for {} is n-to-m (aggregation with {} writers)\n", filename, + nstreams); + } - internal->io.DefineVariable("nprocs"); - internal->io.DefineVariable("ncolumns"); + internal->io.DefineVariable("ntimestep"); + internal->io.DefineVariable("natoms"); - internal->io.DefineVariable("boxxlo"); - internal->io.DefineVariable("boxxhi"); - internal->io.DefineVariable("boxylo"); - internal->io.DefineVariable("boxyhi"); - internal->io.DefineVariable("boxzlo"); - internal->io.DefineVariable("boxzhi"); + internal->io.DefineVariable("nprocs"); + internal->io.DefineVariable("ncolumns"); - internal->io.DefineVariable("boxxy"); - internal->io.DefineVariable("boxxz"); - internal->io.DefineVariable("boxyz"); + internal->io.DefineVariable("boxxlo"); + internal->io.DefineVariable("boxxhi"); + internal->io.DefineVariable("boxylo"); + internal->io.DefineVariable("boxyhi"); + internal->io.DefineVariable("boxzlo"); + internal->io.DefineVariable("boxzhi"); - internal->io.DefineAttribute("triclinic", domain->triclinic); - internal->io.DefineAttribute("scaled", scale_flag); - internal->io.DefineAttribute("image", image_flag); + internal->io.DefineVariable("boxxy"); + internal->io.DefineVariable("boxxz"); + internal->io.DefineVariable("boxyz"); - int *boundaryptr = reinterpret_cast(domain->boundary); - internal->io.DefineAttribute("boundary", boundaryptr, 6); + internal->io.DefineAttribute("triclinic", domain->triclinic); + internal->io.DefineAttribute("scaled", scale_flag); + internal->io.DefineAttribute("image", image_flag); - auto nColumns = static_cast(size_one); - internal->io.DefineAttribute("columns", columnNames.data(), nColumns); - internal->io.DefineAttribute("columnstr", columns); - internal->io.DefineAttribute("boundarystr", boundstr); - internal->io.DefineAttribute("LAMMPS/dump_style", "atom"); - internal->io.DefineAttribute("LAMMPS/version", lmp->version); - internal->io.DefineAttribute("LAMMPS/num_ver", std::to_string(lmp->num_ver)); + int *boundaryptr = reinterpret_cast(domain->boundary); + internal->io.DefineAttribute("boundary", boundaryptr, 6); - // local dimension variables - internal->io.DefineVariable("nme", {adios2::LocalValueDim}); - internal->io.DefineVariable("offset", {adios2::LocalValueDim}); + auto nColumns = static_cast(size_one); + internal->io.DefineAttribute("columns", columnNames.data(), nColumns); + internal->io.DefineAttribute("columnstr", columns); + internal->io.DefineAttribute("boundarystr", boundstr); + internal->io.DefineAttribute("LAMMPS/dump_style", "atom"); + internal->io.DefineAttribute("LAMMPS/version", lmp->version); + internal->io.DefineAttribute("LAMMPS/num_ver", std::to_string(lmp->num_ver)); - // atom table size is not known at the moment - // it will be correctly defined at the moment of write - size_t UnknownSizeYet = 1; - internal->varAtoms = internal->io.DefineVariable( + // local dimension variables + internal->io.DefineVariable("nme", {adios2::LocalValueDim}); + internal->io.DefineVariable("offset", {adios2::LocalValueDim}); + + // atom table size is not known at the moment + // it will be correctly defined at the moment of write + size_t UnknownSizeYet = 1; + internal->varAtoms = internal->io.DefineVariable( "atoms", {UnknownSizeYet, nColumns}, {UnknownSizeYet, 0}, {UnknownSizeYet, nColumns}); + } } diff --git a/src/ADIOS/dump_custom_adios.cpp b/src/ADIOS/dump_custom_adios.cpp index 194141c15a..be827f507d 100644 --- a/src/ADIOS/dump_custom_adios.cpp +++ b/src/ADIOS/dump_custom_adios.cpp @@ -70,7 +70,11 @@ DumpCustomADIOS::DumpCustomADIOS(LAMMPS *lmp, int narg, char **arg) : DumpCustom internal = new DumpCustomADIOSInternal(); try { +#if defined(MPI_STUBS) + internal->ad = new adios2::ADIOS("adios2_config.xml", adios2::DebugON); +#else internal->ad = new adios2::ADIOS("adios2_config.xml", world, adios2::DebugON); +#endif } catch (std::ios_base::failure &e) { error->all(FLERR, "ADIOS initialization failed with error: {}", e.what()); } @@ -96,11 +100,19 @@ void DumpCustomADIOS::openfile() if (multifile) { // if one file per timestep, replace '*' with current timestep auto filecurrent = utils::star_subst(filename, update->ntimestep, padflag); +#if defined(MPI_STUBS) + internal->fh = internal->io.Open(filecurrent, adios2::Mode::Write); +#else internal->fh = internal->io.Open(filecurrent, adios2::Mode::Write, world); +#endif if (!internal->fh) error->one(FLERR, "Cannot open dump file {}", filecurrent); } else { if (!singlefile_opened) { +#if defined(MPI_STUBS) + internal->fh = internal->io.Open(filename, adios2::Mode::Write); +#else internal->fh = internal->io.Open(filename, adios2::Mode::Write, world); +#endif if (!internal->fh) error->one(FLERR, "Cannot open dump file {}", filename); singlefile_opened = 1; } diff --git a/src/ADIOS/reader_adios.cpp b/src/ADIOS/reader_adios.cpp index dd389be702..655df98d26 100644 --- a/src/ADIOS/reader_adios.cpp +++ b/src/ADIOS/reader_adios.cpp @@ -82,7 +82,11 @@ ReaderADIOS::ReaderADIOS(LAMMPS *lmp) : Reader(lmp) internal = new ReadADIOSInternal(); try { +#if defined(MPI_STUBS) + internal->ad = new adios2::ADIOS("adios2_config.xml", adios2::DebugON); +#else internal->ad = new adios2::ADIOS("adios2_config.xml", world, adios2::DebugON); +#endif } catch (std::ios_base::failure &e) { error->one(FLERR, "ADIOS initialization failed with error: {}", e.what()); } @@ -134,7 +138,11 @@ void ReaderADIOS::open_file(const std::string &file) if (internal->fh) internal->fh.Close(); try { +#if defined(MPI_STUBS) + internal->fh = internal->io.Open(file, adios2::Mode::Read); +#else internal->fh = internal->io.Open(file, adios2::Mode::Read, world); +#endif } catch (std::ios_base::failure &e) { error->one(FLERR, "Error opening file {}: {}", file, e.what()); } diff --git a/src/AMOEBA/amoeba_polar.cpp b/src/AMOEBA/amoeba_polar.cpp index a0f4fce301..9342b6213d 100644 --- a/src/AMOEBA/amoeba_polar.cpp +++ b/src/AMOEBA/amoeba_polar.cpp @@ -250,10 +250,12 @@ void PairAmoeba::polar_real() double drc3[3],drc5[3],drc7[3]; double urc3[3],urc5[3],tep[3]; double fix[3],fiy[3],fiz[3]; +#if 0 // for poltyp TCG which is currently not supported double uax[3],uay[3],uaz[3]; double ubx[3],uby[3],ubz[3]; double uaxp[3],uayp[3],uazp[3]; double ubxp[3],ubyp[3],ubzp[3]; +#endif double dmpi[9],dmpk[9]; double dmpik[9]; double bn[5]; @@ -320,6 +322,7 @@ void PairAmoeba::polar_real() uixp = uinp[i][0]; uiyp = uinp[i][1]; uizp = uinp[i][2]; +#if 0 // for poltyp TCG which is currently not supported for (m = 0; m < tcgnab; m++) { uax[m] = uad[m][i][0]; uay[m] = uad[m][i][1]; @@ -334,7 +337,7 @@ void PairAmoeba::polar_real() ubyp[m] = ubp[m][i][1]; ubzp[m] = ubp[m][i][2]; } - +#endif if (amoeba) { pdi = pdamp[itype]; pti = thole[itype]; @@ -1030,7 +1033,8 @@ void PairAmoeba::polar_real() // get the dtau/dr terms used for TCG polarization force } else if (poltyp == TCG) { - /* +#if 0 + // poltyp TCG not yet supported for AMOEBA/HIPPO for (m = 0; m < tcgnab; m++) { ukx = ubd[m][j][0]; uky = ubd[m][j][1]; @@ -1131,8 +1135,7 @@ void PairAmoeba::polar_real() frcx += depx; frcy += depy; frcz += depz; - } - */ +#endif } // increment force-based gradient on the interaction sites diff --git a/src/ASPHERE/pair_ylz.cpp b/src/ASPHERE/pair_ylz.cpp index 2d035b6f5d..afedf635d5 100644 --- a/src/ASPHERE/pair_ylz.cpp +++ b/src/ASPHERE/pair_ylz.cpp @@ -12,7 +12,7 @@ ------------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- - Contributing author: Mehdi Baghaee (SUSTech) + Contributing author: Hongyan Yuan (SUSTech) ------------------------------------------------------------------------- */ #include "pair_ylz.h" diff --git a/src/FEP/compute_fep.cpp b/src/FEP/compute_fep.cpp index bb11a64f54..13fcb98850 100644 --- a/src/FEP/compute_fep.cpp +++ b/src/FEP/compute_fep.cpp @@ -198,8 +198,12 @@ void ComputeFEP::init() if (pert->which == PAIR) { pairflag = 1; + Pair *pair = nullptr; - Pair *pair = force->pair_match(pert->pstyle, 1); + if (lmp->suffix_enable) + pair = force->pair_match(std::string(pert->pstyle)+"/"+lmp->suffix,1); + + if (pair == nullptr) pair = force->pair_match(pert->pstyle,1); if (pair == nullptr) error->all(FLERR, "compute fep pair style " diff --git a/src/PYTHON/python_impl.cpp b/src/PYTHON/python_impl.cpp index 92668674d1..9177381a35 100644 --- a/src/PYTHON/python_impl.cpp +++ b/src/PYTHON/python_impl.cpp @@ -112,37 +112,42 @@ PythonImpl::~PythonImpl() void PythonImpl::command(int narg, char **arg) { - if (narg < 2) error->all(FLERR, "Invalid python command"); + if (narg < 2) utils::missing_cmd_args(FLERR, "python", error); // if invoke is only keyword, invoke the previously defined function if (narg == 2 && strcmp(arg[1], "invoke") == 0) { int ifunc = find(arg[0]); - if (ifunc < 0) error->all(FLERR, "Python invoke of undefined function"); + if (ifunc < 0) error->all(FLERR, "Python invoke of unknown function: {}", arg[0]); char *str = nullptr; if (pfuncs[ifunc].noutput) { str = input->variable->pythonstyle(pfuncs[ifunc].ovarname, pfuncs[ifunc].name); - if (!str) error->all(FLERR, "Python variable does not match Python function"); + if (!str) + error->all(FLERR, + "Python variable {} does not match variable {} " + "registered with Python function {}", + arg[0], pfuncs[ifunc].ovarname, pfuncs[ifunc].name); } invoke_function(ifunc, str); return; } - // if source is only keyword, execute the python code + // if source is only keyword, execute the python code in file - if (narg == 3 && strcmp(arg[1], "source") == 0) { - int err; + if ((narg > 1) && (strcmp(arg[0], "source") == 0)) { + int err = -1; - FILE *fp = fopen(arg[2], "r"); - if (fp == nullptr) + if ((narg > 2) && (strcmp(arg[1], "here") == 0)) { err = execute_string(arg[2]); - else - err = execute_file(arg[2]); - - if (fp) fclose(fp); - if (err) error->all(FLERR, "Could not process Python source command"); + } else { + if (platform::file_is_readable(arg[1])) + err = execute_file(arg[1]); + else + error->all(FLERR, "Could not open python source file {} for processing", arg[1]); + } + if (err) error->all(FLERR, "Failure in python source command"); return; } @@ -162,48 +167,51 @@ void PythonImpl::command(int narg, char **arg) int iarg = 1; while (iarg < narg) { if (strcmp(arg[iarg], "input") == 0) { - if (iarg + 2 > narg) error->all(FLERR, "Invalid python command"); + if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "python input", error); ninput = utils::inumeric(FLERR, arg[iarg + 1], false, lmp); - if (ninput < 0) error->all(FLERR, "Invalid python command"); + if (ninput < 0) error->all(FLERR, "Invalid number of python input arguments: {}", ninput); iarg += 2; delete[] istr; istr = new char *[ninput]; - if (iarg + ninput > narg) error->all(FLERR, "Invalid python command"); + if (iarg + ninput > narg) utils::missing_cmd_args(FLERR, "python input", error); for (int i = 0; i < ninput; i++) istr[i] = arg[iarg + i]; iarg += ninput; } else if (strcmp(arg[iarg], "return") == 0) { - if (iarg + 2 > narg) error->all(FLERR, "Invalid python command"); + if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "python return", error); noutput = 1; ostr = arg[iarg + 1]; iarg += 2; } else if (strcmp(arg[iarg], "format") == 0) { - if (iarg + 2 > narg) error->all(FLERR, "Invalid python command"); + if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "python format", error); format = utils::strdup(arg[iarg + 1]); iarg += 2; } else if (strcmp(arg[iarg], "length") == 0) { - if (iarg + 2 > narg) error->all(FLERR, "Invalid python command"); + if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "python length", error); length_longstr = utils::inumeric(FLERR, arg[iarg + 1], false, lmp); - if (length_longstr <= 0) error->all(FLERR, "Invalid python command"); + if (length_longstr <= 0) error->all(FLERR, "Invalid python return value length"); iarg += 2; } else if (strcmp(arg[iarg], "file") == 0) { - if (iarg + 2 > narg) error->all(FLERR, "Invalid python command"); + if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "python file", error); delete[] pyfile; pyfile = utils::strdup(arg[iarg + 1]); iarg += 2; } else if (strcmp(arg[iarg], "here") == 0) { - if (iarg + 2 > narg) error->all(FLERR, "Invalid python command"); + if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "python here", error); herestr = arg[iarg + 1]; iarg += 2; } else if (strcmp(arg[iarg], "exists") == 0) { existflag = 1; iarg++; } else - error->all(FLERR, "Invalid python command"); + error->all(FLERR, "Unknown python command keyword: {}", arg[iarg]); } - if (pyfile && herestr) error->all(FLERR, "Invalid python command"); - if (pyfile && existflag) error->all(FLERR, "Invalid python command"); - if (herestr && existflag) error->all(FLERR, "Invalid python command"); + if (pyfile && herestr) + error->all(FLERR, "Must not use python 'file' and 'here' keywords at the same time"); + if (pyfile && existflag) + error->all(FLERR, "Must not use python 'file' and 'exists' keywords at the same time"); + if (herestr && existflag) + error->all(FLERR, "Must not use python 'here' and 'exists' keywords at the same time"); // create or overwrite entry in pfuncs vector with name = arg[0] @@ -221,23 +229,21 @@ void PythonImpl::command(int narg, char **arg) if (fp == nullptr) { PyUtils::Print_Errors(); - error->all(FLERR, "Could not open Python file"); + error->all(FLERR, "Could not open Python file: {}", pyfile); } int err = PyRun_SimpleFile(fp, pyfile); - if (err) { PyUtils::Print_Errors(); - error->all(FLERR, "Could not process Python file"); + error->all(FLERR, "Could not process Python file: {}", pyfile); } - fclose(fp); + } else if (herestr) { int err = PyRun_SimpleString(herestr); - if (err) { PyUtils::Print_Errors(); - error->all(FLERR, "Could not process Python string"); + error->all(FLERR, "Could not process Python string: {}", herestr); } } @@ -280,14 +286,17 @@ void PythonImpl::invoke_function(int ifunc, char *result) int ninput = pfuncs[ifunc].ninput; PyObject *pArgs = PyTuple_New(ninput); - if (!pArgs) { error->all(FLERR, "Could not create Python function arguments"); } + if (!pArgs) + error->all(FLERR, "Could not prepare arguments for Python function {}", pfuncs[ifunc].name); for (int i = 0; i < ninput; i++) { int itype = pfuncs[ifunc].itype[i]; if (itype == INT) { if (pfuncs[ifunc].ivarflag[i]) { str = input->variable->retrieve(pfuncs[ifunc].svalue[i]); - if (!str) { error->all(FLERR, "Could not evaluate Python function input variable"); } + if (!str) + error->all(FLERR, "Could not evaluate Python function {} input variable: {}", + pfuncs[ifunc].name, pfuncs[ifunc].svalue[i]); pValue = PY_INT_FROM_LONG(atoi(str)); } else { pValue = PY_INT_FROM_LONG(pfuncs[ifunc].ivalue[i]); @@ -295,7 +304,9 @@ void PythonImpl::invoke_function(int ifunc, char *result) } else if (itype == DOUBLE) { if (pfuncs[ifunc].ivarflag[i]) { str = input->variable->retrieve(pfuncs[ifunc].svalue[i]); - if (!str) { error->all(FLERR, "Could not evaluate Python function input variable"); } + if (!str) + error->all(FLERR, "Could not evaluate Python function {} input variable: {}", + pfuncs[ifunc].name, pfuncs[ifunc].svalue[i]); pValue = PyFloat_FromDouble(atof(str)); } else { pValue = PyFloat_FromDouble(pfuncs[ifunc].dvalue[i]); @@ -303,7 +314,9 @@ void PythonImpl::invoke_function(int ifunc, char *result) } else if (itype == STRING) { if (pfuncs[ifunc].ivarflag[i]) { str = input->variable->retrieve(pfuncs[ifunc].svalue[i]); - if (!str) { error->all(FLERR, "Could not evaluate Python function input variable"); } + if (!str) + error->all(FLERR, "Could not evaluate Python function {} input variable: {}", + pfuncs[ifunc].name, pfuncs[ifunc].svalue[i]); pValue = PY_STRING_FROM_STRING(str); } else { pValue = PY_STRING_FROM_STRING(pfuncs[ifunc].svalue[i]); @@ -311,7 +324,7 @@ void PythonImpl::invoke_function(int ifunc, char *result) } else if (itype == PTR) { pValue = PY_VOID_POINTER(lmp); } else { - error->all(FLERR, "Unsupported variable type"); + error->all(FLERR, "Unsupported variable type: {}", itype); } PyTuple_SetItem(pArgs, i, pValue); } @@ -324,7 +337,7 @@ void PythonImpl::invoke_function(int ifunc, char *result) if (!pValue) { PyUtils::Print_Errors(); - error->one(FLERR, "Python function evaluation failed"); + error->one(FLERR, "Python evaluation of function {} failed", pfuncs[ifunc].name); } // function returned a value @@ -365,9 +378,9 @@ int PythonImpl::variable_match(const char *name, const char *varname, int numeri { int ifunc = find(name); if (ifunc < 0) return -1; - if (pfuncs[ifunc].noutput == 0) return -1; - if (strcmp(pfuncs[ifunc].ovarname, varname) != 0) return -1; - if (numeric && pfuncs[ifunc].otype == STRING) return -1; + if (pfuncs[ifunc].noutput == 0) return -2; + if (strcmp(pfuncs[ifunc].ovarname, varname) != 0) return -3; + if (numeric && pfuncs[ifunc].otype == STRING) return -4; return ifunc; } @@ -400,9 +413,10 @@ int PythonImpl::create_entry(char *name, int ninput, int noutput, int length_lon pfuncs[ifunc].noutput = noutput; if (!format && ninput + noutput) - error->all(FLERR, "Invalid python command"); + error->all(FLERR, "Missing python format keyword"); else if (format && ((int) strlen(format) != ninput + noutput)) - error->all(FLERR, "Invalid python command"); + error->all(FLERR, "Input/output arguments ({}) and format characters ({}) are inconsistent", + (ninput + noutput), strlen(format)); // process inputs as values or variables @@ -448,7 +462,7 @@ int PythonImpl::create_entry(char *name, int ninput, int noutput, int length_lon if (strcmp(istr[i], "SELF") != 0) error->all(FLERR, "Invalid python command"); } else - error->all(FLERR, "Invalid python command"); + error->all(FLERR, "Invalid python format character: {}", type); } // process output as value or variable @@ -465,7 +479,7 @@ int PythonImpl::create_entry(char *name, int ninput, int noutput, int length_lon else if (type == 's') pfuncs[ifunc].otype = STRING; else - error->all(FLERR, "Invalid python command"); + error->all(FLERR, "Invalid python return format character: {}", type); if (length_longstr) { if (pfuncs[ifunc].otype != STRING) @@ -486,7 +500,9 @@ int PythonImpl::create_entry(char *name, int ninput, int noutput, int length_lon int PythonImpl::execute_string(char *cmd) { PyUtils::GIL lock; - return PyRun_SimpleString(cmd); + int err = PyRun_SimpleString(cmd); + if (err) PyUtils::Print_Errors(); + return err; } /* ---------------------------------------------------------------------- */ @@ -498,6 +514,7 @@ int PythonImpl::execute_file(char *fname) PyUtils::GIL lock; int err = PyRun_SimpleFile(fp, fname); + if (err) PyUtils::Print_Errors(); if (fp) fclose(fp); return err; diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index 5dc5967c29..8857590c4f 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -80,7 +80,7 @@ static const char cite_fix_bond_react[] = #define DELTA 16 #define MAXGUESS 20 // max # of guesses allowed by superimpose algorithm #define MAXCONARGS 14 // max # of arguments for any type of constraint + rxnID -#define NUMVARVALS 4 // max # of keyword values that have variables as input +#define NUMVARVALS 5 // max # of keyword values that have variables as input // various statuses of superimpose algorithm: // ACCEPT: site successfully matched to pre-reacted template @@ -98,7 +98,7 @@ enum{DISTANCE,ANGLE,DIHEDRAL,ARRHENIUS,RMSD,CUSTOM}; enum{ATOM,FRAG}; // keyword values that accept variables as input -enum{NEVERY,RMIN,RMAX,PROB}; +enum{NEVERY,RMIN,RMAX,PROB,NRATE}; // flag for one-proc vs shared reaction sites enum{LOCAL,GLOBAL}; @@ -138,10 +138,15 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : status = PROCEED; // reaction functions used by 'custom' constraint - nrxnfunction = 2; + nrxnfunction = 3; rxnfunclist.resize(nrxnfunction); + peratomflag.resize(nrxnfunction); rxnfunclist[0] = "rxnsum"; + peratomflag[0] = 1; rxnfunclist[1] = "rxnave"; + peratomflag[1] = 1; + rxnfunclist[2] = "rxnbond"; + peratomflag[2] = 0; nvvec = 0; ncustomvars = 0; vvec = nullptr; @@ -224,8 +229,10 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : memory->create(nghostlyskips,nreacts,"bond/react:nghostlyskips"); memory->create(seed,nreacts,"bond/react:seed"); memory->create(limit_duration,nreacts,"bond/react:limit_duration"); + memory->create(rate_limit,3,nreacts,"bond/react:rate_limit"); memory->create(stabilize_steps_flag,nreacts,"bond/react:stabilize_steps_flag"); memory->create(custom_charges_fragid,nreacts,"bond/react:custom_charges_fragid"); + memory->create(rescale_charges_flag,nreacts,"bond/react:rescale_charges_flag"); memory->create(create_atoms_flag,nreacts,"bond/react:create_atoms_flag"); memory->create(modify_create_fragid,nreacts,"bond/react:modify_create_fragid"); memory->create(overlapsq,nreacts,"bond/react:overlapsq"); @@ -249,8 +256,11 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : fraction[i] = 1; seed[i] = 12345; max_rxn[i] = INT_MAX; + for (int j = 0; j < 3; j++) + rate_limit[j][i] = 0; stabilize_steps_flag[i] = 0; custom_charges_fragid[i] = -1; + rescale_charges_flag[i] = 0; create_atoms_flag[i] = 0; modify_create_fragid[i] = -1; overlapsq[i] = 0; @@ -286,55 +296,31 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : if (groupid == -1) error->all(FLERR,"Could not find fix group ID"); groupbits[rxn] = group->bitmask[groupid]; - if (strncmp(arg[iarg],"v_",2) == 0) { - const char *str = &arg[iarg][2]; - var_id[NEVERY][rxn] = input->variable->find(str); - if (var_id[NEVERY][rxn] < 0) - error->all(FLERR,"Fix bond/react: Variable name {} does not exist", str); - if (!input->variable->equalstyle(var_id[NEVERY][rxn])) - error->all(FLERR,"Fix bond/react: Variable {} is not equal-style", str); - var_flag[NEVERY][rxn] = 1; - } else { + if (strncmp(arg[iarg],"v_",2) == 0) read_variable_keyword(&arg[iarg][2],NEVERY,rxn); + else { nevery[rxn] = utils::inumeric(FLERR,arg[iarg],false,lmp); if (nevery[rxn] <= 0) error->all(FLERR,"Illegal fix bond/react command: " "'Nevery' must be a positive integer"); } iarg++; + double cutoff; if (strncmp(arg[iarg],"v_",2) == 0) { - const char *str = &arg[iarg][2]; - var_id[RMIN][rxn] = input->variable->find(str); - if (var_id[RMIN][rxn] < 0) - error->all(FLERR,"Fix bond/react: Variable name {} does not exist", str); - if (!input->variable->equalstyle(var_id[RMIN][rxn])) - error->all(FLERR,"Fix bond/react: Variable {} is not equal-style", str); - double cutoff = input->variable->compute_equal(var_id[RMIN][rxn]); - cutsq[rxn][0] = cutoff*cutoff; - var_flag[RMIN][rxn] = 1; - } else { - double cutoff = utils::numeric(FLERR,arg[iarg],false,lmp); + read_variable_keyword(&arg[iarg][2],RMIN,rxn); + cutoff = input->variable->compute_equal(var_id[RMIN][rxn]); + } else cutoff = utils::numeric(FLERR,arg[iarg],false,lmp); if (cutoff < 0.0) error->all(FLERR,"Illegal fix bond/react command: " "'Rmin' cannot be negative"); cutsq[rxn][0] = cutoff*cutoff; - } iarg++; if (strncmp(arg[iarg],"v_",2) == 0) { - const char *str = &arg[iarg][2]; - var_id[RMAX][rxn] = input->variable->find(str); - if (var_id[RMAX][rxn] < 0) - error->all(FLERR,"Fix bond/react: Variable name {} does not exist", str); - if (!input->variable->equalstyle(var_id[RMAX][rxn])) - error->all(FLERR,"Fix bond/react: Variable is {} not equal-style", str); - double cutoff = input->variable->compute_equal(var_id[RMAX][rxn]); - cutsq[rxn][1] = cutoff*cutoff; - var_flag[RMAX][rxn] = 1; - } else { - double cutoff = utils::numeric(FLERR,arg[iarg],false,lmp); + read_variable_keyword(&arg[iarg][2],RMAX,rxn); + cutoff = input->variable->compute_equal(var_id[RMAX][rxn]); + } else cutoff = utils::numeric(FLERR,arg[iarg],false,lmp); if (cutoff < 0.0) error->all(FLERR,"Illegal fix bond/react command:" "'Rmax' cannot be negative"); cutsq[rxn][1] = cutoff*cutoff; - } iarg++; unreacted_mol[rxn] = atom->find_molecule(arg[iarg++]); @@ -344,7 +330,7 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : if (reacted_mol[rxn] == -1) error->all(FLERR,"Reacted molecule template ID for " "fix bond/react does not exist"); - // read superimpose file + //read map file files[rxn] = utils::strdup(arg[iarg]); iarg++; @@ -354,14 +340,8 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : "'prob' keyword has too few arguments"); // check if probability is a variable if (strncmp(arg[iarg+1],"v_",2) == 0) { - const char *str = &arg[iarg+1][2]; - var_id[PROB][rxn] = input->variable->find(str); - if (var_id[PROB][rxn] < 0) - error->all(FLERR,"Fix bond/react: Variable name {} does not exist", str); - if (!input->variable->equalstyle(var_id[PROB][rxn])) - error->all(FLERR,"Fix bond/react: Variable {} is not equal-style", str); + read_variable_keyword(&arg[iarg+1][2],PROB,rxn); fraction[rxn] = input->variable->compute_equal(var_id[PROB][rxn]); - var_flag[PROB][rxn] = 1; } else { // otherwise probability should be a number fraction[rxn] = utils::numeric(FLERR,arg[iarg+1],false,lmp); @@ -380,6 +360,14 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : if (max_rxn[rxn] < 0) error->all(FLERR,"Illegal fix bond/react command: " "'max_rxn' cannot be negative"); iarg += 2; + } else if (strcmp(arg[iarg],"rate_limit") == 0) { + if (iarg+3 > narg) error->all(FLERR,"Illegal fix bond/react command: " + "'rate_limit' has too few arguments"); + rate_limit[0][rxn] = 1; // serves as flag for rate_limit keyword + if (strncmp(arg[iarg+1],"v_",2) == 0) read_variable_keyword(&arg[iarg+1][2],NRATE,rxn); + else rate_limit[1][rxn] = utils::numeric(FLERR,arg[iarg+1],false,lmp); + rate_limit[2][rxn] = utils::numeric(FLERR,arg[iarg+2],false,lmp); + iarg += 3; } else if (strcmp(arg[iarg],"stabilize_steps") == 0) { if (stabilization_flag == 0) error->all(FLERR,"Stabilize_steps keyword " "used without stabilization keyword"); @@ -398,6 +386,13 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : "'custom_charges' keyword does not exist"); } iarg += 2; + } else if (strcmp(arg[iarg],"rescale_charges") == 0) { + if (iarg+2 > narg) error->all(FLERR,"Illegal fix bond/react command: " + "'rescale_charges' has too few arguments"); + if (strcmp(arg[iarg+1],"no") == 0) rescale_charges_flag[rxn] = 0; //default + else if (strcmp(arg[iarg+1],"yes") == 0) rescale_charges_flag[rxn] = 1; + else error->one(FLERR,"Bond/react: Illegal option for 'rescale_charges' keyword"); + iarg += 2; } else if (strcmp(arg[iarg],"molecule") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix bond/react command: " "'molecule' has too few arguments"); @@ -433,21 +428,24 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : } max_natoms = 0; // the number of atoms in largest molecule template + max_rate_limit_steps = 0; for (int myrxn = 0; myrxn < nreacts; myrxn++) { twomol = atom->molecules[reacted_mol[myrxn]]; max_natoms = MAX(max_natoms,twomol->natoms); + max_rate_limit_steps = MAX(max_rate_limit_steps,rate_limit[2][myrxn]); } memory->create(equivalences,max_natoms,2,nreacts,"bond/react:equivalences"); memory->create(reverse_equiv,max_natoms,2,nreacts,"bond/react:reverse_equiv"); memory->create(edge,max_natoms,nreacts,"bond/react:edge"); memory->create(landlocked_atoms,max_natoms,nreacts,"bond/react:landlocked_atoms"); + memory->create(store_rxn_count,max_rate_limit_steps,nreacts,"bond/react:store_rxn_count"); memory->create(custom_charges,max_natoms,nreacts,"bond/react:custom_charges"); memory->create(delete_atoms,max_natoms,nreacts,"bond/react:delete_atoms"); memory->create(create_atoms,max_natoms,nreacts,"bond/react:create_atoms"); memory->create(chiral_atoms,max_natoms,6,nreacts,"bond/react:chiral_atoms"); - for (int j = 0; j < nreacts; j++) + for (int j = 0; j < nreacts; j++) { for (int i = 0; i < max_natoms; i++) { edge[i][j] = 0; custom_charges[i][j] = 1; // update all partial charges by default @@ -462,6 +460,10 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : equivalences[i][m][j] = i+1; } } + for (int i = 0; i < max_rate_limit_steps; i++) { + store_rxn_count[i][j] = -1; + } + } // read all map files afterward for (int i = 0; i < nreacts; i++) { @@ -471,7 +473,7 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : onemol->check_attributes(); twomol->check_attributes(); get_molxspecials(); - read(i); + read_map_file(i); fclose(fp); if (ncreate == 0 && onemol->natoms != twomol->natoms) error->all(FLERR,"Fix bond/react: Reaction templates must contain the same number of atoms"); @@ -598,6 +600,7 @@ FixBondReact::~FixBondReact() memory->destroy(equivalences); memory->destroy(reverse_equiv); memory->destroy(landlocked_atoms); + memory->destroy(store_rxn_count); memory->destroy(custom_charges); memory->destroy(delete_atoms); memory->destroy(create_atoms); @@ -617,8 +620,10 @@ FixBondReact::~FixBondReact() memory->destroy(limit_duration); memory->destroy(var_flag); memory->destroy(var_id); + memory->destroy(rate_limit); memory->destroy(stabilize_steps_flag); memory->destroy(custom_charges_fragid); + memory->destroy(rescale_charges_flag); memory->destroy(molecule_keyword); memory->destroy(nconstraints); memory->destroy(constraintstr); @@ -819,6 +824,16 @@ void FixBondReact::init_list(int /*id*/, NeighList *ptr) void FixBondReact::post_integrate() { + // update store_rxn_count on every step + for (int myrxn = 0; myrxn < nreacts; myrxn++) { + if (rate_limit[0][myrxn] == 1) { + for (int i = rate_limit[2][myrxn]-1; i > 0; i--) { + store_rxn_count[i][myrxn] = store_rxn_count[i-1][myrxn]; + } + store_rxn_count[0][myrxn] = reaction_count_total[myrxn]; + } + } + // check if any reactions could occur on this timestep int nevery_check = 1; for (int i = 0; i < nreacts; i++) { @@ -873,6 +888,8 @@ void FixBondReact::post_integrate() for (int i = 0; i < nreacts; i++) { nattempt[i] = 0; } + // reset per-bond compute map flag + atoms2bondflag = 0; int nlocal = atom->nlocal; int nall = atom->nlocal + atom->nghost; @@ -911,8 +928,22 @@ void FixBondReact::post_integrate() int j; for (rxnID = 0; rxnID < nreacts; rxnID++) { + int rate_limit_flag = 1; + if (rate_limit[0][rxnID] == 1) { + int myrxn_count = store_rxn_count[rate_limit[2][rxnID]-1][rxnID]; + if (myrxn_count == -1) rate_limit_flag = 0; + else { + int nrxns_delta = reaction_count_total[rxnID] - myrxn_count; + int my_nrate; + if (var_flag[NRATE][rxnID] == 1) { + my_nrate = input->variable->compute_equal(var_id[NRATE][rxnID]); + } else my_nrate = rate_limit[1][rxnID]; + if (nrxns_delta > my_nrate) rate_limit_flag = 0; + } + } if ((update->ntimestep % nevery[rxnID]) || - (max_rxn[rxnID] <= reaction_count_total[rxnID])) continue; + (max_rxn[rxnID] <= reaction_count_total[rxnID]) || + (rate_limit_flag == 0)) continue; for (int ii = 0; ii < nall; ii++) { partner[ii] = 0; finalpartner[ii] = 0; @@ -1385,9 +1416,25 @@ void FixBondReact::superimpose_algorithm() std::random_device rnd; std::minstd_rand park_rng(rnd()); - // check if we overstepped our reaction limit + // check if we overstepped our reaction limit, via either max_rxn or rate_limit for (int i = 0; i < nreacts; i++) { - if (reaction_count_total[i] > max_rxn[i]) { + int overstep = 0; + int max_rxn_overstep = reaction_count_total[i] - max_rxn[i]; + overstep = MAX(overstep,max_rxn_overstep); + if (rate_limit[0][i] == 1) { + int myrxn_count = store_rxn_count[rate_limit[2][i]-1][i]; + if (myrxn_count != -1) { + int nrxn_delta = reaction_count_total[i] - myrxn_count; + int my_nrate; + if (var_flag[NRATE][i] == 1) { + my_nrate = input->variable->compute_equal(var_id[NRATE][i]); + } else my_nrate = rate_limit[1][i]; + int rate_limit_overstep = nrxn_delta - my_nrate; + overstep = MAX(overstep,rate_limit_overstep); + } + } + + if (overstep > 0) { // let's randomly choose rxns to skip, unbiasedly from local and ghostly int *local_rxncounts; int *all_localskips; @@ -1395,8 +1442,11 @@ void FixBondReact::superimpose_algorithm() memory->create(all_localskips,nprocs,"bond/react:all_localskips"); MPI_Gather(&local_rxn_count[i],1,MPI_INT,local_rxncounts,1,MPI_INT,0,world); if (me == 0) { - int overstep = reaction_count_total[i] - max_rxn[i]; int delta_rxn = reaction_count[i] + ghostly_rxn_count[i]; + // when using variable input for rate_limit, rate_limit_overstep could be > delta_rxn (below) + // we need to limit overstep to the number of reactions on this timestep + // essentially skipping all reactions, would be more efficient to use a skip_all flag + if (overstep > delta_rxn) overstep = delta_rxn; int *rxn_by_proc; memory->create(rxn_by_proc,delta_rxn,"bond/react:rxn_by_proc"); for (int j = 0; j < delta_rxn; j++) @@ -1414,14 +1464,15 @@ void FixBondReact::superimpose_algorithm() else all_localskips[rxn_by_proc[j]]++; } memory->destroy(rxn_by_proc); + reaction_count_total[i] -= overstep; } - reaction_count_total[i] = max_rxn[i]; MPI_Scatter(&all_localskips[0],1,MPI_INT,&nlocalskips[i],1,MPI_INT,0,world); MPI_Bcast(&nghostlyskips[i],1,MPI_INT,0,world); memory->destroy(local_rxncounts); memory->destroy(all_localskips); } } + MPI_Bcast(&reaction_count_total[0], nreacts, MPI_INT, 0, world); // this updates topology next step next_reneighbor = update->ntimestep; @@ -2118,7 +2169,7 @@ double FixBondReact::get_temperature(tagint **myglove, int row_offset, int col) } /* ---------------------------------------------------------------------- -get per-atom variable names used by custom constraint +get per-atom variable names used by custom constraint ------------------------------------------------------------------------- */ void FixBondReact::customvarnames() @@ -2140,6 +2191,7 @@ void FixBondReact::customvarnames() // find next reaction special function occurrence pos1 = std::string::npos; for (int i = 0; i < nrxnfunction; i++) { + if (peratomflag[i] == 0) continue; pos = varstr.find(rxnfunclist[i],prev3+1); if (pos == std::string::npos) continue; if (pos < pos1) pos1 = pos; @@ -2261,12 +2313,67 @@ double FixBondReact::custom_constraint(const std::string& varstr) } /* ---------------------------------------------------------------------- -currently two 'rxn' functions: rxnsum and rxnave +currently three 'rxn' functions: rxnsum, rxnave, and rxnbond ------------------------------------------------------------------------- */ double FixBondReact::rxnfunction(const std::string& rxnfunc, const std::string& varid, const std::string& fragid) { + int ifrag = -1; + if (fragid != "all") { + ifrag = onemol->findfragment(fragid.c_str()); + if (ifrag < 0) error->one(FLERR,"Bond/react: Molecule fragment " + "in reaction special function does not exist"); + } + + // start with 'rxnbond' per-bond function + // for 'rxnbond', varid corresponds to 'compute bond/local' name, + // and fragid is a pre-reaction fragment containing the two atoms in the bond + if (rxnfunc == "rxnbond") { + int icompute,ibond,nsum; + double perbondval; + std::set aset; + std::string computeid = varid; + std::map,int>::iterator it; + + if (computeid.substr(0,2) != "c_") error->one(FLERR,"Bond/react: Reaction special function compute " + "name should begin with 'c_'"); + computeid = computeid.substr(2); + icompute = modify->find_compute(computeid); + if (icompute < 0) error->one(FLERR,"Bond/react: Reaction special function compute name does not exist"); + cperbond = modify->compute[icompute]; + std::string compute_style = cperbond->style; + if (compute_style != "bond/local") error->one(FLERR,"Bond/react: Compute used by reaction " + "special function 'rxnbond' must be of style 'bond/local'"); + if (cperbond->size_local_cols > 0) error->one(FLERR,"Bond/react: 'Compute bond/local' used by reaction " + "special function 'rxnbond' must compute one value"); + + if (atoms2bondflag == 0) { + atoms2bondflag = 1; + get_atoms2bond(cperbond->groupbit); + } + + nsum = 0; + for (int i = 0; i < onemol->natoms; i++) { + if (onemol->fragmentmask[ifrag][i]) { + aset.insert(glove[i][1]); + nsum++; + } + } + if (nsum != 2) error->one(FLERR,"Bond/react: Molecule fragment of reaction special function 'rxnbond' " + "must contain exactly two atoms"); + + if (cperbond->invoked_local != lmp->update->ntimestep) + cperbond->compute_local(); + + it = atoms2bond.find(aset); + if (it == atoms2bond.end()) error->one(FLERR,"Bond/react: Unable to locate bond referenced by " + "reaction special function 'rxnbond'"); + ibond = it->second; + perbondval = cperbond->vector_local[ibond]; + return perbondval; + } + int ivar = -1; for (int i = 0; i < ncustomvars; i++) { if (varid == customvarstrs[i]) { @@ -2280,13 +2387,6 @@ double FixBondReact::rxnfunction(const std::string& rxnfunc, const std::string& error->one(FLERR,"Fix bond/react: Reaction special function variable " "name does not exist"); - int ifrag = -1; - if (fragid != "all") { - ifrag = onemol->findfragment(fragid.c_str()); - if (ifrag < 0) error->one(FLERR,"Fix bond/react: Molecule fragment " - "in reaction special function does not exist"); - } - int iatom; int nsum = 0; double sumvvec = 0; @@ -2313,6 +2413,39 @@ double FixBondReact::rxnfunction(const std::string& rxnfunc, const std::string& return 0.0; } +/* ---------------------------------------------------------------------- +populate map to get bond index from atom IDs +------------------------------------------------------------------------- */ + +void FixBondReact::get_atoms2bond(int cgroupbit) +{ + int i,m,atom1,atom2,btype,nb; + std::set aset; + + int nlocal = atom->nlocal; + tagint *tag = atom->tag; + int *num_bond = atom->num_bond; + tagint **bond_atom = atom->bond_atom; + int **bond_type = atom->bond_type; + int *mask = atom->mask; + + m = 0; + atoms2bond.clear(); + for (atom1 = 0; atom1 < nlocal; atom1++) { + if (!(mask[atom1] & cgroupbit)) continue; + nb = num_bond[atom1]; + for (i = 0; i < nb; i++) { + btype = bond_type[atom1][i]; + atom2 = atom->map(bond_atom[atom1][i]); + if (atom2 < 0 || !(mask[atom2] & cgroupbit)) continue; + if (newton_bond == 0 && tag[atom1] > tag[atom2]) continue; + if (btype == 0) continue; + aset = {tag[atom1], tag[atom2]}; + atoms2bond.insert(std::make_pair(aset,m++)); + } + } +} + /* ---------------------------------------------------------------------- return handedness (1 or -1) of a chiral center, given ordered set of coordinates ------------------------------------------------------------------------- */ @@ -2985,6 +3118,8 @@ void FixBondReact::update_everything() } delete [] iskip; + if (update_num_mega == 0) continue; + // if inserted atoms and global map exists, reset map now instead // of waiting for comm since other pre-exchange fixes may use it // invoke map_init() b/c atom count has grown @@ -3012,6 +3147,33 @@ void FixBondReact::update_everything() } } + // get charge rescale delta + double charge_rescale_addend = 0; + if (rescale_charges_flag[rxnID] == 1) { + double sim_total_charge = 0; + double mol_total_charge = 0; + int n_custom_charge = 0; + for (int i = 0; i < update_num_mega; i++) { + rxnID = update_mega_glove[0][i]; + twomol = atom->molecules[reacted_mol[rxnID]]; + for (int j = 0; j < twomol->natoms; j++) { + int jj = equivalences[j][1][rxnID]-1; + if (atom->map(update_mega_glove[jj+1][i]) >= 0 && + atom->map(update_mega_glove[jj+1][i]) < nlocal) { + if (landlocked_atoms[j][rxnID] == 1) + type[atom->map(update_mega_glove[jj+1][i])] = twomol->type[j]; + if (twomol->qflag && atom->q_flag && custom_charges[jj][rxnID] == 1) { + double *q = atom->q; + sim_total_charge += q[atom->map(update_mega_glove[jj+1][i])]; + mol_total_charge += twomol->q[j]; + n_custom_charge++; + } + } + } + } + charge_rescale_addend = (sim_total_charge-mol_total_charge)/n_custom_charge; + } + // update charges and types of landlocked atoms for (int i = 0; i < update_num_mega; i++) { rxnID = update_mega_glove[0][i]; @@ -3024,7 +3186,7 @@ void FixBondReact::update_everything() type[atom->map(update_mega_glove[jj+1][i])] = twomol->type[j]; if (twomol->qflag && atom->q_flag && custom_charges[jj][rxnID] == 1) { double *q = atom->q; - q[atom->map(update_mega_glove[jj+1][i])] = twomol->q[j]; + q[atom->map(update_mega_glove[jj+1][i])] = twomol->q[j]+charge_rescale_addend; } } } @@ -3795,10 +3957,24 @@ int FixBondReact::insert_atoms(tagint **my_mega_glove, int iupdate) } /* ---------------------------------------------------------------------- -read superimpose file +add equal-style variable to keyword argument list ------------------------------------------------------------------------- */ -void FixBondReact::read(int myrxn) +void FixBondReact::read_variable_keyword(const char *myarg, int keyword, int myrxn) +{ + var_id[keyword][myrxn] = input->variable->find(myarg); + if (var_id[keyword][myrxn] < 0) + error->all(FLERR,"Fix bond/react: Variable name {} does not exist",myarg); + if (!input->variable->equalstyle(var_id[keyword][myrxn])) + error->all(FLERR,"Fix bond/react: Variable {} is not equal-style",myarg); + var_flag[keyword][myrxn] = 1; +} + +/* ---------------------------------------------------------------------- +read map file +------------------------------------------------------------------------- */ + +void FixBondReact::read_map_file(int myrxn) { char line[MAXLINE],keyword[MAXLINE]; char *eof,*ptr; @@ -4299,34 +4475,71 @@ void FixBondReact::unpack_reverse_comm(int n, int *list, double *buf) void FixBondReact::write_restart(FILE *fp) { + int revision = 1; set[0].nreacts = nreacts; + set[0].max_rate_limit_steps = max_rate_limit_steps; + for (int i = 0; i < nreacts; i++) { set[i].reaction_count_total = reaction_count_total[i]; + strncpy(set[i].rxn_name,rxn_name[i],MAXLINE-1); set[i].rxn_name[MAXLINE-1] = '\0'; } - if (me == 0) { - int size = nreacts*sizeof(Set); - fwrite(&size,sizeof(int),1,fp); - fwrite(set,sizeof(Set),nreacts,fp); + int rbufcount = max_rate_limit_steps*nreacts; + int *rbuf; + if (rbufcount) { + memory->create(rbuf,rbufcount,"bond/react:rbuf"); + memcpy(rbuf,&store_rxn_count[0][0],sizeof(int)*rbufcount); } + + if (me == 0) { + int size = nreacts*sizeof(Set)+(rbufcount+1)*sizeof(int); + fwrite(&size,sizeof(int),1,fp); + fwrite(&revision,sizeof(int),1,fp); + fwrite(set,sizeof(Set),nreacts,fp); + if (rbufcount) fwrite(rbuf,sizeof(int),rbufcount,fp); + } + if (rbufcount) memory->destroy(rbuf); } /* ---------------------------------------------------------------------- use selected state info from restart file to restart the Fix + bond/react restart revisions numbers added after LAMMPS version 3 Nov 2022 ------------------------------------------------------------------------- */ void FixBondReact::restart(char *buf) { - Set *set_restart = (Set *) buf; - for (int i = 0; i < set_restart[0].nreacts; i++) { + int n,revision,r_nreacts,r_max_rate_limit_steps,ibufcount,n2cpy; + int **ibuf; + + n = 0; + if (lmp->restart_ver > utils::date2num("3 Nov 2022")) revision = buf[n++]; + else revision = 0; + + Set *set_restart = (Set *) &buf[n*sizeof(int)]; + r_nreacts = set_restart[0].nreacts; + + if (revision > 0) { + r_max_rate_limit_steps = set_restart[0].max_rate_limit_steps; + ibufcount = r_max_rate_limit_steps*r_nreacts; + memory->create(ibuf,r_max_rate_limit_steps,r_nreacts,"bond/react:ibuf"); + memcpy(&ibuf[0][0],&buf[sizeof(int)+r_nreacts*sizeof(Set)],sizeof(int)*ibufcount); + n2cpy = r_max_rate_limit_steps; + } else n2cpy = 0; + + if (max_rate_limit_steps < n2cpy) n2cpy = max_rate_limit_steps; + for (int i = 0; i < r_nreacts; i++) { for (int j = 0; j < nreacts; j++) { if (strcmp(set_restart[i].rxn_name,rxn_name[j]) == 0) { reaction_count_total[j] = set_restart[i].reaction_count_total; + // read rate_limit restart information + for (int k = 0; k < n2cpy; k++) + store_rxn_count[k][j] = ibuf[k][i]; } } } + if (revision > 0) memory->destroy(ibuf); } /* ---------------------------------------------------------------------- diff --git a/src/REACTION/fix_bond_react.h b/src/REACTION/fix_bond_react.h index e9ca7ec362..b55608b07e 100644 --- a/src/REACTION/fix_bond_react.h +++ b/src/REACTION/fix_bond_react.h @@ -25,6 +25,10 @@ FixStyle(bond/react,FixBondReact); #define LMP_FIX_BOND_REACT_H #include "fix.h" +#include "compute.h" + +#include +#include namespace LAMMPS_NS { @@ -64,8 +68,11 @@ class FixBondReact : public Fix { int stabilization_flag; int reset_mol_ids_flag; int custom_exclude_flag; + int **rate_limit; + int **store_rxn_count; int *stabilize_steps_flag; int *custom_charges_fragid; + int *rescale_charges_flag; int *create_atoms_flag; int *modify_create_fragid; double *overlapsq; @@ -74,9 +81,11 @@ class FixBondReact : public Fix { int *nconstraints; char **constraintstr; int nrxnfunction; - std::vector rxnfunclist; + std::vector rxnfunclist; // lists current special rxn function + std::vector peratomflag; // 1 if special rxn function uses per-atom variable (vs. per-bond) + int atoms2bondflag; // 1 if atoms2bond map has been populated on this timestep int narrhenius; - int **var_flag, **var_id; // for keyword values with variable inputs + int **var_flag, **var_id; // for keyword values with variable inputs int status; int *groupbits; @@ -86,6 +95,7 @@ class FixBondReact : public Fix { int *reaction_count_total; int nmax; // max num local atoms int max_natoms; // max natoms in a molecule template + int max_rate_limit_steps; // max rate limit interval tagint *partner, *finalpartner; double **distsq; int *nattempt; @@ -160,7 +170,8 @@ class FixBondReact : public Fix { // but whose first neighbors haven't int glove_counter; // used to determine when to terminate Superimpose Algorithm - void read(int); + void read_variable_keyword(const char *, int, int); + void read_map_file(int); void EdgeIDs(char *, int); void Equivalences(char *, int); void DeleteAtoms(char *, int); @@ -184,6 +195,7 @@ class FixBondReact : public Fix { double custom_constraint(const std::string &); // evaulate expression for custom constraint double rxnfunction(const std::string &, const std::string &, const std::string &); // eval rxn_sum and rxn_ave + void get_atoms2bond(int); int get_chirality(double[12]); // get handedness given an ordered set of coordinates void open(char *); @@ -198,16 +210,18 @@ class FixBondReact : public Fix { void ghost_glovecast(); void update_everything(); int insert_atoms(tagint **, int); - void unlimit_bond(); + void unlimit_bond(); // removes atoms from stabilization, and other post-reaction every-step operations void limit_bond(int); void dedup_mega_gloves(int); //dedup global mega_glove void write_restart(FILE *) override; void restart(char *buf) override; + // store restart data struct Set { int nreacts; char rxn_name[MAXLINE]; int reaction_count_total; + int max_rate_limit_steps; }; Set *set; @@ -221,7 +235,9 @@ class FixBondReact : public Fix { int ncustomvars; std::vector customvarstrs; int nvvec; - double **vvec; // per-atom vector to store variable constraint atom-style variable values + double **vvec; // per-atom vector to store custom constraint atom-style variable values + Compute *cperbond; // pointer to 'compute bond/local' used by custom constraint ('rxnbond' function) + std::map, int> atoms2bond; // maps atom pair to index of local bond array std::vector> constraints; // DEBUG diff --git a/src/lammps.cpp b/src/lammps.cpp index 9b9dda4ac8..5e9fed61b7 100644 --- a/src/lammps.cpp +++ b/src/lammps.cpp @@ -130,6 +130,7 @@ LAMMPS::LAMMPS(int narg, char **arg, MPI_Comm communicator) : version = (const char *) LAMMPS_VERSION; num_ver = utils::date2num(version); + restart_ver = -1; external_comm = 0; mdicomm = nullptr; @@ -993,6 +994,8 @@ void LAMMPS::destroy() delete python; python = nullptr; + + restart_ver = -1; // reset last restart version id } /* ---------------------------------------------------------------------- diff --git a/src/lammps.h b/src/lammps.h index abe08b45f8..03b019ed43 100644 --- a/src/lammps.h +++ b/src/lammps.h @@ -49,6 +49,8 @@ class LAMMPS { // that is constructed so that will be greater // for newer versions in numeric or string // value comparisons + int restart_ver; // -1 or numeric version id of LAMMPS version in restart + // file, in case LAMMPS was initialized from a restart // MPI_Comm world; // MPI communicator FILE *infile; // infile diff --git a/src/read_restart.cpp b/src/read_restart.cpp index 35d4629291..d4b455aaf9 100644 --- a/src/read_restart.cpp +++ b/src/read_restart.cpp @@ -599,6 +599,8 @@ void ReadRestart::header() if (flag == VERSION) { char *version = read_string(); + lmp->restart_ver = utils::date2num(version); + if (me == 0) utils::logmesg(lmp," restart file = {}, LAMMPS = {}\n", version, lmp->version); delete[] version; diff --git a/src/respa.cpp b/src/respa.cpp index 3e90a8baf9..71878a5af3 100644 --- a/src/respa.cpp +++ b/src/respa.cpp @@ -372,7 +372,7 @@ void Respa::setup(int flag) mesg += fmt::format(" {}:{}", ilevel + 1, step[ilevel]); mesg += "\n r-RESPA fixes :"; - for (int l = 0; l < modify->n_post_force_respa_any; ++l) { + for (int l = 0; l < modify->n_post_force_respa; ++l) { Fix *f = modify->get_fix_by_index(modify->list_post_force_respa[l]); if (f->respa_level >= 0) mesg += fmt::format(" {}:{}[{}]", MIN(f->respa_level + 1, nlevels), f->style, f->id); diff --git a/src/variable.cpp b/src/variable.cpp index 819f130a02..bb45649208 100644 --- a/src/variable.cpp +++ b/src/variable.cpp @@ -985,8 +985,21 @@ char *Variable::retrieve(const char *name) str = data[ivar][1] = utils::strdup(result); } else if (style[ivar] == PYTHON) { int ifunc = python->variable_match(data[ivar][0],name,0); - if (ifunc < 0) - error->all(FLERR,"Python variable {} does not match Python function {}", name, data[ivar][0]); + if (ifunc < 0) { + if (ifunc == -1) { + error->all(FLERR, "Could not find Python function {} linked to variable {}", + data[ivar][0], name); + } else if (ifunc == -2) { + error->all(FLERR, "Python function {} for variable {} does not have a return value", + data[ivar][0], name); + } else if (ifunc == -3) { + error->all(FLERR,"Python variable {} does not match variable name registered with " + "Python function {}", name, data[ivar][0]); + } else { + error->all(FLERR, "Unknown error verifying function {} linked to python style variable {}", + data[ivar][0],name); + } + } python->invoke_function(ifunc,data[ivar][1]); str = data[ivar][1]; // if Python func returns a string longer than VALUELENGTH diff --git a/tools/singularity/ubuntu18.04.def b/tools/singularity/ubuntu18.04.def index 4f8dd22c43..02351d9ecb 100644 --- a/tools/singularity/ubuntu18.04.def +++ b/tools/singularity/ubuntu18.04.def @@ -5,7 +5,7 @@ From: ubuntu:18.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/tools/singularity/ubuntu18.04_amd_rocm.def b/tools/singularity/ubuntu18.04_amd_rocm.def index 721591c589..38eaa6e322 100644 --- a/tools/singularity/ubuntu18.04_amd_rocm.def +++ b/tools/singularity/ubuntu18.04_amd_rocm.def @@ -45,7 +45,7 @@ From: ubuntu:18.04 ########################################################################### # Common Software ########################################################################### - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get install --no-install-recommends -y \ bc \ diff --git a/tools/singularity/ubuntu18.04_gpu.def b/tools/singularity/ubuntu18.04_gpu.def index d6947bb742..cab62c623f 100644 --- a/tools/singularity/ubuntu18.04_gpu.def +++ b/tools/singularity/ubuntu18.04_gpu.def @@ -48,7 +48,7 @@ From: ubuntu:18.04 ########################################################################### # Common Software ########################################################################### - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get install --no-install-recommends -y \ bc \ diff --git a/tools/singularity/ubuntu18.04_intel_opencl.def b/tools/singularity/ubuntu18.04_intel_opencl.def index a1df05fe0e..2d562771bb 100644 --- a/tools/singularity/ubuntu18.04_intel_opencl.def +++ b/tools/singularity/ubuntu18.04_intel_opencl.def @@ -5,7 +5,7 @@ From: ubuntu:18.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/tools/singularity/ubuntu18.04_nvidia.def b/tools/singularity/ubuntu18.04_nvidia.def index 7c0819a7ab..2a3a34b1b2 100644 --- a/tools/singularity/ubuntu18.04_nvidia.def +++ b/tools/singularity/ubuntu18.04_nvidia.def @@ -5,7 +5,7 @@ From: nvidia/cuda:11.6.2-devel-ubuntu18.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/tools/singularity/ubuntu20.04.def b/tools/singularity/ubuntu20.04.def index 8b3c92ada0..2a2a1dd660 100644 --- a/tools/singularity/ubuntu20.04.def +++ b/tools/singularity/ubuntu20.04.def @@ -5,7 +5,7 @@ From: ubuntu:20.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/tools/singularity/ubuntu20.04_amd_rocm.def b/tools/singularity/ubuntu20.04_amd_rocm.def index 3e018363b2..f947de9ee9 100644 --- a/tools/singularity/ubuntu20.04_amd_rocm.def +++ b/tools/singularity/ubuntu20.04_amd_rocm.def @@ -36,7 +36,7 @@ From: ubuntu:20.04 ########################################################################### # Common Software ########################################################################### - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get install --no-install-recommends -y \ bc \ diff --git a/tools/singularity/ubuntu20.04_gpu.def b/tools/singularity/ubuntu20.04_gpu.def index 9c9dea571d..f84c1f8926 100644 --- a/tools/singularity/ubuntu20.04_gpu.def +++ b/tools/singularity/ubuntu20.04_gpu.def @@ -39,7 +39,7 @@ From: ubuntu:20.04 # Common Software ########################################################################### apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get install --no-install-recommends -y \ bc \ diff --git a/tools/singularity/ubuntu20.04_intel_opencl.def b/tools/singularity/ubuntu20.04_intel_opencl.def index 96769516d4..cc547fef29 100644 --- a/tools/singularity/ubuntu20.04_intel_opencl.def +++ b/tools/singularity/ubuntu20.04_intel_opencl.def @@ -5,7 +5,7 @@ From: ubuntu:20.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/tools/singularity/ubuntu20.04_nvidia.def b/tools/singularity/ubuntu20.04_nvidia.def index 654b69f0b9..8ec334ad8b 100644 --- a/tools/singularity/ubuntu20.04_nvidia.def +++ b/tools/singularity/ubuntu20.04_nvidia.def @@ -5,7 +5,7 @@ From: nvidia/cuda:11.6.2-devel-ubuntu20.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/tools/singularity/ubuntu20.04_oneapi.def b/tools/singularity/ubuntu20.04_oneapi.def index 8b66d8d441..25c913f392 100644 --- a/tools/singularity/ubuntu20.04_oneapi.def +++ b/tools/singularity/ubuntu20.04_oneapi.def @@ -5,7 +5,7 @@ From: ubuntu:20.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/tools/singularity/ubuntu22.04.def b/tools/singularity/ubuntu22.04.def index a1a2a45b49..8a4f262f75 100644 --- a/tools/singularity/ubuntu22.04.def +++ b/tools/singularity/ubuntu22.04.def @@ -5,7 +5,7 @@ From: ubuntu:22.04 export DEBIAN_FRONTEND=noninteractive apt-get update apt-get install --no-install-recommends -y software-properties-common gpg gpg-agent - add-apt-repository ppa:openkim/latest + add-apt-repository ppa:openkim/latest -y apt-get update apt-get upgrade --no-install-recommends -y apt-get install --no-install-recommends -y \ diff --git a/unittest/formats/test_file_operations.cpp b/unittest/formats/test_file_operations.cpp index 6149b0c99c..9c75029130 100644 --- a/unittest/formats/test_file_operations.cpp +++ b/unittest/formats/test_file_operations.cpp @@ -303,6 +303,7 @@ TEST_F(FileOperationsTest, error_all_one) TEST_F(FileOperationsTest, write_restart) { + ASSERT_EQ(lmp->restart_ver, -1); BEGIN_HIDE_OUTPUT(); command("echo none"); END_HIDE_OUTPUT(); @@ -356,6 +357,7 @@ TEST_F(FileOperationsTest, write_restart) BEGIN_HIDE_OUTPUT(); command("clear"); END_HIDE_OUTPUT(); + ASSERT_EQ(lmp->restart_ver, -1); ASSERT_EQ(lmp->atom->natoms, 0); ASSERT_EQ(lmp->update->ntimestep, 0); ASSERT_EQ(lmp->domain->triclinic, 0); @@ -369,18 +371,21 @@ TEST_F(FileOperationsTest, write_restart) command("change_box all triclinic"); command("write_restart triclinic.restart"); END_HIDE_OUTPUT(); + ASSERT_EQ(lmp->restart_ver, lmp->num_ver); ASSERT_EQ(lmp->atom->natoms, 1); ASSERT_EQ(lmp->update->ntimestep, 333); ASSERT_EQ(lmp->domain->triclinic, 1); BEGIN_HIDE_OUTPUT(); command("clear"); END_HIDE_OUTPUT(); + ASSERT_EQ(lmp->restart_ver, -1); ASSERT_EQ(lmp->atom->natoms, 0); ASSERT_EQ(lmp->update->ntimestep, 0); ASSERT_EQ(lmp->domain->triclinic, 0); BEGIN_HIDE_OUTPUT(); command("read_restart triclinic.restart"); END_HIDE_OUTPUT(); + ASSERT_EQ(lmp->restart_ver, lmp->num_ver); ASSERT_EQ(lmp->atom->natoms, 1); ASSERT_EQ(lmp->update->ntimestep, 333); ASSERT_EQ(lmp->domain->triclinic, 1); diff --git a/unittest/python/test_python_package.cpp b/unittest/python/test_python_package.cpp index 4e5aa53b0c..f998ebbcd8 100644 --- a/unittest/python/test_python_package.cpp +++ b/unittest/python/test_python_package.cpp @@ -276,7 +276,7 @@ TEST_F(PythonPackageTest, RunSource) { // execute python script from file auto output = CAPTURE_OUTPUT([&] { - command("python xyz source ${input_dir}/run.py"); + command("python source ${input_dir}/run.py"); }); ASSERT_THAT(output, HasSubstr(LOREM_IPSUM)); @@ -286,7 +286,7 @@ TEST_F(PythonPackageTest, RunSourceInline) { // execute inline python script auto output = CAPTURE_OUTPUT([&] { - command("python xyz source \"\"\"\n" + command("python source here \"\"\"\n" "from __future__ import print_function\n" "print(2+2)\n" "\"\"\"");