ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge branch 'master' into feature/user-pace

Browse Source
This commit is contained in:
Axel Kohlmeyer
2021-04-07 16:26:26 -04:00
parent feda2dc08d 6569c2e434
commit 0496626533
469 changed files with 4650 additions and 2609 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
.github/workflows/unittest-macos.yml vendored
View File

@ -24,7 +24,9 @@ jobs:
shell: bash
working-directory: ${{github.workspace}}/build
run: |
cmake -C $GITHUB_WORKSPACE/cmake/presets/most.cmake $GITHUB_WORKSPACE/cmake \
cmake -C $GITHUB_WORKSPACE/cmake/presets/clang.cmake \
-C $GITHUB_WORKSPACE/cmake/presets/most.cmake \
$GITHUB_WORKSPACE/cmake \
-DENABLE_TESTING=ON -DBUILD_SHARED_LIBS=ON -DLAMMPS_EXCEPTIONS=ON
cmake --build . --parallel 2

5
cmake/CMakeLists.txt
View File

@ -22,6 +22,11 @@ set(LAMMPS_TOOLS_DIR ${LAMMPS_DIR}/tools)
set(LAMMPS_PYTHON_DIR ${LAMMPS_DIR}/python)
set(LAMMPS_POTENTIALS_DIR ${LAMMPS_DIR}/potentials)
set(LAMMPS_DOWNLOADS_URL "https://download.lammps.org" CACHE STRING "Base URL for LAMMPS downloads")
set(LAMMPS_POTENTIALS_URL "${LAMMPS_DOWNLOADS_URL}/potentials")
set(LAMMPS_THIRDPARTY_URL "${LAMMPS_DOWNLOADS_URL}/thirdparty")
mark_as_advanced(LAMMPS_DOWNLOADS_URL)
find_package(Git)
# by default, install into $HOME/.local (not /usr/local), so that no root access (and sudo!!) is needed

8
cmake/Modules/Documentation.cmake
View File

@ -55,11 +55,15 @@ if(BUILD_DOC)
COMMAND ${DOCENV_BINARY_DIR}/pip $ENV{PIP_OPTIONS} install -r ${DOC_BUILD_DIR}/requirements.txt --upgrade
)
set(MATHJAX_URL "https://github.com/mathjax/MathJax/archive/3.1.2.tar.gz" CACHE STRING "URL for MathJax tarball")
set(MATHJAX_MD5 "a4a6a093a89bc2ccab1452d766b98e53" CACHE STRING "MD5 checksum of MathJax tarball")
mark_as_advanced(MATHJAX_URL)
# download mathjax distribution and unpack to folder "mathjax"
if(NOT EXISTS ${DOC_BUILD_STATIC_DIR}/mathjax/es5)
file(DOWNLOAD "https://github.com/mathjax/MathJax/archive/3.1.2.tar.gz"
file(DOWNLOAD ${MATHJAX_URL}
"${CMAKE_CURRENT_BINARY_DIR}/mathjax.tar.gz"
EXPECTED_MD5 a4a6a093a89bc2ccab1452d766b98e53)
EXPECTED_MD5 ${MATHJAX_MD5})
execute_process(COMMAND ${CMAKE_COMMAND} -E tar xzf mathjax.tar.gz WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})
file(GLOB MATHJAX_VERSION_DIR ${CMAKE_CURRENT_BINARY_DIR}/MathJax-*)
execute_process(COMMAND ${CMAKE_COMMAND} -E rename ${MATHJAX_VERSION_DIR} ${DOC_BUILD_STATIC_DIR}/mathjax)

6
cmake/Modules/GTest.cmake
View File

@ -8,10 +8,12 @@ endif()
include(ExternalProject)
set(GTEST_URL "https://github.com/google/googletest/archive/release-1.10.0.tar.gz" CACHE STRING "URL for GTest tarball")
set(GTEST_MD5 "ecd1fa65e7de707cd5c00bdac56022cd" CACHE STRING "MD5 checksum of GTest tarball")
mark_as_advanced(GTEST_URL)
mark_as_advanced(GTEST_MD5)
ExternalProject_Add(googletest
URL ${GTEST_URL}
URL_MD5 ecd1fa65e7de707cd5c00bdac56022cd
URL ${GTEST_URL}
URL_MD5 ${GTEST_MD5}
SOURCE_DIR "${CMAKE_BINARY_DIR}/gtest-src"
BINARY_DIR "${CMAKE_BINARY_DIR}/gtest-build"
CMAKE_ARGS ${CMAKE_REQUEST_PIC} ${CMAKE_EXTRA_GTEST_OPTS}

1
cmake/Modules/LAMMPSUtils.cmake
View File

@ -86,7 +86,6 @@ endfunction(GenerateBinaryHeader)
# fetch missing potential files
function(FetchPotentials pkgfolder potfolder)
if (EXISTS "${pkgfolder}/potentials.txt")
set(LAMMPS_POTENTIALS_URL "https://download.lammps.org/potentials")
file(STRINGS "${pkgfolder}/potentials.txt" linelist REGEX "^[^#].")
foreach(line ${linelist})
string(FIND ${line} " " blank)

17
cmake/Modules/MPI4WIN.cmake
View File

@ -1,16 +1,25 @@
# 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 https://download.lammps.org/thirdparty/mpich2-win64-devel.tar.gz
URL_MD5 4939fdb59d13182fd5dd65211e469f14
URL ${MPICH2_WIN64_DEVEL_URL}
URL_MD5 ${MPICH2_WIN64_DEVEL_MD5}
CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND ""
BUILD_BYPRODUCTS <SOURCE_DIR>/lib/libmpi.a)
else()
ExternalProject_Add(mpi4win_build
URL https://download.lammps.org/thirdparty/mpich2-win32-devel.tar.gz
URL_MD5 a61d153500dce44e21b755ee7257e031
URL ${MPICH2_WIN32_DEVEL_URL}
URL_MD5 ${MPICH2_WIN32_DEVEL_MD5}
CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND ""
BUILD_BYPRODUCTS <SOURCE_DIR>/lib/libmpi.a)
endif()

10
cmake/Modules/OpenCLLoader.cmake
View File

@ -1,11 +1,13 @@
message(STATUS "Downloading and building OpenCL loader library")
set(OPENCL_LOADER_URL "${LAMMPS_THIRDPARTY_URL}/opencl-loader-2020.12.18.tar.gz" CACHE STRING "URL for OpenCL loader tarball")
set(OPENCL_LOADER_MD5 "011cdcbd41030be94f3fced6d763a52a" CACHE STRING "MD5 checksum of OpenCL loader tarball")
mark_as_advanced(OPENCL_LOADER_URL)
mark_as_advanced(OPENCL_LOADER_MD5)
include(ExternalProject)
set(OPENCL_LOADER_URL "https://download.lammps.org/thirdparty/opencl-loader-2020.12.18.tar.gz" CACHE STRING "URL for OpenCL loader tarball")
mark_as_advanced(OPENCL_LOADER_URL)
ExternalProject_Add(opencl_loader
URL ${OPENCL_LOADER_URL}
URL_MD5 011cdcbd41030be94f3fced6d763a52a
URL ${OPENCL_LOADER_URL}
URL_MD5 ${OPENCL_LOADER_MD5}
SOURCE_DIR "${CMAKE_BINARY_DIR}/opencl_loader-src"
BINARY_DIR "${CMAKE_BINARY_DIR}/opencl_loader-build"
CMAKE_ARGS ${CMAKE_REQUEST_PIC} ${CMAKE_EXTRA_OPENCL_LOADER_OPTS}

25
cmake/Modules/Packages/GPU.cmake
View File

@ -131,7 +131,7 @@ if(GPU_API STREQUAL "CUDA")
add_library(gpu STATIC ${GPU_LIB_SOURCES} ${GPU_LIB_CUDPP_SOURCES} ${GPU_OBJS})
target_link_libraries(gpu PRIVATE ${CUDA_LIBRARIES} ${CUDA_CUDA_LIBRARY})
target_include_directories(gpu PRIVATE ${LAMMPS_LIB_BINARY_DIR}/gpu ${CUDA_INCLUDE_DIRS})
target_compile_definitions(gpu PRIVATE -D_${GPU_PREC_SETTING} -DMPI_GERYON -DUCL_NO_EXIT ${GPU_CUDA_MPS_FLAGS})
target_compile_definitions(gpu PRIVATE -DUSE_CUDA -D_${GPU_PREC_SETTING} -DMPI_GERYON -DUCL_NO_EXIT ${GPU_CUDA_MPS_FLAGS})
if(CUDPP_OPT)
target_include_directories(gpu PRIVATE ${LAMMPS_LIB_SOURCE_DIR}/gpu/cudpp_mini)
target_compile_definitions(gpu PRIVATE -DUSE_CUDPP)
@ -218,7 +218,7 @@ elseif(GPU_API STREQUAL "HIP")
if(NOT DEFINED HIP_PLATFORM)
if(NOT DEFINED ENV{HIP_PLATFORM})
set(HIP_PLATFORM "hcc" CACHE PATH "HIP Platform to be used during compilation")
set(HIP_PLATFORM "amd" CACHE PATH "HIP Platform to be used during compilation")
else()
set(HIP_PLATFORM $ENV{HIP_PLATFORM} CACHE PATH "HIP Platform used during compilation")
endif()
@ -226,7 +226,7 @@ elseif(GPU_API STREQUAL "HIP")
set(ENV{HIP_PLATFORM} ${HIP_PLATFORM})
if(HIP_PLATFORM STREQUAL "hcc")
if(HIP_PLATFORM STREQUAL "hcc" OR HIP_PLATFORM STREQUAL "amd")
set(HIP_ARCH "gfx906" CACHE STRING "HIP target architecture")
elseif(HIP_PLATFORM STREQUAL "nvcc")
find_package(CUDA REQUIRED)
@ -284,7 +284,7 @@ elseif(GPU_API STREQUAL "HIP")
set(CUBIN_FILE "${LAMMPS_LIB_BINARY_DIR}/gpu/${CU_NAME}.cubin")
set(CUBIN_H_FILE "${LAMMPS_LIB_BINARY_DIR}/gpu/${CU_NAME}_cubin.h")
if(HIP_PLATFORM STREQUAL "hcc")
if(HIP_PLATFORM STREQUAL "hcc" OR HIP_PLATFORM STREQUAL "amd")
configure_file(${CU_FILE} ${CU_CPP_FILE} COPYONLY)
if(HIP_COMPILER STREQUAL "clang")
@ -338,11 +338,16 @@ elseif(GPU_API STREQUAL "HIP")
if(DOWNLOAD_CUB)
message(STATUS "CUB download requested")
set(CUB_URL "https://github.com/NVlabs/cub/archive/1.12.0.tar.gz" CACHE STRING "URL for CUB tarball")
set(CUB_MD5 "1cf595beacafff104700921bac8519f3" CACHE STRING "MD5 checksum of CUB tarball")
mark_as_advanced(CUB_URL)
mark_as_advanced(CUB_MD5)
include(ExternalProject)
ExternalProject_Add(CUB
GIT_REPOSITORY https://github.com/NVlabs/cub
TIMEOUT 5
URL ${CUB_URL}
URL_MD5 ${CUB_MD5}
PREFIX "${CMAKE_CURRENT_BINARY_DIR}"
CONFIGURE_COMMAND ""
BUILD_COMMAND ""
@ -354,7 +359,7 @@ elseif(GPU_API STREQUAL "HIP")
else()
find_package(CUB)
if(NOT CUB_FOUND)
message(FATAL_ERROR "CUB library not found. Help CMake to find it by setting CUB_INCLUDE_DIR, or set DOWNLOAD_VORO=ON to download it")
message(FATAL_ERROR "CUB library not found. Help CMake to find it by setting CUB_INCLUDE_DIR, or set DOWNLOAD_CUB=ON to download it")
endif()
endif()
@ -381,6 +386,12 @@ elseif(GPU_API STREQUAL "HIP")
target_compile_definitions(hip_get_devices PRIVATE -D__HIP_PLATFORM_HCC__)
target_include_directories(hip_get_devices PRIVATE ${HIP_ROOT_DIR}/../include)
elseif(HIP_PLATFORM STREQUAL "amd")
target_compile_definitions(gpu PRIVATE -D__HIP_PLATFORM_AMD__)
target_include_directories(gpu PRIVATE ${HIP_ROOT_DIR}/../include)
target_compile_definitions(hip_get_devices PRIVATE -D__HIP_PLATFORM_AMD__)
target_include_directories(hip_get_devices PRIVATE ${HIP_ROOT_DIR}/../include)
endif()
target_link_libraries(lammps PRIVATE gpu)

8
cmake/Modules/Packages/KIM.cmake
View File

@ -35,9 +35,13 @@ if(DOWNLOAD_KIM)
include(ExternalProject)
enable_language(C)
enable_language(Fortran)
set(KIM_URL "https://s3.openkim.org/kim-api/kim-api-2.2.1.txz" CACHE STRING "URL for KIM tarball")
set(KIM_MD5 "ae1ddda2ef7017ea07934e519d023dca" CACHE STRING "MD5 checksum of KIM tarball")
mark_as_advanced(KIM_URL)
mark_as_advanced(KIM_MD5)
ExternalProject_Add(kim_build
URL https://s3.openkim.org/kim-api/kim-api-2.2.1.txz
URL_MD5 ae1ddda2ef7017ea07934e519d023dca
URL ${KIM_URL}
URL_MD5 ${KIM_MD5}
BINARY_DIR build
CMAKE_ARGS ${CMAKE_REQUEST_PIC}
-DCMAKE_C_COMPILER=${CMAKE_C_COMPILER}

8
cmake/Modules/Packages/KOKKOS.cmake
View File

@ -37,9 +37,13 @@ if(DOWNLOAD_KOKKOS)
list(APPEND KOKKOS_LIB_BUILD_ARGS "-DCMAKE_CXX_EXTENSIONS=${CMAKE_CXX_EXTENSIONS}")
list(APPEND KOKKOS_LIB_BUILD_ARGS "-DCMAKE_TOOLCHAIN_FILE=${CMAKE_TOOLCHAIN_FILE}")
include(ExternalProject)
set(KOKKOS_URL "https://github.com/kokkos/kokkos/archive/3.3.01.tar.gz" CACHE STRING "URL for KOKKOS tarball")
set(KOKKOS_MD5 "08201d1c7cf5bc458ce0f5b44a629d5a" CACHE STRING "MD5 checksum of KOKKOS tarball")
mark_as_advanced(KOKKOS_URL)
mark_as_advanced(KOKKOS_MD5)
ExternalProject_Add(kokkos_build
URL https://github.com/kokkos/kokkos/archive/3.3.01.tar.gz
URL_MD5 08201d1c7cf5bc458ce0f5b44a629d5a
URL ${KOKKOS_URL}
URL_MD5 ${KOKKOS_MD5}
CMAKE_ARGS ${KOKKOS_LIB_BUILD_ARGS}
BUILD_BYPRODUCTS <INSTALL_DIR>/lib/libkokkoscore.a
)

8
cmake/Modules/Packages/LATTE.cmake
View File

@ -15,10 +15,14 @@ endif()
option(DOWNLOAD_LATTE "Download the LATTE library instead of using an already installed one" ${DOWNLOAD_LATTE_DEFAULT})
if(DOWNLOAD_LATTE)
message(STATUS "LATTE download requested - we will build our own")
set(LATTE_URL "https://github.com/lanl/LATTE/archive/v1.2.2.tar.gz" CACHE STRING "URL for LATTE tarball")
set(LATTE_MD5 "820e73a457ced178c08c71389a385de7" CACHE STRING "MD5 checksum of LATTE tarball")
mark_as_advanced(LATTE_URL)
mark_as_advanced(LATTE_MD5)
include(ExternalProject)
ExternalProject_Add(latte_build
URL https://github.com/lanl/LATTE/archive/v1.2.2.tar.gz
URL_MD5 820e73a457ced178c08c71389a385de7
URL ${LATTE_URL}
URL_MD5 ${LATTE_MD5}
SOURCE_SUBDIR cmake
CMAKE_ARGS -DCMAKE_INSTALL_PREFIX=<INSTALL_DIR> ${CMAKE_REQUEST_PIC} -DCMAKE_INSTALL_LIBDIR=lib
-DBLAS_LIBRARIES=${BLAS_LIBRARIES} -DLAPACK_LIBRARIES=${LAPACK_LIBRARIES}

9
cmake/Modules/Packages/MSCG.cmake
View File

@ -7,10 +7,15 @@ else()
endif()
option(DOWNLOAD_MSCG "Download MSCG library instead of using an already installed one)" ${DOWNLOAD_MSCG_DEFAULT})
if(DOWNLOAD_MSCG)
set(MSCG_URL "https://github.com/uchicago-voth/MSCG-release/archive/1.7.3.1.tar.gz" CACHE STRING "URL for MSCG tarball")
set(MSCG_MD5 "8c45e269ee13f60b303edd7823866a91" CACHE STRING "MD5 checksum of MSCG tarball")
mark_as_advanced(MSCG_URL)
mark_as_advanced(MSCG_MD5)
include(ExternalProject)
ExternalProject_Add(mscg_build
URL https://github.com/uchicago-voth/MSCG-release/archive/1.7.3.1.tar.gz
URL_MD5 8c45e269ee13f60b303edd7823866a91
URL ${MSCG_URL}
URL_MD5 ${MSCG_MD5}
SOURCE_SUBDIR src/CMake
CMAKE_ARGS ${CMAKE_REQUEST_PIC} ${EXTRA_MSCG_OPTS}
-DCMAKE_C_COMPILER=${CMAKE_C_COMPILER}

10
cmake/Modules/Packages/USER-PLUMED.cmake
View File

@ -53,10 +53,16 @@ if(DOWNLOAD_PLUMED)
elseif(PLUMED_MODE STREQUAL "RUNTIME")
set(PLUMED_BUILD_BYPRODUCTS "<INSTALL_DIR>/lib/libplumedWrapper.a")
endif()
set(PLUMED_URL "https://github.com/plumed/plumed2/releases/download/v2.7.0/plumed-src-2.7.0.tgz" CACHE STRING "URL for PLUMED tarball")
set(PLUMED_MD5 "95f29dd0c067577f11972ff90dfc7d12" CACHE STRING "MD5 checksum of PLUMED tarball")
mark_as_advanced(PLUMED_URL)
mark_as_advanced(PLUMED_MD5)
include(ExternalProject)
ExternalProject_Add(plumed_build
URL https://github.com/plumed/plumed2/releases/download/v2.7.0/plumed-src-2.7.0.tgz
URL_MD5 95f29dd0c067577f11972ff90dfc7d12
URL ${PLUMED_URL}
URL_MD5 ${PLUMED_MD5}
BUILD_IN_SOURCE 1
CONFIGURE_COMMAND <SOURCE_DIR>/configure --prefix=<INSTALL_DIR>
${CONFIGURE_REQUEST_PIC}

10
cmake/Modules/Packages/USER-SCAFACOS.cmake
View File

@ -14,15 +14,19 @@ endif()
option(DOWNLOAD_SCAFACOS "Download ScaFaCoS library instead of using an already installed one" ${DOWNLOAD_SCAFACOS_DEFAULT})
if(DOWNLOAD_SCAFACOS)
message(STATUS "ScaFaCoS download requested - we will build our own")
set(SCAFACOS_URL "https://github.com/scafacos/scafacos/releases/download/v1.0.1/scafacos-1.0.1.tar.gz" CACHE STRING "URL for SCAFACOS tarball")
set(SCAFACOS_MD5 "bd46d74e3296bd8a444d731bb10c1738" CACHE STRING "MD5 checksum of SCAFACOS tarball")
mark_as_advanced(SCAFACOS_URL)
mark_as_advanced(SCAFACOS_MD5)
# version 1.0.1 needs a patch to compile and linke cleanly with GCC 10 and later.
file(DOWNLOAD https://download.lammps.org/thirdparty/scafacos-1.0.1-fix.diff ${CMAKE_CURRENT_BINARY_DIR}/scafacos-1.0.1.fix.diff
file(DOWNLOAD ${LAMMPS_THIRDPARTY_URL}/scafacos-1.0.1-fix.diff ${CMAKE_CURRENT_BINARY_DIR}/scafacos-1.0.1.fix.diff
EXPECTED_HASH MD5=4baa1333bb28fcce102d505e1992d032)
include(ExternalProject)
ExternalProject_Add(scafacos_build
URL https://github.com/scafacos/scafacos/releases/download/v1.0.1/scafacos-1.0.1.tar.gz
URL_MD5 bd46d74e3296bd8a444d731bb10c1738
URL ${SCAFACOS_URL}
URL_MD5 ${SCAFACOS_MD5}
PATCH_COMMAND patch -p1 < ${CMAKE_CURRENT_BINARY_DIR}/scafacos-1.0.1.fix.diff
CONFIGURE_COMMAND <SOURCE_DIR>/configure --prefix=<INSTALL_DIR> --disable-doc
--enable-fcs-solvers=fmm,p2nfft,direct,ewald,p3m

8
cmake/Modules/Packages/USER-SMD.cmake
View File

@ -7,10 +7,14 @@ endif()
option(DOWNLOAD_EIGEN3 "Download Eigen3 instead of using an already installed one)" ${DOWNLOAD_EIGEN3_DEFAULT})
if(DOWNLOAD_EIGEN3)
message(STATUS "Eigen3 download requested - we will build our own")
set(EIGEN3_URL "https://gitlab.com/libeigen/eigen/-/archive/3.3.7/eigen-3.3.7.tar.gz" CACHE STRING "URL for Eigen3 tarball")
set(EIGEN3_MD5 "9e30f67e8531477de4117506fe44669b" CACHE STRING "MD5 checksum of Eigen3 tarball")
mark_as_advanced(EIGEN3_URL)
mark_as_advanced(EIGEN3_MD5)
include(ExternalProject)
ExternalProject_Add(Eigen3_build
URL https://gitlab.com/libeigen/eigen/-/archive/3.3.7/eigen-3.3.7.tar.gz
URL_MD5 9e30f67e8531477de4117506fe44669b
URL ${EIGEN3_URL}
URL_MD5 ${EIGEN3_MD5}
CONFIGURE_COMMAND "" BUILD_COMMAND "" INSTALL_COMMAND ""
)
ExternalProject_get_property(Eigen3_build SOURCE_DIR)

9
cmake/Modules/Packages/VORONOI.cmake
View File

@ -7,6 +7,11 @@ endif()
option(DOWNLOAD_VORO "Download and compile the Voro++ library instead of using an already installed one" ${DOWNLOAD_VORO_DEFAULT})
if(DOWNLOAD_VORO)
message(STATUS "Voro++ download requested - we will build our own")
set(VORO_URL "${LAMMPS_THIRDPARTY_URL}/voro++-0.4.6.tar.gz" CACHE STRING "URL for Voro++ tarball")
set(VORO_MD5 "2338b824c3b7b25590e18e8df5d68af9" CACHE STRING "MD5 checksum for Voro++ tarball")
mark_as_advanced(VORO_URL)
mark_as_advanced(VORO_MD5)
include(ExternalProject)
if(BUILD_SHARED_LIBS)
@ -22,8 +27,8 @@ if(DOWNLOAD_VORO)
endif()
ExternalProject_Add(voro_build
URL https://download.lammps.org/thirdparty/voro++-0.4.6.tar.gz
URL_MD5 2338b824c3b7b25590e18e8df5d68af9
URL ${VORO_URL}
URL_MD5 ${VORO_MD5}
PATCH_COMMAND patch -b -p0 < ${LAMMPS_LIB_SOURCE_DIR}/voronoi/voro-make.patch
CONFIGURE_COMMAND ""
BUILD_COMMAND make ${VORO_BUILD_OPTIONS}

5
cmake/Modules/YAML.cmake
View File

@ -2,10 +2,13 @@ message(STATUS "Downloading and building YAML library")
include(ExternalProject)
set(YAML_URL "https://pyyaml.org/download/libyaml/yaml-0.2.5.tar.gz" CACHE STRING "URL for libyaml tarball")
set(YAML_MD5 "bb15429d8fb787e7d3f1c83ae129a999" CACHE STRING "MD5 checksum of libyaml tarball")
mark_as_advanced(YAML_URL)
mark_as_advanced(YAML_MD5)
ExternalProject_Add(libyaml
URL ${YAML_URL}
URL_MD5 bb15429d8fb787e7d3f1c83ae129a999
URL_MD5 ${YAML_MD5}
SOURCE_DIR "${CMAKE_BINARY_DIR}/yaml-src"
BINARY_DIR "${CMAKE_BINARY_DIR}/yaml-build"
CONFIGURE_COMMAND <SOURCE_DIR>/configure ${CONFIGURE_REQUEST_PIC}

13
doc/src/Build_extras.rst
View File

@ -125,7 +125,7 @@ CMake build
# default is sm_50
-D HIP_ARCH=value # primary GPU hardware choice for GPU_API=hip
# value depends on selected HIP_PLATFORM
# default is 'gfx906' for HIP_PLATFORM=hcc and 'sm_50' for HIP_PLATFORM=nvcc
# default is 'gfx906' for HIP_PLATFORM=amd and 'sm_50' for HIP_PLATFORM=nvcc
-D HIP_USE_DEVICE_SORT=value # enables GPU sorting
# value = yes (default) or no
-D CUDPP_OPT=value # use GPU binning on with CUDA (should be off for modern GPUs)
@ -169,17 +169,24 @@ desired, you can set :code:`USE_STATIC_OPENCL_LOADER` to :code:`no`.
If you are compiling with HIP, note that before running CMake you will have to
set appropriate environment variables. Some variables such as
:code:`HCC_AMDGPU_TARGET` or :code:`CUDA_PATH` are necessary for :code:`hipcc`
:code:`HCC_AMDGPU_TARGET` (for ROCm <= 4.0) or :code:`CUDA_PATH` are necessary for :code:`hipcc`
and the linker to work correctly.
.. code:: bash
# AMDGPU target
# AMDGPU target (ROCm <= 4.0)
export HIP_PLATFORM=hcc
export HCC_AMDGPU_TARGET=gfx906
cmake -D PKG_GPU=on -D GPU_API=HIP -D HIP_ARCH=gfx906 -D CMAKE_CXX_COMPILER=hipcc ..
make -j 4
.. code:: bash
# AMDGPU target (ROCm >= 4.1)
export HIP_PLATFORM=amd
cmake -D PKG_GPU=on -D GPU_API=HIP -D HIP_ARCH=gfx906 -D CMAKE_CXX_COMPILER=hipcc ..
make -j 4
.. code:: bash
# CUDA target (not recommended, use GPU_ARCH=cuda)

45
doc/src/Developer_unittest.rst
View File

@ -4,10 +4,10 @@ Adding tests for unit testing
This section discusses adding or expanding tests for the unit test
infrastructure included into the LAMMPS source code distribution.
Unlike example inputs, unit tests focus on testing the "local" behavior
of individual features, tend to run very fast, and should be set up to
cover as much of the added code as possible. When contributing code to
the distribution, the LAMMPS developers will appreciate if additions
to the integrated unit test facility are included.
of individual features, tend to run fast, and should be set up to cover
as much of the added code as possible. When contributing code to the
distribution, the LAMMPS developers will appreciate if additions to the
integrated unit test facility are included.
Given the complex nature of MD simulations where many operations can
only be performed when suitable "real" simulation environment has been
@ -50,6 +50,9 @@ available:
* - File name:
- Test name:
- Description:
* - ``test_argutils.cpp``
- ArgInfo
- Tests for ``ArgInfo`` class used by LAMMPS
* - ``test_fmtlib.cpp``
- FmtLib
- Tests for ``fmtlib::`` functions used by LAMMPS
@ -155,23 +158,27 @@ have the desired effect:
{
ASSERT_EQ(lmp->update->ntimestep, 0);
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("reset_timestep 10");
if (!verbose) ::testing::internal::GetCapturedStdout();
BEGIN_HIDE_OUTPUT();
command("reset_timestep 10");
END_HIDE_OUTPUT();
ASSERT_EQ(lmp->update->ntimestep, 10);
if (!verbose) ::testing::internal::CaptureStdout();
lmp->input->one("reset_timestep 0");
if (!verbose) ::testing::internal::GetCapturedStdout();
BEGIN_HIDE_OUTPUT();
command("reset_timestep 0");
END_HIDE_OUTPUT();
ASSERT_EQ(lmp->update->ntimestep, 0);
TEST_FAILURE(".*ERROR: Timestep must be >= 0.*", command("reset_timestep -10"););
TEST_FAILURE(".*ERROR: Illegal reset_timestep .*", command("reset_timestep"););
TEST_FAILURE(".*ERROR: Illegal reset_timestep .*", command("reset_timestep 10 10"););
TEST_FAILURE(".*ERROR: Expected integer .*", command("reset_timestep xxx"););
}
Please note the use of the (global) verbose variable to control whether
the LAMMPS command will be silent by capturing the output or not. In
the default case, verbose == false, the test output will be compact and
not mixed with LAMMPS output. However setting the verbose flag (via
setting the ``TEST_ARGS`` environment variable, ``TEST_ARGS=-v``) can be
helpful to understand why tests fail unexpectedly.
Please note the use of the ``BEGIN_HIDE_OUTPUT`` and ``END_HIDE_OUTPUT``
functions that will capture output from running LAMMPS. This is normally
discarded but by setting the verbose flag (via setting the ``TEST_ARGS``
environment variable, ``TEST_ARGS=-v``) it can be printed and used to
understand why tests fail unexpectedly.
Another complexity of these tests stems from the need to capture
situations where LAMMPS will stop with an error, i.e. handle so-called
@ -210,6 +217,12 @@ The following test programs are currently available:
* - ``test_lattice_region.cpp``
- LatticeRegion
- Tests to validate the :doc:`lattice <lattice>` and :doc:`region <region>` commands
* - ``test_groups.cpp``
- GroupTest
- Tests to validate the :doc:`group <group>` command
* - ``test_variables.cpp``
- VariableTest
- Tests to validate the :doc:`variable <variable>` command
* - ``test_kim_commands.cpp``
- KimCommands
- Tests for several commands from the :ref:`KIM package <PKG-KIM>`

2
doc/src/Library_scatter.rst
View File

@ -76,7 +76,7 @@ It documents the following functions:
-----------------------
.. doxygenfunction:: lammps_create_atoms(void *handle, int n, int *id, int *type, double *x, double *v, int *image, int bexpand)
.. doxygenfunction:: lammps_create_atoms(void *handle, int n, const int *id, const int *type, const double *x, const double *v, const int *image, int bexpand)
:project: progguide

2
doc/src/Packages_standard.rst
View File

@ -71,6 +71,8 @@ package:
+----------------------------------+--------------------------------------+----------------------------------------------------+------------------------------------------------------+---------+
| :ref:`PERI <PKG-PERI>` | Peridynamics models | :doc:`pair_style peri <pair_peri>` | peri | no |
+----------------------------------+--------------------------------------+----------------------------------------------------+------------------------------------------------------+---------+
| :ref:`PLUGIN <PKG-PLUGIN>` | Plugin loader command | :doc:`plugin <plugin>` | plugins | no |
+----------------------------------+--------------------------------------+----------------------------------------------------+------------------------------------------------------+---------+
| :ref:`POEMS <PKG-POEMS>` | coupled rigid body motion | :doc:`fix poems <fix_poems>` | rigid | int |
+----------------------------------+--------------------------------------+----------------------------------------------------+------------------------------------------------------+---------+
| :ref:`PYTHON <PKG-PYTHON>` | embed Python code in an input script | :doc:`python <python>` | python | sys |

2
doc/src/Python_module.rst
View File

@ -142,7 +142,7 @@ Style Constants
Type Constants
--------------
.. py:data:: LMP_TYPE_SCALAR, LMP_TYLE_VECTOR, LMP_TYPE_ARRAY, LMP_SIZE_VECTOR, LMP_SIZE_ROWS, LMP_SIZE_COLS
.. py:data:: LMP_TYPE_SCALAR, LMP_TYPE_VECTOR, LMP_TYPE_ARRAY, LMP_SIZE_VECTOR, LMP_SIZE_ROWS, LMP_SIZE_COLS
:type: int
Constants in the :py:mod:`lammps` module to select what type of data

37
doc/src/Python_neighbor.rst
View File

@ -1,6 +1,43 @@
Neighbor list access
====================
Access to neighbor lists is handled through a couple of wrapper classes
that allows to treat it like either a python list or a NumPy array. The
access procedure is similar to that of the C-library interface: use one
of the "find" functions to look up the index of the neighbor list in the
global table of neighbor lists and then get access to the neighbor list
data. The code sample below demonstrates reading the neighbor list data
using the NumPy access method.
.. code-block:: python
from lammps import lammps
import numpy as np
lmp = lammps()
lmp.commands_string("""
region box block -2 2 -2 2 -2 2
lattice fcc 1.0
create_box 1 box
create_atoms 1 box
mass 1 1.0
pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0
run 0 post no""")
# look up the neighbor list
nlidx = lmp.find_pair_neighlist('lj/cut')
nl = lmp.numpy.get_neighlist(nlidx)
tags = lmp.extract_atom('id')
print("half neighbor list with {} entries".format(nl.size))
# print neighbor list contents
for i in range(0,nl.size):
idx, nlist = nl.get(i)
print("\natom {} with ID {} has {} neighbors:".format(idx,tags[idx],nlist.size))
if nlist.size > 0:
for n in np.nditer(nlist):
print(" atom {} with ID {}".format(n,tags[n]))
**Methods:**
* :py:meth:`lammps.get_neighlist() <lammps.lammps.get_neighlist()>`: Get neighbor list for given index

4
doc/src/fix_ave_histo.rst
View File

@ -138,8 +138,8 @@ vector or columns of the array had been listed one by one. E.g. these
.. code-block:: LAMMPS
compute myCOM all com/chunk
fix 1 all ave/histo 100 1 100 c_myCOM[*] file tmp1.com mode vector
fix 2 all ave/histo 100 1 100 c_myCOM[1] c_myCOM[2] c_myCOM[3] file tmp2.com mode vector
fix 1 all ave/histo 100 1 100 -10.0 10.0 100 c_myCOM[*] file tmp1.com mode vector
fix 2 all ave/histo 100 1 100 -10.0 10.0 100 c_myCOM[1] c_myCOM[2] c_myCOM[3] file tmp2.com mode vector
If the fix ave/histo/weight command is used, exactly two values must
be specified. If the values are vectors, they must be the same

4
doc/src/fix_wall_gran.rst
View File

@ -29,6 +29,7 @@ Syntax
gamma_t = damping coefficient for collisions in tangential direction (1/time units or 1/time-distance units - see discussion below)
xmu = static yield criterion (unitless value between 0.0 and 1.0e4)
dampflag = 0 or 1 if tangential damping force is excluded or included
optional keyword = *limit_damping*, limit damping to prevent attractive interaction
.. parsed-literal::
@ -95,7 +96,8 @@ Specifically, delta = radius - r = overlap of particle with wall, m_eff
= mass of particle, and the effective radius of contact = RiRj/Ri+Rj is
set to the radius of the particle.
The parameters *Kn*\ , *Kt*\ , *gamma_n*, *gamma_t*, *xmu* and *dampflag*
The parameters *Kn*\ , *Kt*\ , *gamma_n*, *gamma_t*, *xmu*, *dampflag*,
and the optional keyword *limit_damping*
have the same meaning and units as those specified with the
:doc:`pair_style gran/\* <pair_gran>` commands. This means a NULL can be
used for either *Kt* or *gamma_t* as described on that page. If a

3
doc/src/fix_wall_gran_region.rst
View File

@ -181,7 +181,8 @@ radius - r = overlap of particle with wall, m_eff = mass of particle,
and the effective radius of contact is just the radius of the
particle.
The parameters *Kn*\ , *Kt*\ , *gamma_n*, *gamma_t*, *xmu* and *dampflag*
The parameters *Kn*\ , *Kt*\ , *gamma_n*, *gamma_t*, *xmu*, *dampflag*,
and the optional keyword *limit_damping*
have the same meaning and units as those specified with the
:doc:`pair_style gran/\* <pair_gran>` commands. This means a NULL can be
used for either *Kt* or *gamma_t* as described on that page. If a

17
doc/src/pair_gran.rst
View File

@ -26,7 +26,7 @@ Syntax
.. code-block:: LAMMPS
pair_style style Kn Kt gamma_n gamma_t xmu dampflag
pair_style style Kn Kt gamma_n gamma_t xmu dampflag keyword
* style = *gran/hooke* or *gran/hooke/history* or *gran/hertz/history*
* Kn = elastic constant for normal particle repulsion (force/distance units or pressure units - see discussion below)
@ -36,6 +36,13 @@ Syntax
* xmu = static yield criterion (unitless value between 0.0 and 1.0e4)
* dampflag = 0 or 1 if tangential damping force is excluded or included
* keyword = *limit_damping*
.. parsed-literal::
*limit_damping* value = none
limit damping to prevent attractive interaction
.. note::
Versions of LAMMPS before 9Jan09 had different style names for
@ -54,6 +61,8 @@ Examples
pair_style gran/hooke/history 200000.0 NULL 50.0 NULL 0.5 1
pair_style gran/hooke 200000.0 70000.0 50.0 30.0 0.5 0
pair_style gran/hooke 200000.0 70000.0 50.0 30.0 0.5 0 limit_damping
Description
"""""""""""
@ -208,6 +217,12 @@ potential is used as a sub-style of :doc:`pair_style hybrid <pair_hybrid>`, then
pair_coeff command to determine which atoms interact via a granular
potential.
If two particles are moving away from each other while in contact, there
is a possibility that the particles could experience an effective attractive
force due to damping. If the *limit_damping* keyword is used, this option
will zero out the normal component of the force if there is an effective
attractive force.
----------
.. include:: accel_styles.rst

16
doc/src/pair_granular.rst
View File

@ -24,7 +24,7 @@ Examples
pair_coeff * * hooke 1000.0 50.0 tangential linear_history 500.0 1.0 0.4 damping mass_velocity
pair_style granular
pair_coeff * * hertz 1000.0 50.0 tangential mindlin 1000.0 1.0 0.4
pair_coeff * * hertz 1000.0 50.0 tangential mindlin 1000.0 1.0 0.4 limit_damping
pair_style granular
pair_coeff * * hertz/material 1e8 0.3 0.3 tangential mindlin_rescale NULL 1.0 0.4 damping tsuji
@ -623,6 +623,14 @@ Finally, the twisting torque on each particle is given by:
----------
If two particles are moving away from each other while in contact, there
is a possibility that the particles could experience an effective attractive
force due to damping. If the optional *limit_damping* keyword is used, this option
will zero out the normal component of the force if there is an effective
attractive force. This keyword cannot be used with the JKR or DMT models.
----------
The *granular* pair style can reproduce the behavior of the
*pair gran/\** styles with the appropriate settings (some very
minor differences can be expected due to corrections in
@ -657,6 +665,12 @@ then LAMMPS will use that cutoff for the specified atom type
combination, and automatically set pairwise cutoffs for the remaining
atom types.
If two particles are moving away from each other while in contact, there
is a possibility that the particles could experience an effective attractive
force due to damping. If the *limit_damping* keyword is used, this option
will zero out the normal component of the force if there is an effective
attractive force. This keyword cannot be used with the JKR or DMT models.
----------
.. include:: accel_styles.rst

9
doc/src/rerun.rst
View File

@ -15,7 +15,7 @@ Syntax
.. parsed-literal::
keyword = *first* or *last* or *every* or *skip* or *start* or *stop* or *dump*
keyword = *first* or *last* or *every* or *skip* or *start* or *stop* or *post* or *dump*
*first* args = Nfirst
Nfirst = dump timestep to start on
*last* args = Nlast
@ -28,6 +28,7 @@ Syntax
Nstart = timestep on which pseudo run will start
*stop* args = Nstop
Nstop = timestep to which pseudo run will end
*post* value = *yes* or *no*
*dump* args = same as :doc:`read_dump <read_dump>` command starting with its field arguments
Examples
@ -154,6 +155,10 @@ Also note that an error will occur if you read a snapshot from the
dump file with a timestep value larger than the *stop* setting you
have specified.
The *post* keyword can be used to minimize the output to the screen that
happens after a *rerun* command, similar to the post keyword of the
:doc:`run command <run>`. It is set to *no* by default.
The *dump* keyword is required and must be the last keyword specified.
Its arguments are passed internally to the :doc:`read_dump <read_dump>`
command. The first argument following the *dump* keyword should be
@ -226,4 +231,4 @@ Default
The option defaults are first = 0, last = a huge value (effectively
infinity), start = same as first, stop = same as last, every = 0, skip
= 1;
= 1, post = no;

1
doc/utils/sphinx-config/false_positives.txt
View File

@ -2693,6 +2693,7 @@ representable
Reproducibility
reproducibility
repuls
reqid
rescale
rescaled
rescales

56
examples/deposit/in.deposit.molecule.rigid-nve-small Normal file
View File

@ -0,0 +1,56 @@
# sample surface deposition script for molecules
units lj
atom_style bond
boundary p p f
lattice fcc 1.0
region box block 0 5 0 5 0 10
create_box 3 box bond/types 1 extra/bond/per/atom 1
region substrate block INF INF INF INF INF 3
create_atoms 1 region substrate
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0
pair_coeff 1 2 1.0 1.0 5.0
mass * 1.0
bond_style harmonic
bond_coeff 1 5.0 1.0
neigh_modify delay 0
molecule dimer molecule.dimer
region slab block 0 5 0 5 8 9
group addatoms empty
region mobile block 0 5 0 5 2 INF
group mobile region mobile
compute add addatoms temp
compute_modify add dynamic/dof yes extra/dof 0
fix 1 addatoms rigid/nve/small molecule mol dimer
fix 2 mobile langevin 0.1 0.1 0.1 587283
fix 3 mobile nve
fix 4 addatoms deposit 100 0 100 12345 region slab near 1.0 &
mol dimer vz -1.0 -1.0 rigid 1
fix 5 addatoms wall/reflect zhi EDGE
thermo_style custom step atoms temp epair etotal press
thermo 100
thermo_modify temp add lost/bond ignore lost warn
#dump 1 all atom 50 dump.deposit.atom
#dump 2 all image 50 image.*.jpg type type &
# axes yes 0.8 0.02 view 80 -30
#dump_modify 2 pad 5
#dump 3 all movie 50 tmp.mpg type type &
# axes yes 0.8 0.02 view 80 -30
#dump_modify 3 pad 5
run 10000

56
examples/deposit/in.deposit.molecule.rigid-nvt-small Normal file
View File

@ -0,0 +1,56 @@
# sample surface deposition script for molecules
units lj
atom_style bond
boundary p p f
lattice fcc 1.0
region box block 0 5 0 5 0 10
create_box 3 box bond/types 1 extra/bond/per/atom 1
region substrate block INF INF INF INF INF 3
create_atoms 1 region substrate
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0
pair_coeff 1 2 1.0 1.0 5.0
mass * 1.0
bond_style harmonic
bond_coeff 1 5.0 1.0
neigh_modify delay 0
molecule dimer molecule.dimer
region slab block 0 5 0 5 8 9
group addatoms empty
region mobile block 0 5 0 5 2 INF
group mobile region mobile
compute add addatoms temp
compute_modify add dynamic/dof yes extra/dof 0
fix 1 addatoms rigid/nvt/small molecule temp 0.1 0.1 0.1 mol dimer
fix 2 mobile langevin 0.1 0.1 0.1 587283
fix 3 mobile nve
fix 4 addatoms deposit 100 0 100 12345 region slab near 1.0 &
mol dimer vz -1.0 -1.0 rigid 1
fix 5 addatoms wall/reflect zhi EDGE
thermo_style custom step atoms temp epair etotal press
thermo 100
thermo_modify temp add lost/bond ignore lost warn
#dump 1 all atom 50 dump.deposit.atom
#dump 2 all image 50 image.*.jpg type type &
# axes yes 0.8 0.02 view 80 -30
#dump_modify 2 pad 5
#dump 3 all movie 50 tmp.mpg type type &
# axes yes 0.8 0.02 view 80 -30
#dump_modify 3 pad 5
run 10000

218
examples/deposit/log.10Mar21.deposit.molecule.rigid-nve-small.g++.1 Normal file
View File

@ -0,0 +1,218 @@
LAMMPS (10 Mar 2021)
using 1 OpenMP thread(s) per MPI task
# sample surface deposition script for molecules
units lj
atom_style bond
boundary p p f
lattice fcc 1.0
Lattice spacing in x,y,z = 1.5874011 1.5874011 1.5874011
region box block 0 5 0 5 0 10
create_box 3 box bond/types 1 extra/bond/per/atom 1
Created orthogonal box = (0.0000000 0.0000000 0.0000000) to (7.9370053 7.9370053 15.874011)
1 by 1 by 1 MPI processor grid
region substrate block INF INF INF INF INF 3
create_atoms 1 region substrate
Created 350 atoms
create_atoms CPU = 0.001 seconds
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0
pair_coeff 1 2 1.0 1.0 5.0
mass * 1.0
bond_style harmonic
bond_coeff 1 5.0 1.0
neigh_modify delay 0
molecule dimer molecule.dimer
Read molecule template dimer:
1 molecules
2 atoms with max type 3
1 bonds with max type 1
0 angles with max type 0
0 dihedrals with max type 0
0 impropers with max type 0
region slab block 0 5 0 5 8 9
group addatoms empty
0 atoms in group addatoms
region mobile block 0 5 0 5 2 INF
group mobile region mobile
150 atoms in group mobile
compute add addatoms temp
compute_modify add dynamic/dof yes extra/dof 0
fix 1 addatoms rigid/nve/small molecule mol dimer
create bodies CPU = 0.000 seconds
0 rigid bodies with 0 atoms
1.0000000 = max distance from body owner to body atom
fix 2 mobile langevin 0.1 0.1 0.1 587283
fix 3 mobile nve
fix 4 addatoms deposit 100 0 100 12345 region slab near 1.0 mol dimer vz -1.0 -1.0 rigid 1
fix 5 addatoms wall/reflect zhi EDGE
thermo_style custom step atoms temp epair etotal press
thermo 100
thermo_modify temp add lost/bond ignore lost warn
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:468)
#dump 1 all atom 50 dump.deposit.atom
#dump 2 all image 50 image.*.jpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 2 pad 5
#dump 3 all movie 50 tmp.mpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 3 pad 5
run 10000
WARNING: Should not allow rigid bodies to bounce off relecting walls (src/fix_wall_reflect.cpp:182)
Neighbor list info ...
update every 1 steps, delay 0 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 5.3
ghost atom cutoff = 5.3
binsize = 2.65, bins = 3 3 6
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair lj/cut, perpetual
attributes: half, newton on
pair build: half/bin/newton
stencil: half/bin/3d/newton
bin: standard
Per MPI rank memory allocation (min/avg/max) = 5.929 | 5.929 | 5.929 Mbytes
Step Atoms Temp E_pair TotEng Press
0 350 0 -6.9215833 -6.9215833 -1.0052629
100 352 1.0079368 -6.8875167 -6.8803581 -0.73353914
200 354 1.0081552 -6.8594643 -6.8452248 -0.70421276
300 356 1.0085803 -6.8171524 -6.7959042 -0.6917826
400 358 1.0099188 -6.7852701 -6.7570601 -0.70371699
500 360 1.0140221 -6.7493429 -6.7141338 -0.68415307
600 362 1.026148 -6.7105231 -6.6680032 -0.68314418
700 364 1.0683344 -6.6725162 -6.621154 -0.65747369
800 366 1.0958953 -6.6412275 -6.5813425 -0.68789615
900 368 1.1250034 -6.6101882 -6.541404 -0.66674346
1000 370 1.2326365 -6.5993719 -6.5160856 -0.69688655
1100 372 1.1397444 -6.5912865 -6.5070312 -0.6333041
1200 374 1.0514447 -6.5905753 -6.5062348 -0.71020272
1300 376 1.0033284 -6.5747792 -6.4880554 -0.65459641
1400 378 0.82996548 -6.5681806 -6.4913319 -0.60438595
1500 380 0.90239071 -6.5752982 -6.4862465 -0.66528503
1600 382 0.86403681 -6.5692212 -6.4787461 -0.65027922
1700 384 0.64748046 -6.5644242 -6.4927629 -0.63046709
1800 386 0.74288742 -6.5515741 -6.4649681 -0.67770665
1900 388 0.72584537 -6.5565195 -6.4676596 -0.66175025
2000 390 0.73351281 -6.5631154 -6.4690753 -0.64686647
2100 392 0.76490681 -6.5573734 -6.4549305 -0.68861526
2200 394 0.65926638 -6.5511574 -6.4591279 -0.71726284
2300 396 0.70267414 -6.5728555 -6.4708258 -0.64360436
2400 398 0.60523691 -6.5845973 -6.4933555 -0.63956839
2500 400 0.5103586 -6.5812006 -6.5014571 -0.68698448
2600 402 0.52401744 -6.6003358 -6.5156066 -0.68987222
2700 404 0.46421291 -6.5662747 -6.4887143 -0.72872856
2800 406 0.48254258 -6.5724266 -6.4892296 -0.75447872
2900 408 0.53073083 -6.5809842 -6.4866754 -0.67592283
3000 410 0.55314547 -6.5922077 -6.4910226 -0.74646953
3100 412 0.47150308 -6.5907273 -6.5020344 -0.64994935
3200 414 0.43047803 -6.5836315 -6.5004473 -0.79764713
3300 416 0.46289353 -6.5739439 -6.4821441 -0.76441674
3400 418 0.59724106 -6.5980575 -6.4766089 -0.73735626
3500 420 0.43571285 -6.5955972 -6.5048237 -0.64145941
3600 422 0.42461639 -6.6060271 -6.5154691 -0.70124484
3700 424 0.44323254 -6.6059723 -6.5092766 -0.74498147
3800 426 0.4037907 -6.592043 -6.5019958 -0.72171799
3900 428 0.41668443 -6.5975302 -6.5026079 -0.68327878
4000 430 0.44895183 -6.5958671 -6.4914597 -0.73939433
4100 432 0.39798052 -6.5952115 -6.5007833 -0.74871822
4200 434 0.45156734 -6.6237274 -6.5144772 -0.75111778
4300 436 0.48297356 -6.6395283 -6.5204465 -0.76105913
4400 438 0.43595234 -6.6347133 -6.5252276 -0.72858275
4500 440 0.44683726 -6.6385398 -6.5242916 -0.74280051
4600 442 0.47875512 -6.6419903 -6.5174274 -0.75579572
4700 444 0.43972605 -6.6406078 -6.5242388 -0.72196509
4800 446 0.43572615 -6.6495404 -6.5323047 -0.7135049
4900 448 0.38063437 -6.6432385 -6.5391588 -0.77087429
5000 450 0.4239223 -6.6617795 -6.5440233 -0.784531
5100 452 0.37201988 -6.6581813 -6.5532421 -0.74611403
5200 454 0.41765777 -6.6661321 -6.5465385 -0.75422239
5300 456 0.38015287 -6.6606624 -6.5502013 -0.78866702
5400 458 0.40549607 -6.6807118 -6.5611878 -0.78932883
5500 460 0.34444407 -6.6720564 -6.5690976 -0.77859171
5600 462 0.36572308 -6.6730078 -6.5621827 -0.85672419
5700 464 0.40055073 -6.6976989 -6.574685 -0.74251563
5800 466 0.36037213 -6.7022014 -6.5900685 -0.66239625
5900 468 0.32810921 -6.6952135 -6.591803 -0.83981757
6000 470 0.35110886 -6.6986862 -6.5866302 -0.82474047
6100 472 0.29965884 -6.6839503 -6.5871326 -0.7864913
6200 474 0.32402637 -6.6902745 -6.5843165 -0.74531083
6300 476 0.35042653 -6.6990084 -6.5830585 -0.74839967
6400 478 0.32695511 -6.6909459 -6.5815048 -0.76549489
6500 480 0.35209088 -6.6902987 -6.5711013 -0.71281516
6600 482 0.354106 -6.6890268 -6.567808 -0.81897158
6700 484 0.3504816 -6.681739 -6.5604463 -0.811609
6800 486 0.37396733 -6.7018369 -6.5710253 -0.80383296
6900 488 0.36435774 -6.7010114 -6.5722169 -0.72063651
7000 490 0.35631012 -6.7089806 -6.581727 -0.74152078
7100 492 0.37646659 -6.719154 -6.5833352 -0.77739001
7200 494 0.36546269 -6.7223269 -6.5891623 -0.8288767
7300 496 0.37688206 -6.7457243 -6.607053 -0.80062121
7400 498 0.30331409 -6.7284953 -6.6158184 -0.79894584
7500 500 0.30382936 -6.7333804 -6.6194444 -0.86924602
7600 502 0.30163143 -6.7294737 -6.6153104 -0.7689497
7700 504 0.30281215 -6.7233976 -6.6077402 -0.8557548
7800 506 0.33378009 -6.7244958 -6.5958651 -0.82584084
7900 508 0.31843128 -6.7250998 -6.6013002 -0.82424684
8000 510 0.35912946 -6.7399052 -6.5990701 -0.81879575
8100 512 0.31405017 -6.733487 -6.6092777 -0.72940457
8200 514 0.30475 -6.7271485 -6.6056043 -0.87958287
8300 516 0.30717434 -6.7277574 -6.6042329 -0.81700937
8400 518 0.33032462 -6.7382165 -6.6043011 -0.75436496
8500 520 0.32218568 -6.7351971 -6.6035347 -0.77661738
8600 522 0.30922651 -6.7275431 -6.6001797 -0.85334327
8700 524 0.30728021 -6.7237477 -6.5962029 -0.94977016
8800 526 0.31061903 -6.7181672 -6.5882505 -0.86132456
8900 528 0.33543344 -6.728929 -6.5875768 -0.87545919
9000 530 0.31887735 -6.7265066 -6.5911341 -0.76892061
9100 532 0.31888326 -6.7216274 -6.5852629 -0.83190298
9200 534 0.33310892 -6.7111349 -6.567661 -0.94282671
9300 536 0.34737171 -6.722515 -6.5718361 -0.95235602
9400 538 0.32752858 -6.7246204 -6.5815549 -0.86227131
9500 540 0.30665764 -6.7225391 -6.5876665 -0.87144326
9600 542 0.32747382 -6.7149245 -6.5699176 -0.86863105
9700 544 0.32463079 -6.7205757 -6.5758643 -0.85393932
9800 546 0.31517825 -6.7178961 -6.5764699 -0.81017759
9900 548 0.33649933 -6.7380644 -6.5860871 -0.80769312
10000 550 0.37394555 -6.7612874 -6.5913121 -0.82102213
Loop time of 16.7275 on 1 procs for 10000 steps with 550 atoms
Performance: 258256.688 tau/day, 597.816 timesteps/s
99.9% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 9.2508 | 9.2508 | 9.2508 | 0.0 | 55.30
Bond | 0.021603 | 0.021603 | 0.021603 | 0.0 | 0.13
Neigh | 4.5325 | 4.5325 | 4.5325 | 0.0 | 27.10
Comm | 0.62777 | 0.62777 | 0.62777 | 0.0 | 3.75
Output | 0.006916 | 0.006916 | 0.006916 | 0.0 | 0.04
Modify | 2.218 | 2.218 | 2.218 | 0.0 | 13.26
Other | | 0.07002 | | | 0.42
Nlocal: 550.000 ave 550 max 550 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 2344.00 ave 2344 max 2344 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 36607.0 ave 36607 max 36607 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 36607
Ave neighs/atom = 66.558182
Ave special neighs/atom = 0.36363636
Neighbor list builds = 847
Dangerous builds = 0
Total wall time: 0:00:16

218
examples/deposit/log.10Mar21.deposit.molecule.rigid-nve-small.g++.4 Normal file
View File

@ -0,0 +1,218 @@
LAMMPS (10 Mar 2021)
using 1 OpenMP thread(s) per MPI task
# sample surface deposition script for molecules
units lj
atom_style bond
boundary p p f
lattice fcc 1.0
Lattice spacing in x,y,z = 1.5874011 1.5874011 1.5874011
region box block 0 5 0 5 0 10
create_box 3 box bond/types 1 extra/bond/per/atom 1
Created orthogonal box = (0.0000000 0.0000000 0.0000000) to (7.9370053 7.9370053 15.874011)
1 by 1 by 4 MPI processor grid
region substrate block INF INF INF INF INF 3
create_atoms 1 region substrate
Created 350 atoms
create_atoms CPU = 0.002 seconds
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0
pair_coeff 1 2 1.0 1.0 5.0
mass * 1.0
bond_style harmonic
bond_coeff 1 5.0 1.0
neigh_modify delay 0
molecule dimer molecule.dimer
Read molecule template dimer:
1 molecules
2 atoms with max type 3
1 bonds with max type 1
0 angles with max type 0
0 dihedrals with max type 0
0 impropers with max type 0
region slab block 0 5 0 5 8 9
group addatoms empty
0 atoms in group addatoms
region mobile block 0 5 0 5 2 INF
group mobile region mobile
150 atoms in group mobile
compute add addatoms temp
compute_modify add dynamic/dof yes extra/dof 0
fix 1 addatoms rigid/nve/small molecule mol dimer
create bodies CPU = 0.000 seconds
0 rigid bodies with 0 atoms
1.0000000 = max distance from body owner to body atom
fix 2 mobile langevin 0.1 0.1 0.1 587283
fix 3 mobile nve
fix 4 addatoms deposit 100 0 100 12345 region slab near 1.0 mol dimer vz -1.0 -1.0 rigid 1
fix 5 addatoms wall/reflect zhi EDGE
thermo_style custom step atoms temp epair etotal press
thermo 100
thermo_modify temp add lost/bond ignore lost warn
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:468)
#dump 1 all atom 50 dump.deposit.atom
#dump 2 all image 50 image.*.jpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 2 pad 5
#dump 3 all movie 50 tmp.mpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 3 pad 5
run 10000
WARNING: Should not allow rigid bodies to bounce off relecting walls (src/fix_wall_reflect.cpp:182)
Neighbor list info ...
update every 1 steps, delay 0 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 5.3
ghost atom cutoff = 5.3
binsize = 2.65, bins = 3 3 6
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair lj/cut, perpetual
attributes: half, newton on
pair build: half/bin/newton
stencil: half/bin/3d/newton
bin: standard
Per MPI rank memory allocation (min/avg/max) = 5.255 | 5.852 | 6.302 Mbytes
Step Atoms Temp E_pair TotEng Press
0 350 0 -6.9215833 -6.9215833 -1.0052629
100 352 1.0079368 -6.8946578 -6.8874992 -0.73775337
200 354 1.0081552 -6.8645575 -6.850318 -0.69629729
300 356 1.0085803 -6.821677 -6.8004288 -0.69532657
400 358 1.0099188 -6.7837923 -6.7555822 -0.68879568
500 360 1.0140221 -6.7446709 -6.7094618 -0.72991641
600 362 1.026146 -6.7129201 -6.6704003 -0.67063836
700 364 1.0683193 -6.6776523 -6.6262908 -0.65572472
800 366 1.0958894 -6.6402029 -6.5803182 -0.66307282
900 368 1.1231769 -6.6050912 -6.5364187 -0.64076928
1000 370 1.1976289 -6.5942508 -6.51333 -0.69249047
1100 372 1.05061 -6.5772306 -6.4995645 -0.6307289
1200 374 0.93723829 -6.5732962 -6.4981166 -0.64963821
1300 376 0.93899122 -6.5578418 -6.476679 -0.65096688
1400 378 0.87050694 -6.546852 -6.4662495 -0.67613401
1500 380 0.84462904 -6.5400883 -6.4567368 -0.64376986
1600 382 0.92417795 -6.5611292 -6.4643567 -0.62092779
1700 384 0.67296541 -6.5589985 -6.4845167 -0.6779111
1800 386 0.87159158 -6.5655058 -6.4638954 -0.6582786
1900 388 0.67504975 -6.5772537 -6.4946123 -0.665098
2000 390 0.65034978 -6.5717854 -6.4884072 -0.73162072
2100 392 0.63389859 -6.5652907 -6.4803935 -0.64937411
2200 394 0.6265637 -6.5582412 -6.4707767 -0.70819854
2300 396 0.85352983 -6.5750751 -6.4511408 -0.69995918
2400 398 0.69792348 -6.5938659 -6.4886513 -0.63755214
2500 400 0.45320495 -6.595022 -6.5242087 -0.70462183
2600 402 0.46824374 -6.5697613 -6.4940502 -0.75497893
2700 404 0.46605975 -6.5599272 -6.4820583 -0.66778506
2800 406 0.50637604 -6.5700492 -6.482743 -0.63812082
2900 408 0.55335921 -6.5873352 -6.4890054 -0.62891652
3000 410 0.51731142 -6.5890566 -6.4944264 -0.64305355
3100 412 0.59439734 -6.6003127 -6.4885025 -0.67140541
3200 414 0.53413632 -6.595564 -6.4923492 -0.76242248
3300 416 0.49221593 -6.5990831 -6.5014681 -0.71330819
3400 418 0.55006126 -6.6015725 -6.4897179 -0.72044309
3500 420 0.49065348 -6.6329864 -6.5307669 -0.65775397
3600 422 0.41907335 -6.6219753 -6.5325995 -0.75936391
3700 424 0.38720116 -6.6153349 -6.5308629 -0.74166397
3800 426 0.40625994 -6.6150209 -6.5244231 -0.7595111
3900 428 0.40460547 -6.6122642 -6.5200936 -0.70484465
4000 430 0.45014991 -6.6254404 -6.5207544 -0.70069933
4100 432 0.44820466 -6.640222 -6.5338771 -0.805972
4200 434 0.39767521 -6.6450316 -6.5488199 -0.72841414
4300 436 0.4155331 -6.6547917 -6.5523381 -0.6894484
4400 438 0.43034353 -6.6615074 -6.5534303 -0.74026062
4500 440 0.38062977 -6.6541618 -6.5568417 -0.85911194
4600 442 0.39357419 -6.6522517 -6.5498512 -0.66232114
4700 444 0.40296801 -6.6647029 -6.5580616 -0.67307577
4800 446 0.38194993 -6.6510248 -6.548258 -0.77887746
4900 448 0.40739835 -6.6601751 -6.5487771 -0.84146416
5000 450 0.38392807 -6.6560665 -6.5494198 -0.77343399
5100 452 0.35209286 -6.6481476 -6.5488294 -0.68065755
5200 454 0.41355143 -6.6615988 -6.543181 -0.81611092
5300 456 0.42200484 -6.6627494 -6.5401273 -0.7774134
5400 458 0.37764326 -6.661431 -6.550117 -0.72187808
5500 460 0.41857405 -6.6674232 -6.542306 -0.74929745
5600 462 0.41682611 -6.6710961 -6.5447852 -0.7557454
5700 464 0.44396148 -6.6924346 -6.5560887 -0.78018312
5800 466 0.37058077 -6.6865013 -6.5711919 -0.74146125
5900 468 0.346812 -6.681215 -6.57191 -0.81184233
6000 470 0.34919462 -6.6827931 -6.571348 -0.87330821
6100 472 0.39360936 -6.6933359 -6.5661634 -0.79237598
6200 474 0.33270778 -6.6847095 -6.5759127 -0.77978526
6300 476 0.35973804 -6.6951254 -6.5760944 -0.80340174
6400 478 0.38124318 -6.7045984 -6.5769857 -0.81628407
6500 480 0.41188302 -6.7133356 -6.573896 -0.7940289
6600 482 0.36998039 -6.7079555 -6.5813025 -0.75055442
6700 484 0.37615026 -6.722917 -6.592741 -0.76534055
6800 486 0.3466597 -6.7188712 -6.5976116 -0.77986211
6900 488 0.32902492 -6.7247054 -6.6084005 -0.75702458
7000 490 0.31856427 -6.7167709 -6.6029979 -0.76014555
7100 492 0.30233891 -6.7144406 -6.6053651 -0.96246708
7200 494 0.32557309 -6.7315347 -6.6129049 -0.8122153
7300 496 0.29919611 -6.7186327 -6.6085455 -0.71917931
7400 498 0.29778169 -6.7259455 -6.6153238 -0.88199391
7500 500 0.35535305 -6.7419073 -6.6086499 -0.90344083
7600 502 0.31979187 -6.7326802 -6.6116434 -0.78324572
7700 504 0.30821359 -6.72895 -6.6112295 -0.81363335
7800 506 0.31374993 -6.7255461 -6.6046346 -0.89904197
7900 508 0.29072338 -6.7177355 -6.6047082 -0.81423073
8000 510 0.30557494 -6.7208283 -6.600995 -0.85853575
8100 512 0.30521237 -6.726874 -6.6061601 -0.75361257
8200 514 0.29622226 -6.7152225 -6.5970794 -0.85766132
8300 516 0.28337698 -6.7023552 -6.5884003 -0.8985415
8400 518 0.32860902 -6.7211074 -6.5878875 -0.89758921
8500 520 0.35483743 -6.7406183 -6.5956126 -0.77602077
8600 522 0.32928486 -6.7326607 -6.5970358 -0.75070137
8700 524 0.36008106 -6.7474714 -6.5980103 -0.87093836
8800 526 0.36082301 -6.743579 -6.5926644 -0.90132107
8900 528 0.35010285 -6.7501553 -6.6026214 -0.85238959
9000 530 0.31513985 -6.7411795 -6.6073937 -0.75884529
9100 532 0.30895083 -6.7421063 -6.6099891 -0.74482692
9200 534 0.33631849 -6.751265 -6.6064087 -0.95632911
9300 536 0.33306096 -6.759187 -6.6147156 -0.76275935
9400 538 0.34582537 -6.7706766 -6.619619 -0.72251598
9500 540 0.32003146 -6.772057 -6.6313024 -0.88195851
9600 542 0.32439637 -6.7743569 -6.6307127 -0.90233104
9700 544 0.33988235 -6.7763721 -6.624862 -0.85185581
9800 546 0.32877587 -6.7744977 -6.62697 -0.8550905
9900 548 0.29570051 -6.7623752 -6.6288243 -1.0371157
10000 550 0.33675914 -6.7757315 -6.6226591 -1.0082157
Loop time of 16.2553 on 4 procs for 10000 steps with 550 atoms
Performance: 265759.101 tau/day, 615.183 timesteps/s
97.0% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.0034051 | 2.6646 | 8.0317 | 200.7 | 16.39
Bond | 0.0022626 | 0.0072971 | 0.019437 | 8.3 | 0.04
Neigh | 0.018657 | 1.218 | 3.7621 | 137.5 | 7.49
Comm | 1.1497 | 6.7957 | 12.192 | 164.2 | 41.81
Output | 0.010692 | 0.01431 | 0.019901 | 3.2 | 0.09
Modify | 3.1502 | 5.3779 | 11.342 | 149.0 | 33.08
Other | | 0.1775 | | | 1.09
Nlocal: 137.500 ave 299 max 2 min
Histogram: 2 0 0 0 0 0 0 1 0 1
Nghost: 1903.75 ave 2686 max 529 min
Histogram: 1 0 0 0 0 0 1 0 0 2
Neighs: 9209.75 ave 23206 max 0 min
Histogram: 2 0 0 0 0 1 0 0 0 1
Total # of neighbors = 36839
Ave neighs/atom = 66.980000
Ave special neighs/atom = 0.36363636
Neighbor list builds = 829
Dangerous builds = 0
Total wall time: 0:00:16

218
examples/deposit/log.10Mar21.deposit.molecule.rigid-nvt-small.g++.1 Normal file
View File

@ -0,0 +1,218 @@
LAMMPS (10 Mar 2021)
using 1 OpenMP thread(s) per MPI task
# sample surface deposition script for molecules
units lj
atom_style bond
boundary p p f
lattice fcc 1.0
Lattice spacing in x,y,z = 1.5874011 1.5874011 1.5874011
region box block 0 5 0 5 0 10
create_box 3 box bond/types 1 extra/bond/per/atom 1
Created orthogonal box = (0.0000000 0.0000000 0.0000000) to (7.9370053 7.9370053 15.874011)
1 by 1 by 1 MPI processor grid
region substrate block INF INF INF INF INF 3
create_atoms 1 region substrate
Created 350 atoms
create_atoms CPU = 0.001 seconds
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0
pair_coeff 1 2 1.0 1.0 5.0
mass * 1.0
bond_style harmonic
bond_coeff 1 5.0 1.0
neigh_modify delay 0
molecule dimer molecule.dimer
Read molecule template dimer:
1 molecules
2 atoms with max type 3
1 bonds with max type 1
0 angles with max type 0
0 dihedrals with max type 0
0 impropers with max type 0
region slab block 0 5 0 5 8 9
group addatoms empty
0 atoms in group addatoms
region mobile block 0 5 0 5 2 INF
group mobile region mobile
150 atoms in group mobile
compute add addatoms temp
compute_modify add dynamic/dof yes extra/dof 0
fix 1 addatoms rigid/nvt/small molecule temp 0.1 0.1 0.1 mol dimer
create bodies CPU = 0.000 seconds
0 rigid bodies with 0 atoms
1.0000000 = max distance from body owner to body atom
fix 2 mobile langevin 0.1 0.1 0.1 587283
fix 3 mobile nve
fix 4 addatoms deposit 100 0 100 12345 region slab near 1.0 mol dimer vz -1.0 -1.0 rigid 1
fix 5 addatoms wall/reflect zhi EDGE
thermo_style custom step atoms temp epair etotal press
thermo 100
thermo_modify temp add lost/bond ignore lost warn
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:468)
#dump 1 all atom 50 dump.deposit.atom
#dump 2 all image 50 image.*.jpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 2 pad 5
#dump 3 all movie 50 tmp.mpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 3 pad 5
run 10000
WARNING: Should not allow rigid bodies to bounce off relecting walls (src/fix_wall_reflect.cpp:182)
Neighbor list info ...
update every 1 steps, delay 0 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 5.3
ghost atom cutoff = 5.3
binsize = 2.65, bins = 3 3 6
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair lj/cut, perpetual
attributes: half, newton on
pair build: half/bin/newton
stencil: half/bin/3d/newton
bin: standard
Per MPI rank memory allocation (min/avg/max) = 5.929 | 5.929 | 5.929 Mbytes
Step Atoms Temp E_pair TotEng Press
0 350 0 -6.9215833 -6.9215833 -1.0052629
100 352 1.0079368 -6.8875167 -6.8803581 -0.73353914
200 354 1.0081552 -6.8594643 -6.8452248 -0.70421276
300 356 1.0085803 -6.8171524 -6.7959042 -0.6917826
400 358 1.0099188 -6.7852701 -6.7570601 -0.70371699
500 360 1.0140221 -6.7493429 -6.7141338 -0.68415307
600 362 1.026148 -6.7105231 -6.6680032 -0.68314418
700 364 1.0683344 -6.6725162 -6.621154 -0.65747369
800 366 1.0958953 -6.6412275 -6.5813425 -0.68789615
900 368 1.1250034 -6.6101882 -6.541404 -0.66674346
1000 370 1.2326365 -6.5993719 -6.5160856 -0.69688655
1100 372 1.1397444 -6.5912865 -6.5070312 -0.6333041
1200 374 1.0514447 -6.5905753 -6.5062348 -0.71020272
1300 376 1.0033284 -6.5747792 -6.4880554 -0.65459641
1400 378 0.82996548 -6.5681806 -6.4913319 -0.60438595
1500 380 0.90239071 -6.5752982 -6.4862465 -0.66528503
1600 382 0.86403681 -6.5692212 -6.4787461 -0.65027922
1700 384 0.64748046 -6.5644242 -6.4927629 -0.63046709
1800 386 0.74288742 -6.5515741 -6.4649681 -0.67770665
1900 388 0.72584537 -6.5565195 -6.4676596 -0.66175025
2000 390 0.73351281 -6.5631154 -6.4690753 -0.64686647
2100 392 0.76490681 -6.5573734 -6.4549305 -0.68861526
2200 394 0.65926638 -6.5511574 -6.4591279 -0.71726284
2300 396 0.70267414 -6.5728555 -6.4708258 -0.64360436
2400 398 0.60523691 -6.5845973 -6.4933555 -0.63956839
2500 400 0.5103586 -6.5812006 -6.5014571 -0.68698448
2600 402 0.52401744 -6.6003358 -6.5156066 -0.68987222
2700 404 0.46421291 -6.5662747 -6.4887143 -0.72872856
2800 406 0.48254258 -6.5724266 -6.4892296 -0.75447872
2900 408 0.53073083 -6.5809842 -6.4866754 -0.67592283
3000 410 0.55314547 -6.5922077 -6.4910226 -0.74646953
3100 412 0.47150308 -6.5907273 -6.5020344 -0.64994935
3200 414 0.43047803 -6.5836315 -6.5004473 -0.79764713
3300 416 0.46289353 -6.5739439 -6.4821441 -0.76441674
3400 418 0.59724106 -6.5980575 -6.4766089 -0.73735626
3500 420 0.43571285 -6.5955972 -6.5048237 -0.64145941
3600 422 0.42461639 -6.6060271 -6.5154691 -0.70124484
3700 424 0.44323254 -6.6059723 -6.5092766 -0.74498147
3800 426 0.4037907 -6.592043 -6.5019958 -0.72171799
3900 428 0.41668443 -6.5975302 -6.5026079 -0.68327878
4000 430 0.44895183 -6.5958671 -6.4914597 -0.73939433
4100 432 0.39798052 -6.5952115 -6.5007833 -0.74871822
4200 434 0.45156734 -6.6237274 -6.5144772 -0.75111778
4300 436 0.48297356 -6.6395283 -6.5204465 -0.76105913
4400 438 0.43595234 -6.6347133 -6.5252276 -0.72858275
4500 440 0.44683726 -6.6385398 -6.5242916 -0.74280051
4600 442 0.47875512 -6.6419903 -6.5174274 -0.75579572
4700 444 0.43972605 -6.6406078 -6.5242388 -0.72196509
4800 446 0.43572615 -6.6495404 -6.5323047 -0.7135049
4900 448 0.38063437 -6.6432385 -6.5391588 -0.77087429
5000 450 0.4239223 -6.6617795 -6.5440233 -0.784531
5100 452 0.37201988 -6.6581813 -6.5532421 -0.74611403
5200 454 0.41765777 -6.6661321 -6.5465385 -0.75422239
5300 456 0.38015287 -6.6606624 -6.5502013 -0.78866702
5400 458 0.40549607 -6.6807118 -6.5611878 -0.78932883
5500 460 0.34444407 -6.6720564 -6.5690976 -0.77859171
5600 462 0.36572308 -6.6730078 -6.5621827 -0.85672419
5700 464 0.40055073 -6.6976989 -6.574685 -0.74251563
5800 466 0.36037213 -6.7022014 -6.5900685 -0.66239625
5900 468 0.32810921 -6.6952135 -6.591803 -0.83981757
6000 470 0.35110886 -6.6986862 -6.5866302 -0.82474047
6100 472 0.29965884 -6.6839503 -6.5871326 -0.7864913
6200 474 0.32402637 -6.6902745 -6.5843165 -0.74531083
6300 476 0.35042653 -6.6990084 -6.5830585 -0.74839967
6400 478 0.32695511 -6.6909459 -6.5815048 -0.76549489
6500 480 0.35209088 -6.6902987 -6.5711013 -0.71281516
6600 482 0.354106 -6.6890268 -6.567808 -0.81897158
6700 484 0.3504816 -6.681739 -6.5604463 -0.811609
6800 486 0.37396733 -6.7018369 -6.5710253 -0.80383296
6900 488 0.36435774 -6.7010114 -6.5722169 -0.72063651
7000 490 0.35631012 -6.7089806 -6.581727 -0.74152078
7100 492 0.37646659 -6.719154 -6.5833352 -0.77739001
7200 494 0.36546269 -6.7223269 -6.5891623 -0.8288767
7300 496 0.37688206 -6.7457243 -6.607053 -0.80062121
7400 498 0.30331409 -6.7284953 -6.6158184 -0.79894584
7500 500 0.30382936 -6.7333804 -6.6194444 -0.86924602
7600 502 0.30163143 -6.7294737 -6.6153104 -0.7689497
7700 504 0.30281215 -6.7233976 -6.6077402 -0.8557548
7800 506 0.33378009 -6.7244958 -6.5958651 -0.82584084
7900 508 0.31843128 -6.7250998 -6.6013002 -0.82424684
8000 510 0.35912946 -6.7399052 -6.5990701 -0.81879575
8100 512 0.31405017 -6.733487 -6.6092777 -0.72940457
8200 514 0.30475 -6.7271485 -6.6056043 -0.87958287
8300 516 0.30717434 -6.7277574 -6.6042329 -0.81700937
8400 518 0.33032462 -6.7382165 -6.6043011 -0.75436496
8500 520 0.32218568 -6.7351971 -6.6035347 -0.77661738
8600 522 0.30922651 -6.7275431 -6.6001797 -0.85334327
8700 524 0.30728021 -6.7237477 -6.5962029 -0.94977016
8800 526 0.31061903 -6.7181672 -6.5882505 -0.86132456
8900 528 0.33543344 -6.728929 -6.5875768 -0.87545919
9000 530 0.31887735 -6.7265066 -6.5911341 -0.76892061
9100 532 0.31888326 -6.7216274 -6.5852629 -0.83190298
9200 534 0.33310892 -6.7111349 -6.567661 -0.94282671
9300 536 0.34737171 -6.722515 -6.5718361 -0.95235602
9400 538 0.32752858 -6.7246204 -6.5815549 -0.86227131
9500 540 0.30665764 -6.7225391 -6.5876665 -0.87144326
9600 542 0.32747382 -6.7149245 -6.5699176 -0.86863105
9700 544 0.32463079 -6.7205757 -6.5758643 -0.85393932
9800 546 0.31517825 -6.7178961 -6.5764699 -0.81017759
9900 548 0.33649933 -6.7380644 -6.5860871 -0.80769312
10000 550 0.37394555 -6.7612874 -6.5913121 -0.82102213
Loop time of 16.6615 on 1 procs for 10000 steps with 550 atoms
Performance: 259279.806 tau/day, 600.185 timesteps/s
99.9% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 9.2639 | 9.2639 | 9.2639 | 0.0 | 55.60
Bond | 0.021559 | 0.021559 | 0.021559 | 0.0 | 0.13
Neigh | 4.544 | 4.544 | 4.544 | 0.0 | 27.27
Comm | 0.62792 | 0.62792 | 0.62792 | 0.0 | 3.77
Output | 0.0069666 | 0.0069666 | 0.0069666 | 0.0 | 0.04
Modify | 2.1282 | 2.1282 | 2.1282 | 0.0 | 12.77
Other | | 0.06897 | | | 0.41
Nlocal: 550.000 ave 550 max 550 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 2344.00 ave 2344 max 2344 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 36607.0 ave 36607 max 36607 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 36607
Ave neighs/atom = 66.558182
Ave special neighs/atom = 0.36363636
Neighbor list builds = 847
Dangerous builds = 0
Total wall time: 0:00:16

218
examples/deposit/log.10Mar21.deposit.molecule.rigid-nvt-small.g++.4 Normal file
View File

@ -0,0 +1,218 @@
LAMMPS (10 Mar 2021)
using 1 OpenMP thread(s) per MPI task
# sample surface deposition script for molecules
units lj
atom_style bond
boundary p p f
lattice fcc 1.0
Lattice spacing in x,y,z = 1.5874011 1.5874011 1.5874011
region box block 0 5 0 5 0 10
create_box 3 box bond/types 1 extra/bond/per/atom 1
Created orthogonal box = (0.0000000 0.0000000 0.0000000) to (7.9370053 7.9370053 15.874011)
1 by 1 by 4 MPI processor grid
region substrate block INF INF INF INF INF 3
create_atoms 1 region substrate
Created 350 atoms
create_atoms CPU = 0.002 seconds
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0
pair_coeff 1 2 1.0 1.0 5.0
mass * 1.0
bond_style harmonic
bond_coeff 1 5.0 1.0
neigh_modify delay 0
molecule dimer molecule.dimer
Read molecule template dimer:
1 molecules
2 atoms with max type 3
1 bonds with max type 1
0 angles with max type 0
0 dihedrals with max type 0
0 impropers with max type 0
region slab block 0 5 0 5 8 9
group addatoms empty
0 atoms in group addatoms
region mobile block 0 5 0 5 2 INF
group mobile region mobile
150 atoms in group mobile
compute add addatoms temp
compute_modify add dynamic/dof yes extra/dof 0
fix 1 addatoms rigid/nvt/small molecule temp 0.1 0.1 0.1 mol dimer
create bodies CPU = 0.000 seconds
0 rigid bodies with 0 atoms
1.0000000 = max distance from body owner to body atom
fix 2 mobile langevin 0.1 0.1 0.1 587283
fix 3 mobile nve
fix 4 addatoms deposit 100 0 100 12345 region slab near 1.0 mol dimer vz -1.0 -1.0 rigid 1
fix 5 addatoms wall/reflect zhi EDGE
thermo_style custom step atoms temp epair etotal press
thermo 100
thermo_modify temp add lost/bond ignore lost warn
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:468)
#dump 1 all atom 50 dump.deposit.atom
#dump 2 all image 50 image.*.jpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 2 pad 5
#dump 3 all movie 50 tmp.mpg type type # axes yes 0.8 0.02 view 80 -30
#dump_modify 3 pad 5
run 10000
WARNING: Should not allow rigid bodies to bounce off relecting walls (src/fix_wall_reflect.cpp:182)
Neighbor list info ...
update every 1 steps, delay 0 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 5.3
ghost atom cutoff = 5.3
binsize = 2.65, bins = 3 3 6
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair lj/cut, perpetual
attributes: half, newton on
pair build: half/bin/newton
stencil: half/bin/3d/newton
bin: standard
Per MPI rank memory allocation (min/avg/max) = 5.255 | 5.852 | 6.302 Mbytes
Step Atoms Temp E_pair TotEng Press
0 350 0 -6.9215833 -6.9215833 -1.0052629
100 352 1.0079368 -6.8946578 -6.8874992 -0.73775337
200 354 1.0081552 -6.8645575 -6.850318 -0.69629729
300 356 1.0085803 -6.821677 -6.8004288 -0.69532657
400 358 1.0099188 -6.7837923 -6.7555822 -0.68879568
500 360 1.0140221 -6.7446709 -6.7094618 -0.72991641
600 362 1.026146 -6.7129201 -6.6704003 -0.67063836
700 364 1.0683193 -6.6776523 -6.6262908 -0.65572472
800 366 1.0958894 -6.6402029 -6.5803182 -0.66307282
900 368 1.1231769 -6.6050912 -6.5364187 -0.64076928
1000 370 1.1976289 -6.5942508 -6.51333 -0.69249047
1100 372 1.05061 -6.5772306 -6.4995645 -0.6307289
1200 374 0.93723829 -6.5732962 -6.4981166 -0.64963821
1300 376 0.93899122 -6.5578418 -6.476679 -0.65096688
1400 378 0.87050694 -6.546852 -6.4662495 -0.67613401
1500 380 0.84462904 -6.5400883 -6.4567368 -0.64376986
1600 382 0.92417795 -6.5611292 -6.4643567 -0.62092779
1700 384 0.67296541 -6.5589985 -6.4845167 -0.6779111
1800 386 0.87159158 -6.5655058 -6.4638954 -0.6582786
1900 388 0.67504975 -6.5772537 -6.4946123 -0.665098
2000 390 0.65034978 -6.5717854 -6.4884072 -0.73162072
2100 392 0.63389859 -6.5652907 -6.4803935 -0.64937411
2200 394 0.6265637 -6.5582412 -6.4707767 -0.70819854
2300 396 0.85352983 -6.5750751 -6.4511408 -0.69995918
2400 398 0.69792348 -6.5938659 -6.4886513 -0.63755214
2500 400 0.45320495 -6.595022 -6.5242087 -0.70462183
2600 402 0.46824374 -6.5697613 -6.4940502 -0.75497893
2700 404 0.46605975 -6.5599272 -6.4820583 -0.66778506
2800 406 0.50637604 -6.5700492 -6.482743 -0.63812082
2900 408 0.55335921 -6.5873352 -6.4890054 -0.62891652
3000 410 0.51731142 -6.5890566 -6.4944264 -0.64305355
3100 412 0.59439734 -6.6003127 -6.4885025 -0.67140541
3200 414 0.53413632 -6.595564 -6.4923492 -0.76242248
3300 416 0.49221593 -6.5990831 -6.5014681 -0.71330819
3400 418 0.55006126 -6.6015725 -6.4897179 -0.72044309
3500 420 0.49065348 -6.6329864 -6.5307669 -0.65775397
3600 422 0.41907335 -6.6219753 -6.5325995 -0.75936391
3700 424 0.38720116 -6.6153349 -6.5308629 -0.74166397
3800 426 0.40625994 -6.6150209 -6.5244231 -0.7595111
3900 428 0.40460547 -6.6122642 -6.5200936 -0.70484465
4000 430 0.45014991 -6.6254404 -6.5207544 -0.70069933
4100 432 0.44820466 -6.640222 -6.5338771 -0.805972
4200 434 0.39767521 -6.6450316 -6.5488199 -0.72841414
4300 436 0.4155331 -6.6547917 -6.5523381 -0.6894484
4400 438 0.43034353 -6.6615074 -6.5534303 -0.74026062
4500 440 0.38062977 -6.6541618 -6.5568417 -0.85911194
4600 442 0.39357419 -6.6522517 -6.5498512 -0.66232114
4700 444 0.40296801 -6.6647029 -6.5580616 -0.67307577
4800 446 0.38194993 -6.6510248 -6.548258 -0.77887746
4900 448 0.40739835 -6.6601751 -6.5487771 -0.84146416
5000 450 0.38392807 -6.6560665 -6.5494198 -0.77343399
5100 452 0.35209286 -6.6481476 -6.5488294 -0.68065755
5200 454 0.41355143 -6.6615988 -6.543181 -0.81611092
5300 456 0.42200484 -6.6627494 -6.5401273 -0.7774134
5400 458 0.37764326 -6.661431 -6.550117 -0.72187808
5500 460 0.41857405 -6.6674232 -6.542306 -0.74929745
5600 462 0.41682611 -6.6710961 -6.5447852 -0.7557454
5700 464 0.44396148 -6.6924346 -6.5560887 -0.78018312
5800 466 0.37058077 -6.6865013 -6.5711919 -0.74146125
5900 468 0.346812 -6.681215 -6.57191 -0.81184233
6000 470 0.34919462 -6.6827931 -6.571348 -0.87330821
6100 472 0.39360936 -6.6933359 -6.5661634 -0.79237598
6200 474 0.33270778 -6.6847095 -6.5759127 -0.77978526
6300 476 0.35973804 -6.6951254 -6.5760944 -0.80340174
6400 478 0.38124318 -6.7045984 -6.5769857 -0.81628407
6500 480 0.41188302 -6.7133356 -6.573896 -0.7940289
6600 482 0.36998039 -6.7079555 -6.5813025 -0.75055442
6700 484 0.37615026 -6.722917 -6.592741 -0.76534055
6800 486 0.3466597 -6.7188712 -6.5976116 -0.77986211
6900 488 0.32902492 -6.7247054 -6.6084005 -0.75702458
7000 490 0.31856427 -6.7167709 -6.6029979 -0.76014555
7100 492 0.30233891 -6.7144406 -6.6053651 -0.96246708
7200 494 0.32557309 -6.7315347 -6.6129049 -0.8122153
7300 496 0.29919611 -6.7186327 -6.6085455 -0.71917931
7400 498 0.29778169 -6.7259455 -6.6153238 -0.88199391
7500 500 0.35535305 -6.7419073 -6.6086499 -0.90344083
7600 502 0.31979187 -6.7326802 -6.6116434 -0.78324572
7700 504 0.30821359 -6.72895 -6.6112295 -0.81363335
7800 506 0.31374993 -6.7255461 -6.6046346 -0.89904197
7900 508 0.29072338 -6.7177355 -6.6047082 -0.81423073
8000 510 0.30557494 -6.7208283 -6.600995 -0.85853575
8100 512 0.30521237 -6.726874 -6.6061601 -0.75361257
8200 514 0.29622226 -6.7152225 -6.5970794 -0.85766132
8300 516 0.28337698 -6.7023552 -6.5884003 -0.8985415
8400 518 0.32860902 -6.7211074 -6.5878875 -0.89758921
8500 520 0.35483743 -6.7406183 -6.5956126 -0.77602077
8600 522 0.32928486 -6.7326607 -6.5970358 -0.75070137
8700 524 0.36008106 -6.7474714 -6.5980103 -0.87093836
8800 526 0.36082301 -6.743579 -6.5926644 -0.90132107
8900 528 0.35010285 -6.7501553 -6.6026214 -0.85238959
9000 530 0.31513985 -6.7411795 -6.6073937 -0.75884529
9100 532 0.30895083 -6.7421063 -6.6099891 -0.74482692
9200 534 0.33631849 -6.751265 -6.6064087 -0.95632911
9300 536 0.33306096 -6.759187 -6.6147156 -0.76275935
9400 538 0.34582537 -6.7706766 -6.619619 -0.72251598
9500 540 0.32003146 -6.772057 -6.6313024 -0.88195851
9600 542 0.32439637 -6.7743569 -6.6307127 -0.90233104
9700 544 0.33988235 -6.7763721 -6.624862 -0.85185581
9800 546 0.32877587 -6.7744977 -6.62697 -0.8550905
9900 548 0.29570051 -6.7623752 -6.6288243 -1.0371157
10000 550 0.33675914 -6.7757315 -6.6226591 -1.0082157
Loop time of 16.1559 on 4 procs for 10000 steps with 550 atoms
Performance: 267395.139 tau/day, 618.970 timesteps/s
96.1% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.0035574 | 2.5467 | 7.5373 | 192.6 | 15.76
Bond | 0.0022929 | 0.0073138 | 0.019471 | 8.3 | 0.05
Neigh | 0.018645 | 1.2357 | 3.8167 | 138.4 | 7.65
Comm | 1.1946 | 6.6335 | 11.87 | 158.8 | 41.06
Output | 0.010714 | 0.014391 | 0.020006 | 3.2 | 0.09
Modify | 3.4474 | 5.6081 | 11.47 | 143.2 | 34.71
Other | | 0.1102 | | | 0.68
Nlocal: 137.500 ave 299 max 2 min
Histogram: 2 0 0 0 0 0 0 1 0 1
Nghost: 1903.75 ave 2686 max 529 min
Histogram: 1 0 0 0 0 0 1 0 0 2
Neighs: 9209.75 ave 23206 max 0 min
Histogram: 2 0 0 0 0 1 0 0 0 1
Total # of neighbors = 36839
Ave neighs/atom = 66.980000
Ave special neighs/atom = 0.36363636
Neighbor list builds = 829
Dangerous builds = 0
Total wall time: 0:00:16

6
lib/gpu/Makefile.hip
View File

@ -1,6 +1,6 @@
# /* ----------------------------------------------------------------------
# Generic Linux Makefile for HIP
# - export HIP_PLATFORM=hcc (or nvcc) before execution
# - export HIP_PLATFORM=amd (or nvcc) before execution
# - change HIP_ARCH for your GPU
# ------------------------------------------------------------------------- */
@ -42,6 +42,10 @@ ifeq (hcc,$(HIP_PLATFORM))
HIP_OPTS += -ffast-math
# possible values: gfx803,gfx900,gfx906
HIP_ARCH = gfx906
else ifeq (amd,$(HIP_PLATFORM))
HIP_OPTS += -ffast-math
# possible values: gfx803,gfx900,gfx906
HIP_ARCH = gfx906
else ifeq (nvcc,$(HIP_PLATFORM))
HIP_OPTS += --use_fast_math
HIP_ARCH = -gencode arch=compute_30,code=[sm_30,compute_30] -gencode arch=compute_32,code=[sm_32,compute_32] -gencode arch=compute_35,code=[sm_35,compute_35] \

5
lib/gpu/README
View File

@ -212,8 +212,8 @@ additionally requires cub (https://nvlabs.github.io/cub). Download and
extract the cub directory to lammps/lib/gpu/ or specify an appropriate
path in lammps/lib/gpu/Makefile.hip.
2. In Makefile.hip it is possible to specify the target platform via
export HIP_PLATFORM=hcc or HIP_PLATFORM=nvcc as well as the target
architecture (gfx803, gfx900, gfx906 etc.)
export HIP_PLATFORM=amd (ROCm >= 4.1), HIP_PLATFORM=hcc (ROCm <= 4.0)
or HIP_PLATFORM=nvcc as well as the target architecture (gfx803, gfx900, gfx906 etc.)
3. If your MPI implementation does not support `mpicxx --showme` command,
it is required to specify the corresponding MPI compiler and linker flags
in lammps/lib/gpu/Makefile.hip and in lammps/src/MAKE/OPTIONS/Makefile.hip.
@ -278,4 +278,3 @@ and
Brown, W.M., Masako, Y. Implementing Molecular Dynamics on Hybrid High
Performance Computers - Three-Body Potentials. Computer Physics Communications.
2013. 184: p. 2785–2793.

2
lib/gpu/lal_device.cpp
View File

@ -1061,7 +1061,7 @@ bool lmp_gpu_config(const std::string &category, const std::string &setting)
return setting == "opencl";
#elif defined(USE_HIP)
return setting == "hip";
#elif defined(USE_CUDA)
#elif defined(USE_CUDA) || defined(USE_CUDART)
return setting == "cuda";
#endif
return false;

8
lib/gpu/lal_lj_cubic.cu
View File

@ -26,10 +26,10 @@ _texture_2d( pos_tex,int4);
// LJ quantities scaled by epsilon and rmin = sigma*2^1/6 (see src/pair_lj_cubic.h)
#define _RT6TWO (numtyp)1.1224621
#define _PHIS (numtyp)-0.7869823 /* energy at s */
#define _DPHIDS (numtyp)2.6899009 /* gradient at s */
#define _A3 (numtyp)27.93357 /* cubic coefficient */
#define _RT6TWO (numtyp)1.1224620483093730 /* 2^1/6 */
#define _PHIS (numtyp)-0.7869822485207097 /* energy at s */
#define _DPHIDS (numtyp)2.6899008972047196 /* gradient at s */
#define _A3 (numtyp)27.9335700460986445 /* cubic coefficient */
__kernel void k_lj_cubic(const __global numtyp4 *restrict x_,
const __global numtyp4 *restrict lj1,

8
lib/gpu/lal_pre_cuda_hip.h
View File

@ -30,7 +30,7 @@
// -------------------------------------------------------------------------
#ifdef __HIP_PLATFORM_HCC__
#if defined(__HIP_PLATFORM_HCC__) || defined(__HIP_PLATFORM_AMD__)
#define CONFIG_ID 303
#define SIMD_SIZE 64
#else
@ -161,7 +161,7 @@
// KERNEL MACROS - TEXTURES
// -------------------------------------------------------------------------
#ifdef __HIP_PLATFORM_HCC__
#if defined(__HIP_PLATFORM_HCC__) || defined(__HIP_PLATFORM_AMD__)
#define _texture(name, type) __device__ type* name
#define _texture_2d(name, type) __device__ type* name
#else
@ -201,7 +201,7 @@
#define mu_tex mu_
#endif
#ifdef __HIP_PLATFORM_HCC__
#if defined(__HIP_PLATFORM_HCC__) || defined(__HIP_PLATFORM_AMD__)
#undef fetch4
#undef fetch
@ -266,7 +266,7 @@ typedef struct _double4 double4;
#endif
#endif
#if defined(CUDA_PRE_NINE) || defined(__HIP_PLATFORM_HCC__)
#if defined(CUDA_PRE_NINE) || defined(__HIP_PLATFORM_HCC__) || defined(__HIP_PLATFORM_AMD__)
#ifdef _SINGLE_SINGLE
#define shfl_down __shfl_down

64
python/lammps/core.py
View File

@ -1747,17 +1747,23 @@ class lammps(object):
# -------------------------------------------------------------------------
def find_pair_neighlist(self, style, exact=True, nsub=0, request=0):
def find_pair_neighlist(self, style, exact=True, nsub=0, reqid=0):
"""Find neighbor list index of pair style neighbor list
Try finding pair instance that matches style. If exact is set, the pair must
match style exactly. If exact is 0, style must only be contained. If pair is
of style pair/hybrid, style is instead matched the nsub-th hybrid sub-style.
Search for a neighbor list requested by a pair style instance that
matches "style". If exact is True, the pair style name must match
exactly. If exact is False, the pair style name is matched against
"style" as regular expression or sub-string. If the pair style is a
hybrid pair style, the style is instead matched against the hybrid
sub-styles. If the same pair style is used as sub-style multiple
types, you must set nsub to a value n > 0 which indicates the nth
instance of that sub-style to be used (same as for the pair_coeff
command). The default value of 0 will fail to match in that case.
Once the pair instance has been identified, multiple neighbor list requests
may be found. Every neighbor list is uniquely identified by its request
index. Thus, providing this request index ensures that the correct neighbor
list index is returned.
Once the pair style instance has been identified, it may have
requested multiple neighbor lists. Those are uniquely identified by
a request ID > 0 as set by the pair style. Otherwise the request
ID is 0.
:param style: name of pair style that should be searched for
:type style: string
@ -1765,44 +1771,58 @@ class lammps(object):
:type exact: bool, optional
:param nsub: match nsub-th hybrid sub-style, defaults to 0
:type nsub: int, optional
:param request: index of neighbor list request, in case there are more than one, defaults to 0
:type request: int, optional
:param reqid: list request id, > 0 in case there are more than one, defaults to 0
:type reqid: int, optional
:return: neighbor list index if found, otherwise -1
:rtype: int
"""
"""
style = style.encode()
exact = int(exact)
idx = self.lib.lammps_find_pair_neighlist(self.lmp, style, exact, nsub, request)
idx = self.lib.lammps_find_pair_neighlist(self.lmp, style, exact, nsub, reqid)
return idx
# -------------------------------------------------------------------------
def find_fix_neighlist(self, fixid, request=0):
def find_fix_neighlist(self, fixid, reqid=0):
"""Find neighbor list index of fix neighbor list
The fix instance requesting the neighbor list is uniquely identified
by the fix ID. In case the fix has requested multiple neighbor
lists, those are uniquely identified by a request ID > 0 as set by
the fix. Otherwise the request ID is 0 (the default).
:param fixid: name of fix
:type fixid: string
:param request: index of neighbor list request, in case there are more than one, defaults to 0
:type request: int, optional
:param reqid: id of neighbor list request, in case there are more than one request, defaults to 0
:type reqid: int, optional
:return: neighbor list index if found, otherwise -1
:rtype: int
"""
"""
fixid = fixid.encode()
idx = self.lib.lammps_find_fix_neighlist(self.lmp, fixid, request)
idx = self.lib.lammps_find_fix_neighlist(self.lmp, fixid, reqid)
return idx
# -------------------------------------------------------------------------
def find_compute_neighlist(self, computeid, request=0):
def find_compute_neighlist(self, computeid, reqid=0):
"""Find neighbor list index of compute neighbor list
The compute instance requesting the neighbor list is uniquely
identified by the compute ID. In case the compute has requested
multiple neighbor lists, those are uniquely identified by a request
ID > 0 as set by the compute. Otherwise the request ID is 0 (the
default).
:param computeid: name of compute
:type computeid: string
:param request: index of neighbor list request, in case there are more than one, defaults to 0
:type request: int, optional
:param reqid: index of neighbor list request, in case there are more than one request, defaults to 0
:type reqid: int, optional
:return: neighbor list index if found, otherwise -1
:rtype: int
"""
"""
computeid = computeid.encode()
idx = self.lib.lammps_find_compute_neighlist(self.lmp, computeid, request)
idx = self.lib.lammps_find_compute_neighlist(self.lmp, computeid, reqid)
return idx

21
python/lammps/data.py
View File

@ -52,7 +52,9 @@ class NeighList:
def get(self, element):
"""
:return: tuple with atom local index, numpy array of neighbor local atom indices
Access a specific neighbor list entry. "element" must be a number from 0 to the size-1 of the list
:return: tuple with atom local index, number of neighbors and ctypes pointer to neighbor's local atom indices
:rtype: (int, int, ctypes.POINTER(c_int))
"""
iatom, numneigh, neighbors = self.lmp.get_neighlist_element_neighbors(self.idx, element)
@ -71,3 +73,20 @@ class NeighList:
for ii in range(inum):
yield self.get(ii)
def find(self, iatom):
"""
Find the neighbor list for a specific (local) atom iatom.
If there is no list for iatom, (-1, None) is returned.
:return: tuple with number of neighbors and ctypes pointer to neighbor's local atom indices
:rtype: (int, ctypes.POINTER(c_int))
"""
inum = self.size
for ii in range(inum):
idx, numneigh, neighbors = self.get(ii)
if idx == iatom:
return numneigh, neighbors
return -1, None

17
python/lammps/numpy_wrapper.py
View File

@ -331,8 +331,25 @@ class NumPyNeighList(NeighList):
def get(self, element):
"""
Access a specific neighbor list entry. "element" must be a number from 0 to the size-1 of the list
:return: tuple with atom local index, numpy array of neighbor local atom indices
:rtype: (int, numpy.array)
"""
iatom, neighbors = self.lmp.numpy.get_neighlist_element_neighbors(self.idx, element)
return iatom, neighbors
def find(self, iatom):
"""
Find the neighbor list for a specific (local) atom iatom.
If there is no list for iatom, None is returned.
:return: numpy array of neighbor local atom indices
:rtype: numpy.array or None
"""
inum = self.size
for ii in range(inum):
idx, neighbors = self.get(ii)
if idx == iatom:
return neighbors
return None

17
src/CLASS2/pair_lj_class2.cpp
View File

@ -11,20 +11,20 @@
#include "pair_lj_class2.h"
#include <cmath>
#include <cstring>
#include "atom.h"
#include "comm.h"
#include "error.h"
#include "force.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "update.h"
#include "respa.h"
#include "math_const.h"
#include "memory.h"
#include "error.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "neighbor.h"
#include "respa.h"
#include "update.h"
#include <cmath>
#include <cstring>
using namespace LAMMPS_NS;
using namespace MathConst;
@ -36,6 +36,7 @@ PairLJClass2::PairLJClass2(LAMMPS *lmp) : Pair(lmp)
respa_enable = 1;
writedata = 1;
centroidstressflag = CENTROID_SAME;
cut_respa = nullptr;
}
/* ---------------------------------------------------------------------- */

32
src/GPU/pair_lj_class2_gpu.cpp
View File

@ -16,25 +16,24 @@
------------------------------------------------------------------------- */
#include "pair_lj_class2_gpu.h"
#include <cmath>
#include <cstdio>
#include <cstring>
#include "atom.h"
#include "atom_vec.h"
#include "comm.h"
#include "domain.h"
#include "error.h"
#include "force.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "gpu_extra.h"
#include "integrate.h"
#include "memory.h"
#include "error.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "neighbor.h"
#include "suffix.h"
#include "universe.h"
#include "update.h"
#include "domain.h"
#include "gpu_extra.h"
#include "suffix.h"
#include <cmath>
using namespace LAMMPS_NS;
@ -46,13 +45,13 @@ int lj96_gpu_init(const int ntypes, double **cutsq, double **host_lj1,
const int nall, const int max_nbors, const int maxspecial,
const double cell_size, int &gpu_mode, FILE *screen);
void lj96_gpu_clear();
int ** lj96_gpu_compute_n(const int ago, const int inum, const int nall,
double **host_x, int *host_type, double *sublo,
double *subhi, tagint *tag, int **nspecial,
tagint **special, const bool eflag, const bool vflag,
const bool eatom, const bool vatom, int &host_start,
int **ilist, int **jnum,
const double cpu_time, bool &success);
int **lj96_gpu_compute_n(const int ago, const int inum, const int nall,
double **host_x, int *host_type, double *sublo,
double *subhi, tagint *tag, int **nspecial,
tagint **special, const bool eflag, const bool vflag,
const bool eatom, const bool vatom, int &host_start,
int **ilist, int **jnum,
const double cpu_time, bool &success);
void lj96_gpu_compute(const int ago, const int inum, const int nall,
double **host_x, int *host_type, int *ilist, int *numj,
int **firstneigh, const bool eflag, const bool vflag,
@ -64,6 +63,7 @@ double lj96_gpu_bytes();
PairLJClass2GPU::PairLJClass2GPU(LAMMPS *lmp) : PairLJClass2(lmp), gpu_mode(GPU_FORCE)
{
respa_enable = 0;
reinitflag = 0;
cpu_time = 0.0;
suffix_flag |= Suffix::GPU;

19
src/GPU/pair_lj_cubic_gpu.cpp
View File

@ -17,25 +17,24 @@
#include "pair_lj_cubic_gpu.h"
#include <cmath>
#include <cstdio>
#include <cstring>
#include "atom.h"
#include "atom_vec.h"
#include "comm.h"
#include "domain.h"
#include "error.h"
#include "force.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "gpu_extra.h"
#include "integrate.h"
#include "memory.h"
#include "error.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "neighbor.h"
#include "suffix.h"
#include "universe.h"
#include "update.h"
#include "domain.h"
#include "gpu_extra.h"
#include "suffix.h"
#include <cmath>
#include <cstring>
#include "pair_lj_cubic_const.h"

238
src/GRANULAR/fix_wall_gran.cpp
View File

@ -54,7 +54,8 @@ enum {XPLANE=0,YPLANE=1,ZPLANE=2,ZCYLINDER,REGION};
enum {NONE,CONSTANT,EQUAL};
enum {DAMPING_NONE, VELOCITY, MASS_VELOCITY, VISCOELASTIC, TSUJI};
enum {TANGENTIAL_NONE, TANGENTIAL_NOHISTORY, TANGENTIAL_HISTORY,
TANGENTIAL_MINDLIN, TANGENTIAL_MINDLIN_RESCALE};
TANGENTIAL_MINDLIN, TANGENTIAL_MINDLIN_RESCALE,
TANGENTIAL_MINDLIN_FORCE, TANGENTIAL_MINDLIN_RESCALE_FORCE};
enum {TWIST_NONE, TWIST_SDS, TWIST_MARSHALL};
enum {ROLL_NONE, ROLL_SDS};
@ -70,6 +71,7 @@ FixWallGran::FixWallGran(LAMMPS *lmp, int narg, char **arg) :
error->all(FLERR,"Fix wall/gran requires atom style sphere");
create_attribute = 1;
limit_damping = 0;
// set interaction style
// disable bonded/history option for now
@ -90,7 +92,6 @@ FixWallGran::FixWallGran(LAMMPS *lmp, int narg, char **arg) :
// wall/particle coefficients
int iarg;
if (pairstyle != GRANULAR) {
size_history = 3;
if (narg < 11) error->all(FLERR,"Illegal fix wall/gran command");
@ -118,6 +119,12 @@ FixWallGran::FixWallGran(LAMMPS *lmp, int narg, char **arg) :
kt /= force->nktv2p;
}
iarg = 10;
if (strcmp(arg[iarg],"limit_damping") == 0) {
limit_damping = 1;
iarg += 1;
}
} else {
iarg = 4;
damping_model = VISCOELASTIC;
@ -211,7 +218,9 @@ FixWallGran::FixWallGran(LAMMPS *lmp, int narg, char **arg) :
iarg += 4;
} else if ((strcmp(arg[iarg+1], "linear_history") == 0) ||
(strcmp(arg[iarg+1], "mindlin") == 0) ||
(strcmp(arg[iarg+1], "mindlin_rescale") == 0)) {
(strcmp(arg[iarg+1], "mindlin_rescale") == 0) ||
(strcmp(arg[iarg+1], "mindlin/force") == 0) ||
(strcmp(arg[iarg+1], "mindlin_rescale/force") == 0)) {
if (iarg + 4 >= narg)
error->all(FLERR,"Illegal pair_coeff command, "
"not enough parameters provided for tangential model");
@ -221,8 +230,14 @@ FixWallGran::FixWallGran(LAMMPS *lmp, int narg, char **arg) :
tangential_model = TANGENTIAL_MINDLIN;
else if (strcmp(arg[iarg+1], "mindlin_rescale") == 0)
tangential_model = TANGENTIAL_MINDLIN_RESCALE;
else if (strcmp(arg[iarg+1], "mindlin/force") == 0)
tangential_model = TANGENTIAL_MINDLIN_FORCE;
else if (strcmp(arg[iarg+1], "mindlin_rescale/force") == 0)
tangential_model = TANGENTIAL_MINDLIN_RESCALE_FORCE;
if ((tangential_model == TANGENTIAL_MINDLIN ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE) &&
tangential_model == TANGENTIAL_MINDLIN_RESCALE ||
tangential_model == TANGENTIAL_MINDLIN_FORCE ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE_FORCE) &&
(strcmp(arg[iarg+2], "NULL") == 0)) {
if (normal_model == NORMAL_HERTZ || normal_model == NORMAL_HOOKE) {
error->all(FLERR, "NULL setting for Mindlin tangential "
@ -293,13 +308,27 @@ FixWallGran::FixWallGran(LAMMPS *lmp, int narg, char **arg) :
strcmp(arg[iarg], "zcylinder") == 0 ||
strcmp(arg[iarg], "region") == 0) {
break;
} else if (strcmp(arg[iarg],"limit_damping") == 0) {
limit_damping = 1;
iarg += 1;
} else {
error->all(FLERR, "Illegal fix wall/gran command");
}
}
size_history = 3*tangential_history + 3*roll_history + twist_history;
//Unlike the pair style, the wall style does not have a 'touch'
//array. Hence, an additional entry in the history is used to
//determine if particles previously contacted for JKR cohesion purposes.
if (normal_model == JKR) size_history += 1;
if (tangential_model == TANGENTIAL_MINDLIN_RESCALE) size_history += 1;
if (tangential_model == TANGENTIAL_MINDLIN_RESCALE ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE_FORCE) size_history += 1;
if (limit_damping && normal_model == JKR)
error->all(FLERR,"Illegal pair_coeff command, "
"cannot limit damping with JRK model");
if (limit_damping && normal_model == DMT)
error->all(FLERR,"Illegal pair_coeff command, "
"Cannot limit damping with DMT model");
}
// wallstyle args
@ -475,36 +504,39 @@ void FixWallGran::init()
if (modify->fix[i]->rigid_flag) break;
if (i < modify->nfix) fix_rigid = modify->fix[i];
tangential_history_index = 0;
if (roll_history) {
if (tangential_history) roll_history_index = 3;
else roll_history_index = 0;
}
if (twist_history) {
if (tangential_history) {
if (roll_history) twist_history_index = 6;
else twist_history_index = 3;
if(pairstyle == GRANULAR) {
tangential_history_index = 0;
if (roll_history) {
if (tangential_history) roll_history_index = 3;
else roll_history_index = 0;
}
else{
if (roll_history) twist_history_index = 3;
else twist_history_index = 0;
if (twist_history) {
if (tangential_history) {
if (roll_history) twist_history_index = 6;
else twist_history_index = 3;
}
else{
if (roll_history) twist_history_index = 3;
else twist_history_index = 0;
}
}
if (normal_model == JKR) {
tangential_history_index += 1;
roll_history_index += 1;
twist_history_index += 1;
}
if (tangential_model == TANGENTIAL_MINDLIN_RESCALE ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE_FORCE) {
roll_history_index += 1;
twist_history_index += 1;
}
}
if (normal_model == JKR) {
tangential_history_index += 1;
roll_history_index += 1;
twist_history_index += 1;
}
if (tangential_model == TANGENTIAL_MINDLIN_RESCALE) {
roll_history_index += 1;
twist_history_index += 1;
}
if (damping_model == TSUJI) {
double cor = normal_coeffs[1];
normal_coeffs[1] = 1.2728-4.2783*cor+11.087*pow(cor,2)-22.348*pow(cor,3)+
27.467*pow(cor,4)-18.022*pow(cor,5)+
4.8218*pow(cor,6);
if (damping_model == TSUJI) {
double cor = normal_coeffs[1];
normal_coeffs[1] = 1.2728-4.2783*cor+11.087*pow(cor,2)-22.348*pow(cor,3)+
27.467*pow(cor,4)-18.022*pow(cor,5)+
4.8218*pow(cor,6);
}
}
}
@ -762,6 +794,7 @@ void FixWallGran::hooke(double rsq, double dx, double dy, double dz,
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radius-r)*rinv - damp;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities
@ -854,6 +887,7 @@ void FixWallGran::hooke_history(double rsq, double dx, double dy, double dz,
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radius-r)*rinv - damp;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities
@ -985,6 +1019,7 @@ void FixWallGran::hertz_history(double rsq, double dx, double dy, double dz,
if (rwall == 0.0) polyhertz = sqrt((radius-r)*radius);
else polyhertz = sqrt((radius-r)*radius*rwall/(rwall+radius));
ccel *= polyhertz;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities
@ -1103,7 +1138,8 @@ void FixWallGran::granular(double rsq, double dx, double dy, double dz,
double signtwist, magtwist, magtortwist, Mtcrit;
double tortwist1, tortwist2, tortwist3;
double shrmag,rsht;
double shrmag,rsht,prjmag;
bool frameupdate;
r = sqrt(rsq);
E = normal_coeffs[0];
@ -1178,11 +1214,18 @@ void FixWallGran::granular(double rsq, double dx, double dy, double dz,
Fdamp = -damp_normal_prefactor*vnnr;
Fntot = Fne + Fdamp;
if (limit_damping && (Fntot < 0.0)) Fntot = 0.0;
//****************************************
// tangential force, including history effects
//****************************************
// For linear, mindlin, mindlin_rescale:
// history = cumulative tangential displacement
//
// For mindlin/force, mindlin_rescale/force:
// history = cumulative tangential elastic force
// tangential component
vt1 = vr1 - vn1;
vt2 = vr2 - vn2;
@ -1224,69 +1267,115 @@ void FixWallGran::granular(double rsq, double dx, double dy, double dz,
int thist2 = thist1 + 1;
if (tangential_history) {
if (tangential_model == TANGENTIAL_MINDLIN) {
if (tangential_model == TANGENTIAL_MINDLIN ||
tangential_model == TANGENTIAL_MINDLIN_FORCE) {
k_tangential *= a;
}
else if (tangential_model == TANGENTIAL_MINDLIN_RESCALE) {
else if (tangential_model ==
TANGENTIAL_MINDLIN_RESCALE ||
tangential_model ==
TANGENTIAL_MINDLIN_RESCALE_FORCE){
k_tangential *= a;
if (a < history[3]) { //On unloading, rescale the shear displacements
// on unloading, rescale the shear displacements/force
if (a < history[thist2+1]) {
double factor = a/history[thist2+1];
history[thist0] *= factor;
history[thist1] *= factor;
history[thist2] *= factor;
}
}
shrmag = sqrt(history[thist0]*history[thist0] +
history[thist1]*history[thist1] +
history[thist2]*history[thist2]);
// rotate and update displacements.
// see e.g. eq. 17 of Luding, Gran. Matter 2008, v10,p235
if (history_update) {
rsht = history[thist0]*nx + history[thist1]*ny + history[thist2]*nz;
if (fabs(rsht) < EPSILON) rsht = 0;
if (rsht > 0) {
// if rhst == shrmag, contacting pair has rotated 90 deg in one step,
// in which case you deserve a crash!
scalefac = shrmag/(shrmag - rsht);
if (tangential_model == TANGENTIAL_MINDLIN_FORCE ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE_FORCE)
frameupdate = fabs(rsht) > EPSILON*Fscrit;
else
frameupdate = fabs(rsht)*k_tangential > EPSILON*Fscrit;
if (frameupdate) {
shrmag = sqrt(history[thist0]*history[thist0] +
history[thist1]*history[thist1] +
history[thist2]*history[thist2]);
// projection
history[thist0] -= rsht*nx;
history[thist1] -= rsht*ny;
history[thist2] -= rsht*nz;
// also rescale to preserve magnitude
prjmag = sqrt(history[thist0]*history[thist0] +
history[thist1]*history[thist1] + history[thist2]*history[thist2]);
if (prjmag > 0) scalefac = shrmag/prjmag;
else scalefac = 0;
history[thist0] *= scalefac;
history[thist1] *= scalefac;
history[thist2] *= scalefac;
}
// update history
if (tangential_model == TANGENTIAL_HISTORY ||
tangential_model == TANGENTIAL_MINDLIN ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE) {
history[thist0] += vtr1*dt;
history[thist1] += vtr2*dt;
history[thist2] += vtr3*dt;
} else{
// tangential force
// see e.g. eq. 18 of Thornton et al, Pow. Tech. 2013, v223,p30-46
history[thist0] -= k_tangential*vtr1*dt;
history[thist1] -= k_tangential*vtr2*dt;
history[thist2] -= k_tangential*vtr3*dt;
}
if (tangential_model == TANGENTIAL_MINDLIN_RESCALE ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE_FORCE)
history[thist2+1] = a;
}
// tangential forces = history + tangential velocity damping
fs1 = -k_tangential*history[thist0] - damp_tangential*vtr1;
fs2 = -k_tangential*history[thist1] - damp_tangential*vtr2;
fs3 = -k_tangential*history[thist2] - damp_tangential*vtr3;
if (tangential_model == TANGENTIAL_HISTORY ||
tangential_model == TANGENTIAL_MINDLIN ||
tangential_model == TANGENTIAL_MINDLIN_RESCALE) {
fs1 = -k_tangential*history[thist0] - damp_tangential*vtr1;
fs2 = -k_tangential*history[thist1] - damp_tangential*vtr2;
fs3 = -k_tangential*history[thist2] - damp_tangential*vtr3;
} else {
fs1 = history[thist0] - damp_tangential*vtr1;
fs2 = history[thist1] - damp_tangential*vtr2;
fs3 = history[thist2] - damp_tangential*vtr3;
}
// rescale frictional displacements and forces if needed
Fscrit = tangential_coeffs[2] * Fncrit;
fs = sqrt(fs1*fs1 + fs2*fs2 + fs3*fs3);
if (fs > Fscrit) {
shrmag = sqrt(history[thist0]*history[thist0] +
history[thist1]*history[thist1] +
history[thist2]*history[thist2]);
if (shrmag != 0.0) {
history[thist0] = -1.0/k_tangential*(Fscrit*fs1/fs +
damp_tangential*vtr1);
history[thist1] = -1.0/k_tangential*(Fscrit*fs2/fs +
damp_tangential*vtr2);
history[thist2] = -1.0/k_tangential*(Fscrit*fs3/fs +
damp_tangential*vtr3);
if (tangential_model == TANGENTIAL_HISTORY ||
tangential_model == TANGENTIAL_MINDLIN ||
tangential_model ==
TANGENTIAL_MINDLIN_RESCALE) {
history[thist0] = -1.0/k_tangential*(Fscrit*fs1/fs +
damp_tangential*vtr1);
history[thist1] = -1.0/k_tangential*(Fscrit*fs2/fs +
damp_tangential*vtr2);
history[thist2] = -1.0/k_tangential*(Fscrit*fs3/fs +
damp_tangential*vtr3);
} else {
history[thist0] = Fscrit*fs1/fs + damp_tangential*vtr1;
history[thist1] = Fscrit*fs2/fs + damp_tangential*vtr2;
history[thist2] = Fscrit*fs3/fs + damp_tangential*vtr3;
}
fs1 *= Fscrit/fs;
fs2 *= Fscrit/fs;
fs3 *= Fscrit/fs;
} else fs1 = fs2 = fs3 = 0.0;
}
} else { // classic pair gran/hooke (no history)
fs = meff*damp_tangential*vrel;
if (vrel != 0.0) Ft = MIN(Fne,fs) / vrel;
fs = damp_tangential*vrel;
if (vrel != 0.0) Ft = MIN(Fscrit,fs) / vrel;
else Ft = 0.0;
fs1 = -Ft*vtr1;
fs2 = -Ft*vtr2;
@ -1297,14 +1386,14 @@ void FixWallGran::granular(double rsq, double dx, double dy, double dz,
// rolling resistance
//****************************************
if (roll_model != ROLL_NONE || twist_model != NONE) {
if (roll_model != ROLL_NONE || twist_model != TWIST_NONE) {
relrot1 = omega[0];
relrot2 = omega[1];
relrot3 = omega[2];
}
if (roll_model != ROLL_NONE) {
// rolling velocity, see eq. 31 of Wang et al, Particuology v 23, p 49 (2015)
// rolling velocity,
// see eq. 31 of Wang et al, Particuology v 23, p 49 (2015)
// This is different from the Marshall papers,
// which use the Bagi/Kuhn formulation
// for rolling velocity (see Wang et al for why the latter is wrong)
@ -1316,21 +1405,29 @@ void FixWallGran::granular(double rsq, double dx, double dy, double dz,
int rhist1 = rhist0 + 1;
int rhist2 = rhist1 + 1;
// rolling displacement
rollmag = sqrt(history[rhist0]*history[rhist0] +
history[rhist1]*history[rhist1] +
history[rhist2]*history[rhist2]);
rolldotn = history[rhist0]*nx + history[rhist1]*ny + history[rhist2]*nz;
k_roll = roll_coeffs[0];
damp_roll = roll_coeffs[1];
Frcrit = roll_coeffs[2] * Fncrit;
if (history_update) {
if (fabs(rolldotn) < EPSILON) rolldotn = 0;
if (rolldotn > 0) { // rotate into tangential plane
scalefac = rollmag/(rollmag - rolldotn);
rolldotn = history[rhist0]*nx + history[rhist1]*ny + history[rhist2]*nz;
frameupdate = fabs(rolldotn)*k_roll > EPSILON*Frcrit;
if (frameupdate) { // rotate into tangential plane
rollmag = sqrt(history[rhist0]*history[rhist0] +
history[rhist1]*history[rhist1] +
history[rhist2]*history[rhist2]);
// projection
history[rhist0] -= rolldotn*nx;
history[rhist1] -= rolldotn*ny;
history[rhist2] -= rolldotn*nz;
// also rescale to preserve magnitude
prjmag = sqrt(history[rhist0]*history[rhist0] +
history[rhist1]*history[rhist1] +
history[rhist2]*history[rhist2]);
if (prjmag > 0) scalefac = rollmag/prjmag;
else scalefac = 0;
history[rhist0] *= scalefac;
history[rhist1] *= scalefac;
history[rhist2] *= scalefac;
@ -1340,17 +1437,16 @@ void FixWallGran::granular(double rsq, double dx, double dy, double dz,
history[rhist2] += vrl3*dt;
}
k_roll = roll_coeffs[0];
damp_roll = roll_coeffs[1];
fr1 = -k_roll*history[rhist0] - damp_roll*vrl1;
fr2 = -k_roll*history[rhist1] - damp_roll*vrl2;
fr3 = -k_roll*history[rhist2] - damp_roll*vrl3;
// rescale frictional displacements and forces if needed
Frcrit = roll_coeffs[2] * Fncrit;
fr = sqrt(fr1*fr1 + fr2*fr2 + fr3*fr3);
if (fr > Frcrit) {
rollmag = sqrt(history[rhist0]*history[rhist0] +
history[rhist1]*history[rhist1] +
history[rhist2]*history[rhist2]);
if (rollmag != 0.0) {
history[rhist0] = -1.0/k_roll*(Frcrit*fr1/fr + damp_roll*vrl1);
history[rhist1] = -1.0/k_roll*(Frcrit*fr2/fr + damp_roll*vrl2);

1
src/GRANULAR/fix_wall_gran.h
View File

@ -73,6 +73,7 @@ class FixWallGran : public Fix {
// for granular model choices
int normal_model, damping_model;
int tangential_model, roll_model, twist_model;
int limit_damping;
// history flags
int normal_history, tangential_history, roll_history, twist_history;

10
src/GRANULAR/pair_gran_hertz_history.cpp
View File

@ -181,6 +181,7 @@ void PairGranHertzHistory::compute(int eflag, int vflag)
ccel = kn*(radsum-r)*rinv - damp;
polyhertz = sqrt((radsum-r)*radi*radj / radsum);
ccel *= polyhertz;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities
@ -277,7 +278,7 @@ void PairGranHertzHistory::compute(int eflag, int vflag)
void PairGranHertzHistory::settings(int narg, char **arg)
{
if (narg != 6) error->all(FLERR,"Illegal pair_style command");
if (narg != 6 && narg != 7) error->all(FLERR,"Illegal pair_style command");
kn = utils::numeric(FLERR,arg[0],false,lmp);
if (strcmp(arg[1],"NULL") == 0) kt = kn * 2.0/7.0;
@ -295,6 +296,12 @@ void PairGranHertzHistory::settings(int narg, char **arg)
xmu < 0.0 || xmu > 10000.0 || dampflag < 0 || dampflag > 1)
error->all(FLERR,"Illegal pair_style command");
limit_damping = 0;
if (narg == 7) {
if (strcmp(arg[6], "limit_damping") == 0) limit_damping = 1;
else error->all(FLERR,"Illegal pair_style command");
}
// convert Kn and Kt from pressure units to force/distance^2
kn /= force->nktv2p;
@ -388,6 +395,7 @@ double PairGranHertzHistory::single(int i, int j, int /*itype*/, int /*jtype*/,
ccel = kn*(radsum-r)*rinv - damp;
polyhertz = sqrt((radsum-r)*radi*radj / radsum);
ccel *= polyhertz;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities

2
src/GRANULAR/pair_gran_hooke.cpp
View File

@ -158,6 +158,7 @@ void PairGranHooke::compute(int eflag, int vflag)
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radsum-r)*rinv - damp;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities
@ -299,6 +300,7 @@ double PairGranHooke::single(int i, int j, int /*itype*/, int /*jtype*/, double
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radsum-r)*rinv - damp;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities

10
src/GRANULAR/pair_gran_hooke_history.cpp
View File

@ -237,6 +237,7 @@ void PairGranHookeHistory::compute(int eflag, int vflag)
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radsum-r)*rinv - damp;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities
@ -356,7 +357,7 @@ void PairGranHookeHistory::allocate()
void PairGranHookeHistory::settings(int narg, char **arg)
{
if (narg != 6) error->all(FLERR,"Illegal pair_style command");
if (narg != 6 && narg != 7) error->all(FLERR,"Illegal pair_style command");
kn = utils::numeric(FLERR,arg[0],false,lmp);
if (strcmp(arg[1],"NULL") == 0) kt = kn * 2.0/7.0;
@ -370,6 +371,12 @@ void PairGranHookeHistory::settings(int narg, char **arg)
dampflag = utils::inumeric(FLERR,arg[5],false,lmp);
if (dampflag == 0) gammat = 0.0;
limit_damping = 0;
if (narg == 7) {
if (strcmp(arg[6], "limit_damping") == 0) limit_damping = 1;
else error->all(FLERR,"Illegal pair_style command");
}
if (kn < 0.0 || kt < 0.0 || gamman < 0.0 || gammat < 0.0 ||
xmu < 0.0 || xmu > 10000.0 || dampflag < 0 || dampflag > 1)
error->all(FLERR,"Illegal pair_style command");
@ -686,6 +693,7 @@ double PairGranHookeHistory::single(int i, int j, int /*itype*/, int /*jtype*/,
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radsum-r)*rinv - damp;
if(limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities

1
src/GRANULAR/pair_gran_hooke_history.h
View File

@ -49,6 +49,7 @@ class PairGranHookeHistory : public Pair {
double dt;
int freeze_group_bit;
int history;
int limit_damping;
int neighprev;
double *onerad_dynamic,*onerad_frozen;

43
src/GRANULAR/pair_granular.cpp
View File

@ -81,6 +81,8 @@ PairGranular::PairGranular(LAMMPS *lmp) : Pair(lmp)
maxrad_dynamic = nullptr;
maxrad_frozen = nullptr;
limit_damping = nullptr;
history_transfer_factors = nullptr;
dt = update->dt;
@ -130,6 +132,7 @@ PairGranular::~PairGranular()
memory->destroy(tangential_model);
memory->destroy(roll_model);
memory->destroy(twist_model);
memory->destroy(limit_damping);
delete [] onerad_dynamic;
delete [] onerad_frozen;
@ -364,12 +367,13 @@ void PairGranular::compute(int eflag, int vflag)
damp_normal = a*meff;
} else if (damping_model[itype][jtype] == TSUJI) {
damp_normal = sqrt(meff*knfac);
}
} else damp_normal = 0.0;
damp_normal_prefactor = normal_coeffs[itype][jtype][1]*damp_normal;
Fdamp = -damp_normal_prefactor*vnnr;
Fntot = Fne + Fdamp;
if (limit_damping[itype][jtype] && (Fntot < 0.0)) Fntot = 0.0;
//****************************************
// tangential force, including history effects
@ -540,20 +544,17 @@ void PairGranular::compute(int eflag, int vflag)
relrot1 = omega[i][0] - omega[j][0];
relrot2 = omega[i][1] - omega[j][1];
relrot3 = omega[i][2] - omega[j][2];
// rolling velocity,
// see eq. 31 of Wang et al, Particuology v 23, p 49 (2015)
// this is different from the Marshall papers,
// which use the Bagi/Kuhn formulation
// for rolling velocity (see Wang et al for why the latter is wrong)
// - 0.5*((radj-radi)/radsum)*vtr1;
// - 0.5*((radj-radi)/radsum)*vtr2;
// - 0.5*((radj-radi)/radsum)*vtr3;
}
//****************************************
// rolling resistance
//****************************************
if (roll_model[itype][jtype] != ROLL_NONE) {
// rolling velocity,
// see eq. 31 of Wang et al, Particuology v 23, p 49 (2015)
// this is different from the Marshall papers,
// which use the Bagi/Kuhn formulation
// for rolling velocity (see Wang et al for why the latter is wrong)
vrl1 = Reff*(relrot2*nz - relrot3*ny);
vrl2 = Reff*(relrot3*nx - relrot1*nz);
vrl3 = Reff*(relrot1*ny - relrot2*nx);
@ -740,6 +741,7 @@ void PairGranular::allocate()
memory->create(tangential_model,n+1,n+1,"pair:tangential_model");
memory->create(roll_model,n+1,n+1,"pair:roll_model");
memory->create(twist_model,n+1,n+1,"pair:twist_model");
memory->create(limit_damping,n+1,n+1,"pair:limit_damping");
onerad_dynamic = new double[n+1];
onerad_frozen = new double[n+1];
@ -793,6 +795,7 @@ void PairGranular::coeff(int narg, char **arg)
roll_model_one = ROLL_NONE;
twist_model_one = TWIST_NONE;
damping_model_one = VISCOELASTIC;
int ld_flag = 0;
int iarg = 2;
while (iarg < narg) {
@ -967,12 +970,15 @@ void PairGranular::coeff(int narg, char **arg)
error->all(FLERR, "Illegal pair_coeff command, not enough parameters");
cutoff_one = utils::numeric(FLERR,arg[iarg+1],false,lmp);
iarg += 2;
} else error->all(FLERR, "Illegal pair coeff command");
} else if (strcmp(arg[iarg], "limit_damping") == 0) {
ld_flag = 1;
iarg += 1;
} else error->all(FLERR, "Illegal pair_coeff command");
}
// error not to specify normal or tangential model
if ((normal_model_one < 0) || (tangential_model_one < 0))
error->all(FLERR, "Illegal pair coeff command, "
error->all(FLERR, "Illegal pair_coeff command, "
"must specify normal or tangential contact model");
int count = 0;
@ -984,6 +990,14 @@ void PairGranular::coeff(int narg, char **arg)
27.467*powint(cor,4)-18.022*powint(cor,5)+4.8218*powint(cor,6);
} else damp = normal_coeffs_one[1];
if (ld_flag && normal_model_one == JKR)
error->all(FLERR,"Illegal pair_coeff command, "
"Cannot limit damping with JKR model");
if (ld_flag && normal_model_one == DMT)
error->all(FLERR,"Illegal pair_coeff command, "
"Cannot limit damping with DMT model");
for (int i = ilo; i <= ihi; i++) {
for (int j = MAX(jlo,i); j <= jhi; j++) {
normal_model[i][j] = normal_model[j][i] = normal_model_one;
@ -1026,11 +1040,14 @@ void PairGranular::coeff(int narg, char **arg)
cutoff_type[i][j] = cutoff_type[j][i] = cutoff_one;
limit_damping[i][j] = limit_damping[j][i] = ld_flag;
setflag[i][j] = 1;
count++;
}
}
if (count == 0) error->all(FLERR,"Incorrect args for pair coefficients");
}
@ -1303,6 +1320,7 @@ void PairGranular::write_restart(FILE *fp)
fwrite(&tangential_model[i][j],sizeof(int),1,fp);
fwrite(&roll_model[i][j],sizeof(int),1,fp);
fwrite(&twist_model[i][j],sizeof(int),1,fp);
fwrite(&limit_damping[i][j],sizeof(int),1,fp);
fwrite(normal_coeffs[i][j],sizeof(double),4,fp);
fwrite(tangential_coeffs[i][j],sizeof(double),3,fp);
fwrite(roll_coeffs[i][j],sizeof(double),3,fp);
@ -1333,6 +1351,7 @@ void PairGranular::read_restart(FILE *fp)
utils::sfread(FLERR,&tangential_model[i][j],sizeof(int),1,fp,nullptr,error);
utils::sfread(FLERR,&roll_model[i][j],sizeof(int),1,fp,nullptr,error);
utils::sfread(FLERR,&twist_model[i][j],sizeof(int),1,fp,nullptr,error);
utils::sfread(FLERR,&limit_damping[i][j],sizeof(int),1,fp,nullptr,error);
utils::sfread(FLERR,normal_coeffs[i][j],sizeof(double),4,fp,nullptr,error);
utils::sfread(FLERR,tangential_coeffs[i][j],sizeof(double),3,fp,nullptr,error);
utils::sfread(FLERR,roll_coeffs[i][j],sizeof(double),3,fp,nullptr,error);
@ -1344,6 +1363,7 @@ void PairGranular::read_restart(FILE *fp)
MPI_Bcast(&tangential_model[i][j],1,MPI_INT,0,world);
MPI_Bcast(&roll_model[i][j],1,MPI_INT,0,world);
MPI_Bcast(&twist_model[i][j],1,MPI_INT,0,world);
MPI_Bcast(&limit_damping[i][j],1,MPI_INT,0,world);
MPI_Bcast(normal_coeffs[i][j],4,MPI_DOUBLE,0,world);
MPI_Bcast(tangential_coeffs[i][j],3,MPI_DOUBLE,0,world);
MPI_Bcast(roll_coeffs[i][j],3,MPI_DOUBLE,0,world);
@ -1529,6 +1549,7 @@ double PairGranular::single(int i, int j, int itype, int jtype,
Fdamp = -damp_normal_prefactor*vnnr;
Fntot = Fne + Fdamp;
if (limit_damping[itype][jtype] && (Fntot < 0.0)) Fntot = 0.0;
jnum = list->numneigh[i];
jlist = list->firstneigh[i];

1
src/GRANULAR/pair_granular.h
View File

@ -70,6 +70,7 @@ class PairGranular : public Pair {
// model choices
int **normal_model, **damping_model;
int **tangential_model, **roll_model, **twist_model;
int **limit_damping;
// history flags
int normal_history, tangential_history, roll_history, twist_history;

1
src/KOKKOS/pair_gran_hooke_history_kokkos.cpp
View File

@ -392,6 +392,7 @@ void PairGranHookeHistoryKokkos<DeviceType>::operator()(TagPairGranHookeHistoryC
F_FLOAT damp = meff*gamman*vnnr*rsqinv;
F_FLOAT ccel = kn*(radsum-r)*rinv - damp;
if(limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities

16
src/KSPACE/pair_lj_cut_coul_long.cpp
View File

@ -17,21 +17,21 @@
#include "pair_lj_cut_coul_long.h"
#include <cmath>
#include <cstring>
#include "atom.h"
#include "comm.h"
#include "error.h"
#include "force.h"
#include "kspace.h"
#include "update.h"
#include "respa.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "math_const.h"
#include "memory.h"
#include "error.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "neighbor.h"
#include "respa.h"
#include "update.h"
#include <cmath>
#include <cstring>
using namespace LAMMPS_NS;
using namespace MathConst;

49
src/PYTHON/fix_python_invoke.cpp
View File

@ -21,6 +21,7 @@
#include "error.h"
#include "lmppython.h"
#include "python_compat.h"
#include "python_utils.h"
#include "update.h"
#include <cstring>
@ -51,12 +52,12 @@ FixPythonInvoke::FixPythonInvoke(LAMMPS *lmp, int narg, char **arg) :
}
// get Python function
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
PyObject *pyMain = PyImport_AddModule("__main__");
if (!pyMain) {
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not initialize embedded Python");
}
@ -64,11 +65,19 @@ FixPythonInvoke::FixPythonInvoke(LAMMPS *lmp, int narg, char **arg) :
pFunc = PyObject_GetAttrString(pyMain, fname);
if (!pFunc) {
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not find Python function");
}
PyGILState_Release(gstate);
lmpPtr = PY_VOID_POINTER(lmp);
}
/* ---------------------------------------------------------------------- */
FixPythonInvoke::~FixPythonInvoke()
{
PyUtils::GIL lock;
Py_CLEAR(lmpPtr);
}
/* ---------------------------------------------------------------------- */
@ -82,18 +91,16 @@ int FixPythonInvoke::setmask()
void FixPythonInvoke::end_of_step()
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
PyObject *ptr = PY_VOID_POINTER(lmp);
PyObject *arglist = Py_BuildValue("(O)", ptr);
PyObject * result = PyObject_CallFunction((PyObject*)pFunc, "O", (PyObject*)lmpPtr);
PyObject *result = PyEval_CallObject((PyObject*)pFunc, arglist);
Py_DECREF(arglist);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result)
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/invoke end_of_step() method failed");
}
Py_CLEAR(result);
}
/* ---------------------------------------------------------------------- */
@ -102,16 +109,14 @@ void FixPythonInvoke::post_force(int vflag)
{
if (update->ntimestep % nevery != 0) return;
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
PyObject *ptr = PY_VOID_POINTER(lmp);
PyObject *arglist = Py_BuildValue("(Oi)", ptr, vflag);
PyObject * result = PyObject_CallFunction((PyObject*)pFunc, "Oi", (PyObject*)lmpPtr, vflag);
PyObject *result = PyEval_CallObject((PyObject*)pFunc, arglist);
Py_DECREF(arglist);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result)
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/invoke post_force() method failed");
}
Py_CLEAR(result);
}

4
src/PYTHON/fix_python_invoke.h
View File

@ -21,6 +21,7 @@ FixStyle(python,FixPythonInvoke)
#ifndef LMP_FIX_PYTHON_INVOKE_H
#define LMP_FIX_PYTHON_INVOKE_H
#include "fix.h"
namespace LAMMPS_NS {
@ -28,12 +29,13 @@ namespace LAMMPS_NS {
class FixPythonInvoke : public Fix {
public:
FixPythonInvoke(class LAMMPS *, int, char **);
virtual ~FixPythonInvoke() {}
virtual ~FixPythonInvoke();
int setmask();
virtual void end_of_step();
virtual void post_force(int);
private:
void * lmpPtr;
void * pFunc;
int selected_callback;
};

177
src/PYTHON/fix_python_move.cpp
View File

@ -21,8 +21,9 @@
#include "error.h"
#include "lmppython.h"
#include "python_compat.h"
#include "python_utils.h"
#include <cstring>
#include <string>
#include <Python.h> // IWYU pragma: export
using namespace LAMMPS_NS;
@ -40,7 +41,7 @@ FixPythonMove::FixPythonMove(LAMMPS *lmp, int narg, char **arg) :
py_move = nullptr;
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
// add current directory to PYTHONPATH
PyObject *py_path = PySys_GetObject((char *)"path");
@ -48,70 +49,50 @@ FixPythonMove::FixPythonMove(LAMMPS *lmp, int narg, char **arg) :
// create integrator instance
char *full_cls_name = arg[3];
char *lastpos = strrchr(full_cls_name, '.');
std::string full_cls_name = arg[3];
size_t lastpos = full_cls_name.rfind(".");
if (lastpos == nullptr) {
if (lastpos == std::string::npos) {
error->all(FLERR,"Fix python/integrate requires fully qualified class name");
}
size_t module_name_length = strlen(full_cls_name) - strlen(lastpos);
size_t cls_name_length = strlen(lastpos)-1;
std::string module_name = full_cls_name.substr(0, lastpos);
std::string cls_name = full_cls_name.substr(lastpos+1);
char *module_name = new char[module_name_length+1];
char *cls_name = new char[cls_name_length+1];
strncpy(module_name, full_cls_name, module_name_length);
module_name[module_name_length] = 0;
strcpy(cls_name, lastpos+1);
PyObject *pModule = PyImport_ImportModule(module_name);
PyObject *pModule = PyImport_ImportModule(module_name.c_str());
if (!pModule) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Loading python integrator module failure");
}
// create LAMMPS atom type to potential file type mapping in python class
// by calling 'lammps_pair_style.map_coeff(name,type)'
PyObject *py_move_type = PyObject_GetAttrString(pModule, cls_name);
PyObject *py_move_type = PyObject_GetAttrString(pModule, cls_name.c_str());
if (!py_move_type) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not find integrator class in module'");
}
delete [] module_name;
delete [] cls_name;
PyObject *ptr = PY_VOID_POINTER(lmp);
PyObject *arglist = Py_BuildValue("(O)", ptr);
PyObject *py_move_obj = PyObject_CallObject(py_move_type, arglist);
Py_DECREF(arglist);
PyObject *py_move_obj = PyObject_CallFunction(py_move_type, "O", ptr);
Py_CLEAR(ptr);
if (!py_move_obj) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not instantiate instance of integrator class'");
}
// check object interface
py_move = (void *) py_move_obj;
PyGILState_Release(gstate);
}
/* ---------------------------------------------------------------------- */
FixPythonMove::~FixPythonMove()
{
PyGILState_STATE gstate = PyGILState_Ensure();
if (py_move) Py_DECREF((PyObject*) py_move);
PyGILState_Release(gstate);
PyUtils::GIL lock;
Py_CLEAR(py_move);
}
/* ---------------------------------------------------------------------- */
@ -130,122 +111,82 @@ int FixPythonMove::setmask()
void FixPythonMove::init()
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject *py_move_obj = (PyObject *) py_move;
PyObject *py_init = PyObject_GetAttrString(py_move_obj,(char *)"init");
if (!py_init) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'init()' method'");
PyUtils::GIL lock;
PyObject * result = PyObject_CallMethod((PyObject *)py_move, "init", nullptr);
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/move init() method failed");
}
PyObject *result = PyEval_CallObject(py_init, nullptr);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result) error->all(FLERR,"Fix python/move init() method failed");
Py_CLEAR(result);
}
/* ---------------------------------------------------------------------- */
void FixPythonMove::initial_integrate(int vflag)
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject *py_move_obj = (PyObject *) py_move;
PyObject *py_initial_integrate = PyObject_GetAttrString(py_move_obj,"initial_integrate");
if (!py_initial_integrate) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'initial_integrate' method'");
PyUtils::GIL lock;
PyObject * result = PyObject_CallMethod((PyObject*)py_move, "initial_integrate", "i", vflag);
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/move initial_integrate() method failed");
}
PyObject *arglist = Py_BuildValue("(i)", vflag);
PyObject *result = PyEval_CallObject(py_initial_integrate, arglist);
Py_DECREF(arglist);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result) error->all(FLERR,"Fix python/move initial_integrate() "
"method failed");
Py_CLEAR(result);
}
/* ---------------------------------------------------------------------- */
void FixPythonMove::final_integrate()
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject *py_move_obj = (PyObject *) py_move;
PyObject *py_final_integrate = PyObject_GetAttrString(py_move_obj,"final_integrate");
if (!py_final_integrate) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'final_integrate' method'");
PyUtils::GIL lock;
PyObject * result = PyObject_CallMethod((PyObject*)py_move, "final_integrate", nullptr);
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/move final_integrate() method failed");
}
PyObject *result = PyEval_CallObject(py_final_integrate, nullptr);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result) error->all(FLERR,"Fix python/move final_integrate() method "
"failed");
Py_CLEAR(result);
}
/* ---------------------------------------------------------------------- */
void FixPythonMove::initial_integrate_respa(int vflag, int ilevel, int iloop)
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject *py_move_obj = (PyObject *) py_move;
PyObject *py_initial_integrate_respa = PyObject_GetAttrString(py_move_obj,"initial_integrate_respa");
if (!py_initial_integrate_respa) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'initial_integrate_respa' method'");
PyUtils::GIL lock;
PyObject * result = PyObject_CallMethod((PyObject*)py_move, "initial_integrate_respa", "iii", vflag, ilevel, iloop);
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/move initial_integrate_respa() method failed");
}
PyObject *arglist = Py_BuildValue("(iii)", vflag, ilevel, iloop);
PyObject *result = PyEval_CallObject(py_initial_integrate_respa, arglist);
Py_DECREF(arglist);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result) error->all(FLERR,"Fix python/move initial_integrate_respa() "
"method failed");
Py_CLEAR(result);
}
/* ---------------------------------------------------------------------- */
void FixPythonMove::final_integrate_respa(int ilevel, int iloop)
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject *py_move_obj = (PyObject *) py_move;
PyObject *py_final_integrate_respa = PyObject_GetAttrString(py_move_obj,"final_integrate_respa");
if (!py_final_integrate_respa) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'final_integrate_respa' method'");
PyUtils::GIL lock;
PyObject * result = PyObject_CallMethod((PyObject*)py_move, "final_integrate_respa", "ii", ilevel, iloop);
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/move final_integrate_respa() method failed");
}
PyObject *arglist = Py_BuildValue("(ii)", ilevel, iloop);
PyObject *result = PyEval_CallObject(py_final_integrate_respa, arglist);
Py_DECREF(arglist);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result) error->all(FLERR,"Fix python/move final_integrate_respa() "
"method failed");
Py_CLEAR(result);
}
/* ---------------------------------------------------------------------- */
void FixPythonMove::reset_dt()
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject *py_move_obj = (PyObject *) py_move;
PyObject *py_reset_dt = PyObject_GetAttrString(py_move_obj,"reset_dt");
if (!py_reset_dt) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'reset_dt' method'");
PyUtils::GIL lock;
PyObject * result = PyObject_CallMethod((PyObject*)py_move, "reset_dt", nullptr);
if (!result) {
PyUtils::Print_Errors();
error->all(FLERR,"Fix python/move reset_dt() method failed");
}
PyObject *result = PyEval_CallObject(py_reset_dt, nullptr);
if (!result && (comm->me == 0)) PyErr_Print();
PyGILState_Release(gstate);
if (!result) error->all(FLERR,"Fix python/move reset_dt() method failed");
Py_CLEAR(result);
}

238
src/PYTHON/pair_python.cpp
View File

@ -24,9 +24,10 @@
#include "memory.h"
#include "neigh_list.h"
#include "python_compat.h"
#include "python_utils.h"
#include "update.h"
#include <cstring>
#include <string>
#include <Python.h> // IWYU pragma: export
using namespace LAMMPS_NS;
@ -50,7 +51,7 @@ PairPython::PairPython(LAMMPS *lmp) : Pair(lmp) {
// add current directory to PYTHONPATH
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
PyObject *py_path = PySys_GetObject((char *)"path");
PyList_Append(py_path, PY_STRING_FROM_STRING("."));
@ -61,14 +62,14 @@ PairPython::PairPython(LAMMPS *lmp) : Pair(lmp) {
if (potentials_path != nullptr) {
PyList_Append(py_path, PY_STRING_FROM_STRING(potentials_path));
}
PyGILState_Release(gstate);
}
/* ---------------------------------------------------------------------- */
PairPython::~PairPython()
{
if (py_potential) Py_DECREF((PyObject*) py_potential);
PyUtils::GIL lock;
Py_CLEAR(py_potential);
delete[] skip_types;
if (allocated) {
@ -103,41 +104,31 @@ void PairPython::compute(int eflag, int vflag)
// prepare access to compute_force and compute_energy functions
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
PyObject *py_pair_instance = (PyObject *) py_potential;
PyObject *py_compute_force = PyObject_GetAttrString(py_pair_instance,"compute_force");
if (!py_compute_force) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not find 'compute_force' method'");
}
if (!PyCallable_Check(py_compute_force)) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Python 'compute_force' is not callable");
}
PyObject *py_compute_energy = PyObject_GetAttrString(py_pair_instance,"compute_energy");
if (!py_compute_energy) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not find 'compute_energy' method'");
}
if (!PyCallable_Check(py_compute_energy)) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Python 'compute_energy' is not callable");
}
PyObject *py_compute_args = PyTuple_New(3);
if (!py_compute_args) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not create tuple for 'compute' function arguments");
}
@ -179,13 +170,11 @@ void PairPython::compute(int eflag, int vflag)
PyTuple_SetItem(py_compute_args,0,py_rsq);
py_value = PyObject_CallObject(py_compute_force,py_compute_args);
if (!py_value) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Calling 'compute_force' function failed");
}
fpair = factor_lj*PyFloat_AsDouble(py_value);
Py_DECREF(py_value);
Py_CLEAR(py_value);
f[i][0] += delx*fpair;
f[i][1] += dely*fpair;
@ -198,8 +187,12 @@ void PairPython::compute(int eflag, int vflag)
if (eflag) {
py_value = PyObject_CallObject(py_compute_energy,py_compute_args);
if (!py_value) {
PyUtils::Print_Errors();
error->all(FLERR,"Calling 'compute_energy' function failed");
}
evdwl = factor_lj*PyFloat_AsDouble(py_value);
Py_DECREF(py_value);
Py_CLEAR(py_value);
} else evdwl = 0.0;
if (evflag) ev_tally(i,j,nlocal,newton_pair,
@ -207,8 +200,7 @@ void PairPython::compute(int eflag, int vflag)
}
}
}
Py_DECREF(py_compute_args);
PyGILState_Release(gstate);
Py_CLEAR(py_compute_args);
if (vflag_fdotr) virial_fdotr_compute();
}
@ -261,112 +253,49 @@ void PairPython::coeff(int narg, char **arg)
error->all(FLERR,"Incorrect args for pair coefficients");
// check if python potential file exists and source it
char * full_cls_name = arg[2];
char * lastpos = strrchr(full_cls_name, '.');
std::string full_cls_name = arg[2];
size_t lastpos = full_cls_name.rfind(".");
if (lastpos == nullptr) {
if (lastpos == std::string::npos) {
error->all(FLERR,"Python pair style requires fully qualified class name");
}
size_t module_name_length = strlen(full_cls_name) - strlen(lastpos);
size_t cls_name_length = strlen(lastpos)-1;
std::string module_name = full_cls_name.substr(0, lastpos);
std::string cls_name = full_cls_name.substr(lastpos+1);
char * module_name = new char[module_name_length+1];
char * cls_name = new char[cls_name_length+1];
strncpy(module_name, full_cls_name, module_name_length);
module_name[module_name_length] = 0;
PyUtils::GIL lock;
strcpy(cls_name, lastpos+1);
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject * pModule = PyImport_ImportModule(module_name);
PyObject * pModule = PyImport_ImportModule(module_name.c_str());
if (!pModule) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Loading python pair style module failure");
}
// create LAMMPS atom type to potential file type mapping in python class
// by calling 'lammps_pair_style.map_coeff(name,type)'
PyObject *py_pair_type = PyObject_GetAttrString(pModule, cls_name);
PyObject *py_pair_type = PyObject_GetAttrString(pModule, cls_name.c_str());
if (!py_pair_type) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not find pair style class in module'");
}
delete [] module_name;
delete [] cls_name;
PyObject * py_pair_instance = PyObject_CallObject(py_pair_type, nullptr);
if (!py_pair_instance) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not instantiate instance of pair style class'");
}
py_potential = (void *) py_pair_instance;
PyObject *py_check_units = PyObject_GetAttrString(py_pair_instance,"check_units");
if (!py_check_units) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'check_units' method'");
}
if (!PyCallable_Check(py_check_units)) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Python 'check_units' is not callable");
}
PyObject *py_units_args = PyTuple_New(1);
if (!py_units_args) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not create tuple for 'check_units' function arguments");
}
PyObject *py_name = PY_STRING_FROM_STRING(update->unit_style);
PyTuple_SetItem(py_units_args,0,py_name);
PyObject *py_value = PyObject_CallObject(py_check_units,py_units_args);
PyObject *py_value = PyObject_CallMethod(py_pair_instance, "check_units", "s", update->unit_style);
if (!py_value) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Calling 'check_units' function failed");
}
Py_DECREF(py_units_args);
Py_CLEAR(py_value);
PyObject *py_map_coeff = PyObject_GetAttrString(py_pair_instance,"map_coeff");
if (!py_map_coeff) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'map_coeff' method'");
}
if (!PyCallable_Check(py_map_coeff)) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Python 'map_coeff' is not callable");
}
PyObject *py_map_args = PyTuple_New(2);
if (!py_map_args) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not create tuple for 'map_coeff' function arguments");
}
delete[] skip_types;
skip_types = new int[ntypes+1];
skip_types[0] = 1;
@ -375,25 +304,20 @@ void PairPython::coeff(int narg, char **arg)
skip_types[i] = 1;
continue;
} else skip_types[i] = 0;
PyObject *py_type = PY_INT_FROM_LONG(i);
py_name = PY_STRING_FROM_STRING(arg[2+i]);
PyTuple_SetItem(py_map_args,0,py_name);
PyTuple_SetItem(py_map_args,1,py_type);
py_value = PyObject_CallObject(py_map_coeff,py_map_args);
const int type = i;
const char * name = arg[2+i];
py_value = PyObject_CallMethod(py_pair_instance, "map_coeff", "si", name, type);
if (!py_value) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Calling 'map_coeff' function failed");
}
Py_CLEAR(py_value);
for (int j = i; j <= ntypes ; j++) {
setflag[i][j] = 1;
cutsq[i][j] = cut_global*cut_global;
}
}
Py_DECREF(py_map_args);
PyGILState_Release(gstate);
}
/* ---------------------------------------------------------------------- */
@ -417,76 +341,52 @@ double PairPython::single(int /* i */, int /* j */, int itype, int jtype,
// prepare access to compute_force and compute_energy functions
PyGILState_STATE gstate = PyGILState_Ensure();
PyObject *py_pair_instance = (PyObject *) py_potential;
PyObject *py_compute_force
= PyObject_GetAttrString(py_pair_instance,"compute_force");
if (!py_compute_force) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'compute_force' method'");
}
if (!PyCallable_Check(py_compute_force)) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Python 'compute_force' is not callable");
}
PyUtils::GIL lock;
PyObject *py_compute_force = (PyObject *) get_member_function("compute_force");
PyObject *py_compute_energy = (PyObject *) get_member_function("compute_energy");
PyObject *py_compute_args = Py_BuildValue("(dii)", rsq, itype, jtype);
PyObject *py_compute_energy
= PyObject_GetAttrString(py_pair_instance,"compute_energy");
if (!py_compute_energy) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Could not find 'compute_energy' method'");
}
if (!PyCallable_Check(py_compute_energy)) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
error->all(FLERR,"Python 'compute_energy' is not callable");
}
PyObject *py_rsq, *py_itype, *py_jtype, *py_value;
PyObject *py_compute_args = PyTuple_New(3);
if (!py_compute_args) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not create tuple for 'compute' function arguments");
}
py_itype = PY_INT_FROM_LONG(itype);
PyTuple_SetItem(py_compute_args,1,py_itype);
py_jtype = PY_INT_FROM_LONG(jtype);
PyTuple_SetItem(py_compute_args,2,py_jtype);
py_rsq = PyFloat_FromDouble(rsq);
PyTuple_SetItem(py_compute_args,0,py_rsq);
py_value = PyObject_CallObject(py_compute_force,py_compute_args);
PyObject * py_value = PyObject_CallObject(py_compute_force, py_compute_args);
if (!py_value) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Calling 'compute_force' function failed");
}
fforce = factor_lj*PyFloat_AsDouble(py_value);
Py_DECREF(py_value);
Py_CLEAR(py_value);
py_value = PyObject_CallObject(py_compute_energy,py_compute_args);
py_value = PyObject_CallObject(py_compute_energy, py_compute_args);
if (!py_value) {
PyErr_Print();
PyErr_Clear();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Calling 'compute_energy' function failed");
}
double evdwl = factor_lj*PyFloat_AsDouble(py_value);
Py_DECREF(py_value);
Py_DECREF(py_compute_args);
PyGILState_Release(gstate);
Py_CLEAR(py_value);
Py_CLEAR(py_compute_args);
return evdwl;
}
/* ---------------------------------------------------------------------- */
void * PairPython::get_member_function(const char * name)
{
PyUtils::GIL lock;
PyObject *py_pair_instance = (PyObject *) py_potential;
PyObject * py_mfunc = PyObject_GetAttrString(py_pair_instance, name);
if (!py_mfunc) {
PyUtils::Print_Errors();
error->all(FLERR, fmt::format("Could not find '{}' method'", name));
}
if (!PyCallable_Check(py_mfunc)) {
PyUtils::Print_Errors();
error->all(FLERR, fmt::format("Python '{}' is not callable", name));
}
return py_mfunc;
}

1
src/PYTHON/pair_python.h
View File

@ -50,6 +50,7 @@ class PairPython : public Pair {
int * skip_types;
virtual void allocate();
void * get_member_function(const char *);
};
}

62
src/PYTHON/python_impl.cpp
View File

@ -21,6 +21,7 @@
#include "input.h"
#include "memory.h"
#include "python_compat.h"
#include "python_utils.h"
#include "variable.h"
#include <cstring>
@ -91,13 +92,12 @@ PythonImpl::PythonImpl(LAMMPS *lmp) : Pointers(lmp)
}
#endif
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
PyObject *pModule = PyImport_AddModule("__main__");
if (!pModule) error->all(FLERR,"Could not initialize embedded Python");
pyMain = (void *) pModule;
PyGILState_Release(gstate);
}
/* ---------------------------------------------------------------------- */
@ -106,23 +106,19 @@ PythonImpl::~PythonImpl()
{
if (pyMain) {
// clean up
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
for (int i = 0; i < nfunc; i++) {
delete [] pfuncs[i].name;
deallocate(i);
PyObject *pFunc = (PyObject *) pfuncs[i].pFunc;
Py_XDECREF(pFunc);
Py_CLEAR(pfuncs[i].pFunc);
}
}
// shutdown Python interpreter
if (!external_interpreter) {
Py_Finalize();
}
else {
PyGILState_Release(gstate);
}
// shutdown Python interpreter
if (!external_interpreter) {
PyGILState_Ensure();
Py_Finalize();
}
memory->sfree(pfuncs);
@ -228,7 +224,7 @@ void PythonImpl::command(int narg, char **arg)
int ifunc = create_entry(arg[0]);
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
// send Python code to Python interpreter
// file: read the file via PyRun_SimpleFile()
@ -239,14 +235,14 @@ void PythonImpl::command(int narg, char **arg)
FILE *fp = fopen(pyfile,"r");
if (fp == nullptr) {
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not open Python file");
}
int err = PyRun_SimpleFile(fp,pyfile);
if (err) {
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not process Python file");
}
@ -255,7 +251,7 @@ void PythonImpl::command(int narg, char **arg)
int err = PyRun_SimpleString(herestr);
if (err) {
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,"Could not process Python string");
}
}
@ -266,13 +262,13 @@ void PythonImpl::command(int narg, char **arg)
PyObject *pFunc = PyObject_GetAttrString(pModule,pfuncs[ifunc].name);
if (!pFunc) {
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,fmt::format("Could not find Python function {}",
pfuncs[ifunc].name));
}
if (!PyCallable_Check(pFunc)) {
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->all(FLERR,fmt::format("Python function {} is not callable",
pfuncs[ifunc].name));
}
@ -284,14 +280,13 @@ void PythonImpl::command(int narg, char **arg)
delete [] istr;
delete [] format;
delete [] pyfile;
PyGILState_Release(gstate);
}
/* ------------------------------------------------------------------ */
void PythonImpl::invoke_function(int ifunc, char *result)
{
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
PyObject *pValue;
char *str;
@ -303,7 +298,6 @@ void PythonImpl::invoke_function(int ifunc, char *result)
PyObject *pArgs = PyTuple_New(ninput);
if (!pArgs) {
PyGILState_Release(gstate);
error->all(FLERR,"Could not create Python function arguments");
}
@ -314,7 +308,6 @@ void PythonImpl::invoke_function(int ifunc, char *result)
str = input->variable->retrieve(pfuncs[ifunc].svalue[i]);
if (!str) {
PyGILState_Release(gstate);
error->all(FLERR,"Could not evaluate Python function input variable");
}
@ -327,7 +320,6 @@ void PythonImpl::invoke_function(int ifunc, char *result)
str = input->variable->retrieve(pfuncs[ifunc].svalue[i]);
if (!str) {
PyGILState_Release(gstate);
error->all(FLERR,"Could not evaluate Python function input variable");
}
@ -339,7 +331,6 @@ void PythonImpl::invoke_function(int ifunc, char *result)
if (pfuncs[ifunc].ivarflag[i]) {
str = input->variable->retrieve(pfuncs[ifunc].svalue[i]);
if (!str) {
PyGILState_Release(gstate);
error->all(FLERR,"Could not evaluate Python function input variable");
}
@ -350,7 +341,6 @@ void PythonImpl::invoke_function(int ifunc, char *result)
} else if (itype == PTR) {
pValue = PY_VOID_POINTER(lmp);
} else {
PyGILState_Release(gstate);
error->all(FLERR,"Unsupported variable type");
}
PyTuple_SetItem(pArgs,i,pValue);
@ -360,15 +350,13 @@ void PythonImpl::invoke_function(int ifunc, char *result)
// error check with one() since only some procs may fail
pValue = PyObject_CallObject(pFunc,pArgs);
Py_CLEAR(pArgs);
if (!pValue) {
PyErr_Print();
PyGILState_Release(gstate);
PyUtils::Print_Errors();
error->one(FLERR,"Python function evaluation failed");
}
Py_DECREF(pArgs);
// function returned a value
// assign it to result string stored by python-style variable
// or if user specified a length, assign it to longstr
@ -385,10 +373,8 @@ void PythonImpl::invoke_function(int ifunc, char *result)
strncpy(pfuncs[ifunc].longstr,pystr,pfuncs[ifunc].length_longstr);
else strncpy(result,pystr,VALUELENGTH-1);
}
Py_DECREF(pValue);
}
PyGILState_Release(gstate);
Py_CLEAR(pValue);
}
/* ------------------------------------------------------------------ */
@ -523,11 +509,8 @@ int PythonImpl::create_entry(char *name)
int PythonImpl::execute_string(char *cmd)
{
PyGILState_STATE gstate = PyGILState_Ensure();
int err = PyRun_SimpleString(cmd);
PyGILState_Release(gstate);
return err;
PyUtils::GIL lock;
return PyRun_SimpleString(cmd);
}
/* ---------------------------------------------------------------------- */
@ -537,9 +520,8 @@ int PythonImpl::execute_file(char *fname)
FILE *fp = fopen(fname,"r");
if (fp == nullptr) return -1;
PyGILState_STATE gstate = PyGILState_Ensure();
PyUtils::GIL lock;
int err = PyRun_SimpleFile(fp,fname);
PyGILState_Release(gstate);
if (fp) fclose(fp);
return err;

42
src/PYTHON/python_utils.h Normal file
View File

@ -0,0 +1,42 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifndef LMP_PYTHON_UTILS_H
#define LMP_PYTHON_UTILS_H
#include <Python.h>
namespace LAMMPS_NS {
namespace PyUtils {
class GIL {
PyGILState_STATE gstate;
public:
GIL() : gstate(PyGILState_Ensure()) {
}
~GIL() {
PyGILState_Release(gstate);
}
};
static void Print_Errors() {
PyErr_Print();
PyErr_Clear();
}
}
}
#endif

6
src/RIGID/fix_rigid_nh.cpp
View File

@ -234,8 +234,6 @@ void FixRigidNH::init()
}
g_f = nf_t + nf_r;
onednft = 1.0 + (double)(dimension) / (double)g_f;
onednfr = (double) (dimension) / (double)g_f;
// see Table 1 in Kamberaj et al
@ -719,6 +717,8 @@ void FixRigidNH::final_integrate()
void FixRigidNH::nhc_temp_integrate()
{
if (g_f == 0) return;
int i,j,k;
double kt,gfkt_t,gfkt_r,tmp,ms,s,s2;
@ -1063,6 +1063,8 @@ void FixRigidNH::compute_press_target()
void FixRigidNH::nh_epsilon_dot()
{
if (g_f == 0) return;
int i;
double volume,scale,f_epsilon;

2
src/RIGID/fix_rigid_nh.h
View File

@ -38,8 +38,6 @@ class FixRigidNH : public FixRigid {
double boltz,nktv2p,mvv2e; // boltzman constant, conversion factors
int nf_t,nf_r; // trans/rot degrees of freedom
double onednft,onednfr; // factors 1 + dimension/trans(rot)
// degrees of freedom
double *w,*wdti1,*wdti2,*wdti4; // Yoshida-Suzuki coefficients
double *q_t,*q_r; // trans/rot thermostat masses
double *eta_t,*eta_r; // trans/rot thermostat positions

19
src/RIGID/fix_rigid_nh_small.cpp
View File

@ -754,6 +754,8 @@ void FixRigidNHSmall::final_integrate()
void FixRigidNHSmall::nhc_temp_integrate()
{
if (g_f == 0) return;
int i,j,k;
double kt,gfkt_t,gfkt_r,tmp,ms,s,s2;
@ -1148,6 +1150,8 @@ void FixRigidNHSmall::compute_press_target()
void FixRigidNHSmall::nh_epsilon_dot()
{
if (g_f == 0) return;
int i;
double volume,scale,f_epsilon;
@ -1204,8 +1208,6 @@ void FixRigidNHSmall::compute_dof()
nf_r = nfall[1];
g_f = nf_t + nf_r;
onednft = 1.0 + (double)(dimension) / (double)g_f;
onednfr = (double) (dimension) / (double)g_f;
}
/* ----------------------------------------------------------------------
@ -1367,19 +1369,6 @@ int FixRigidNHSmall::modify_param(int narg, char **arg)
return FixRigidSmall::modify_param(narg,arg);
}
/* ----------------------------------------------------------------------
disallow using fix rigid/n??/small fixes with fix deposit
we would need custom functionality to update data structures
used by all fixes derived from this class but not fix rigid/small
------------------------------------------------------------------------- */
void FixRigidNHSmall::set_molecule(int, tagint, int,
double *, double *, double *)
{
error->all(FLERR,fmt::format("Molecule update not (yet) implemented for "
"fix {}", style));
}
/* ---------------------------------------------------------------------- */
void FixRigidNHSmall::allocate_chain()

2
src/RIGID/fix_rigid_nh_small.h
View File

@ -38,7 +38,6 @@ class FixRigidNHSmall : public FixRigidSmall {
int dimension; // # of dimensions
int nf_t,nf_r; // trans/rot degrees of freedom
double onednft,onednfr; // factors 1 + dimension/trans(rot) degrees of freedom
double *w,*wdti1,*wdti2,*wdti4; // Yoshida-Suzuki coefficients
double *q_t,*q_r; // trans/rot thermostat masses
double *eta_t,*eta_r; // trans/rot thermostat positions
@ -80,7 +79,6 @@ class FixRigidNHSmall : public FixRigidSmall {
void nh_epsilon_dot();
void compute_dof();
void set_molecule(int, tagint, int, double *, double *, double *);
void allocate_chain();
void allocate_order();
void deallocate_chain();

1
src/USER-OMP/pair_gran_hertz_history_omp.cpp
View File

@ -224,6 +224,7 @@ void PairGranHertzHistoryOMP::eval(int iifrom, int iito, ThrData * const thr)
ccel = kn*(radsum-r)*rinv - damp;
polyhertz = sqrt((radsum-r)*radi*radj / radsum);
ccel *= polyhertz;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities

1
src/USER-OMP/pair_gran_hooke_history_omp.cpp
View File

@ -223,6 +223,7 @@ void PairGranHookeHistoryOMP::eval(int iifrom, int iito, ThrData * const thr)
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radsum-r)*rinv - damp;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities

1
src/USER-OMP/pair_gran_hooke_omp.cpp
View File

@ -197,6 +197,7 @@ void PairGranHookeOMP::eval(int iifrom, int iito, ThrData * const thr)
damp = meff*gamman*vnnr*rsqinv;
ccel = kn*(radsum-r)*rinv - damp;
if (limit_damping && (ccel < 0.0)) ccel = 0.0;
// relative velocities

1
src/USER-QUIP/pair_quip.cpp
View File

@ -318,7 +318,6 @@ void PairQUIP::init_style()
// Initialise neighbor list
int irequest_full = neighbor->request(this);
neighbor->requests[irequest_full]->id = 1;
neighbor->requests[irequest_full]->half = 0;
neighbor->requests[irequest_full]->full = 1;
}

45
src/USER-REAXC/reaxc_tool_box.cpp
View File

@ -25,10 +25,11 @@
----------------------------------------------------------------------*/
#include "reaxc_tool_box.h"
#include "reaxc_defs.h"
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include "reaxc_defs.h"
#include "error.h"
@ -49,62 +50,58 @@ int Tokenize( char* s, char*** tok )
return count;
}
/* safe malloc */
void *smalloc( LAMMPS_NS::Error *error_ptr, rc_bigint n, const char *name )
{
void *ptr;
char errmsg[256];
if (n <= 0) {
snprintf(errmsg, 256, "Trying to allocate %ld bytes for array %s. "
"returning NULL.", n, name);
auto errmsg = fmt::format("Trying to allocate {} bytes for array {}. "
"returning NULL.", n, name);
if (error_ptr) error_ptr->one(FLERR,errmsg);
else fputs(errmsg,stderr);
else fputs(errmsg.c_str(),stderr);
return nullptr;
}
ptr = malloc( n );
if (ptr == nullptr) {
snprintf(errmsg, 256, "Failed to allocate %ld bytes for array %s", n, name);
auto errmsg = fmt::format("Failed to allocate {} bytes for array {}",
n, name);
if (error_ptr) error_ptr->one(FLERR,errmsg);
else fputs(errmsg,stderr);
else fputs(errmsg.c_str(),stderr);
}
return ptr;
}
/* safe calloc */
void *scalloc( LAMMPS_NS::Error *error_ptr, rc_bigint n, rc_bigint size, const char *name )
{
void *ptr;
char errmsg[256];
if (n <= 0) {
snprintf(errmsg, 256, "Trying to allocate %ld elements for array %s. "
"returning NULL.\n", n, name );
auto errmsg = fmt::format("Trying to allocate {} elements for array {}. "
"returning NULL.\n", n, name);
if (error_ptr) error_ptr->one(FLERR,errmsg);
else fputs(errmsg,stderr);
else fputs(errmsg.c_str(),stderr);
return nullptr;
}
if (size <= 0) {
snprintf(errmsg, 256, "Elements size for array %s is %ld. "
"returning NULL", name, size );
auto errmsg = fmt::format("Elements size for array {} is {}. "
"returning NULL", name, size);
if (error_ptr) error_ptr->one(FLERR,errmsg);
else fputs(errmsg,stderr);
else fputs(errmsg.c_str(),stderr);
return nullptr;
}
ptr = calloc( n, size );
if (ptr == nullptr) {
char errmsg[256];
snprintf(errmsg, 256, "Failed to allocate %ld bytes for array %s", n*size, name);
auto errmsg = fmt::format("Failed to allocate {} bytes for array {}",
n*size, name);
if (error_ptr) error_ptr->one(FLERR,errmsg);
else fputs(errmsg,stderr);
else fputs(errmsg.c_str(),stderr);
}
return ptr;
@ -115,14 +112,14 @@ void *scalloc( LAMMPS_NS::Error *error_ptr, rc_bigint n, rc_bigint size, const c
void sfree( LAMMPS_NS::Error* error_ptr, void *ptr, const char *name )
{
if (ptr == nullptr) {
char errmsg[256];
snprintf(errmsg, 256, "Trying to free the already NULL pointer %s", name );
auto errmsg = fmt::format("Trying to free the already free()'d pointer {}",
name);
if (error_ptr) error_ptr->one(FLERR,errmsg);
else fputs(errmsg,stderr);
else fputs(errmsg.c_str(),stderr);
return;
}
free( ptr );
free(ptr);
ptr = nullptr;
}

16
src/USER-REAXC/reaxc_traj.h
View File

@ -36,14 +36,18 @@
#define HEADER_LINE_LEN 62
#define STR_LINE "%-37s%-24s\n"
#define INT_LINE "%-37s%-24d\n"
#define BIGINT_LINE "%-37s%-24ld\n"
#if defined(LAMMPS_SMALLSMALL)
#define BIGINT_LINE "%-37s%-24d\n"
#else
#define BIGINT_LINE "%-37s%-24lld\n"
#endif
#define INT2_LINE "%-36s%-12d,%-12d\n"
#define REAL_LINE "%-37s%-24.3f\n"
#define SCI_LINE "%-37s%-24g\n"
#define REAL3_LINE "%-32s%9.3f,%9.3f,%9.3f\n"
#if defined(LAMMPS_BIGBIG)
#define INIT_DESC "%9ld%3d%9s%10.3f\n" //AtomID - AtomType, AtomName, AtomMass
#define INIT_DESC "%9lld%3d%9s%10.3f\n" //AtomID - AtomType, AtomName, AtomMass
#else
#define INIT_DESC "%9d%3d%9s%10.3f\n" //AtomID - AtomType, AtomName, AtomMass
#endif
@ -56,7 +60,7 @@
#define SIZE_INFO_LEN3 31
#if defined(LAMMPS_BIGBIG)
#define ATOM_BASIC "%9ld%10.3f%10.3f%10.3f%10.3f\n" //AtomID AtomType (X Y Z) Charge
#define ATOM_BASIC "%9lld%10.3f%10.3f%10.3f%10.3f\n" //AtomID AtomType (X Y Z) Charge
#define ATOM_wV "%9ld%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f\n" //AtomID (X Y Z) (Vx Vy Vz) Charge
#define ATOM_wF "%9ld%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f\n" //AtomID (X Y Z) (Fx Fy Fz) Charge
#define ATOM_FULL "%9ld%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f\n" //AtomID (X Y Z) (Vx Vy Vz) (Fx Fy Fz) Charge
@ -72,9 +76,9 @@
#define ATOM_FULL_LEN 110
#if defined(LAMMPS_BIGBIG)
#define BOND_BASIC "%9ld%9ld%10.3f%10.3f\n" // Atom1 Atom2 Dist Total_BO
#define BOND_FULL "%9ld%9ld%10.3f%10.3f%10.3f%10.3f%10.3f\n" // Atom1 Atom2 Dist Total_BO BOs BOpi BOpi2
#define ANGLE_BASIC "%9ld%9ld%9ld%10.3f\n" // Atom1 Atom2 Atom3 Theta
#define BOND_BASIC "%9lld%9lld%10.3f%10.3f\n" // Atom1 Atom2 Dist Total_BO
#define BOND_FULL "%9lld%9lld%10.3f%10.3f%10.3f%10.3f%10.3f\n" // Atom1 Atom2 Dist Total_BO BOs BOpi BOpi2
#define ANGLE_BASIC "%9lld%9lld%9lld%10.3f\n" // Atom1 Atom2 Atom3 Theta
#else
#define BOND_BASIC "%9d%9d%10.3f%10.3f\n" // Atom1 Atom2 Dist Total_BO
#define BOND_FULL "%9d%9d%10.3f%10.3f%10.3f%10.3f%10.3f\n" // Atom1 Atom2 Dist Total_BO BOs BOpi BOpi2

22
src/USER-YAFF/pair_lj_switch3_coulgauss_long.cpp
View File

@ -106,6 +106,7 @@ void PairLJSwitch3CoulGaussLong::compute(int eflag, int vflag)
firstneigh = list->firstneigh;
// loop over neighbors of my atoms
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
@ -193,7 +194,7 @@ void PairLJSwitch3CoulGaussLong::compute(int eflag, int vflag)
offset[itype][jtype];
} else evdwl = 0.0;
// Truncation, see Yaff Switch33
// Truncation, see Yaff Switch3
if (truncw>0) {
if (rsq < cut_ljsq[itype][jtype]) {
if (r>cut_lj[itype][jtype]-truncw) {
@ -262,19 +263,18 @@ void PairLJSwitch3CoulGaussLong::allocate()
void PairLJSwitch3CoulGaussLong::settings(int narg, char **arg)
{
if (narg < 2 || narg > 3) error->all(FLERR,"Illegal pair_style command");
if (narg < 2 || narg > 3) error->all(FLERR,"Illegal pair_style command");
cut_lj_global = utils::numeric(FLERR,arg[0],false,lmp);
if (narg == 2) {
cut_coul = cut_lj_global;
truncw = utils::numeric(FLERR,arg[1],false,lmp);
}
else {
} else {
cut_coul = utils::numeric(FLERR,arg[1],false,lmp);
truncw = utils::numeric(FLERR,arg[2],false,lmp);
}
if (truncw>0.0) truncwi = 1.0/truncw;
else truncwi = 0.0;
// reset cutoffs that have been explicitly set
if (allocated) {
@ -332,6 +332,9 @@ void PairLJSwitch3CoulGaussLong::init_style()
cut_coulsq = cut_coul * cut_coul;
if (truncw>0.0) truncwi = 1.0/truncw;
else truncwi = 0.0;
// insure use of KSpace long-range solver, set g_ewald
if (force->kspace == nullptr)
@ -375,8 +378,7 @@ double PairLJSwitch3CoulGaussLong::init_one(int i, int j)
double r6inv = r2inv*r2inv*r2inv;
double r12inv = r6inv*r6inv;
offset[i][j] = lj3[i][j]*r12inv-lj4[i][j]*r6inv;
}
else {offset[i][j] = 0.0;}
} else {offset[i][j] = 0.0;}
} else offset[i][j] = 0.0;
cut_ljsq[j][i] = cut_ljsq[i][j];
@ -431,8 +433,7 @@ double PairLJSwitch3CoulGaussLong::init_one(int i, int j)
double t71 = -0.4e1 * cg * (0.2e1 * t10 * t11 - 0.2e1 * t10 * t14 + (cg5 - 0.2e1 * cg1) * t58 * cg5) * t26 / t4 / t41 / t9;
etail_ij = 2.0*MY_PI*all[0]*all[1]*t71;
ptail_ij = 2.0*MY_PI*all[0]*all[1]*t71;
}
else {
} else {
double t1 = pow(cg3, 0.2e1);
double t2 = t1 * t1;
double t3 = t2 * t1;
@ -618,6 +619,7 @@ double PairLJSwitch3CoulGaussLong::single(int i, int j, int itype, int jtype,
expn2 = 0.0;
erfc2 = 0.0;
forcecoul2 = 0.0;
prefactor2 = 0.0;
} else {
r = sqrt(rsq);
rrij = lj2[itype][jtype]*r;

20
src/USER-YAFF/pair_mm3_switch3_coulgauss_long.cpp
View File

@ -106,6 +106,7 @@ void PairMM3Switch3CoulGaussLong::compute(int eflag, int vflag)
firstneigh = list->firstneigh;
// loop over neighbors of my atoms
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
@ -170,6 +171,7 @@ void PairMM3Switch3CoulGaussLong::compute(int eflag, int vflag)
expn2 = 0.0;
erfc2 = 0.0;
forcecoul2 = 0.0;
prefactor2 = 0.0;
} else {
rrij = lj2[itype][jtype]*r;
expn2 = exp(-rrij*rrij);
@ -263,19 +265,18 @@ void PairMM3Switch3CoulGaussLong::allocate()
void PairMM3Switch3CoulGaussLong::settings(int narg, char **arg)
{
if (narg < 2 || narg > 3) error->all(FLERR,"Illegal pair_style command");
if (narg < 2 || narg > 3) error->all(FLERR,"Illegal pair_style command");
cut_lj_global = utils::numeric(FLERR,arg[0],false,lmp);
if (narg == 2) {
cut_coul = cut_lj_global;
truncw = utils::numeric(FLERR,arg[1],false,lmp);
}
else {
} else {
cut_coul = utils::numeric(FLERR,arg[1],false,lmp);
truncw = utils::numeric(FLERR,arg[2],false,lmp);
}
if (truncw>0.0) truncwi = 1.0/truncw;
else truncwi = 0.0;
// reset cutoffs that have been explicitly set
if (allocated) {
@ -333,6 +334,9 @@ void PairMM3Switch3CoulGaussLong::init_style()
cut_coulsq = cut_coul * cut_coul;
if (truncw>0.0) truncwi = 1.0/truncw;
else truncwi = 0.0;
// insure use of KSpace long-range solver, set g_ewald
if (force->kspace == nullptr)
@ -375,8 +379,7 @@ double PairMM3Switch3CoulGaussLong::init_one(int i, int j)
double r6inv = r2inv*r2inv*r2inv;
double expb = lj3[i][j]*exp(-lj1[i][j]*r);
offset[i][j] = expb-lj4[i][j]*r6inv;
}
else {offset[i][j] = 0.0;}
} else {offset[i][j] = 0.0;}
} else offset[i][j] = 0.0;
cut_ljsq[j][i] = cut_ljsq[i][j];
@ -426,8 +429,7 @@ double PairMM3Switch3CoulGaussLong::init_one(int i, int j)
double t64 = cg * (0.6388888889e3 * ((-t3 + (0.7e1 / 0.36e2 * cg5 - t5) * t1 - 0.2e1 / 0.3e1 * t8 * (cg5 - cg1 / 0.4e1) * cg3 + cg5 * t14) * t20 + t3 + (cg5 / 0.12e2 + t5) * t1 + (cg5 + cg1 / 0.3e1) * cg1 * cg3 / 0.2e1 + t30 * cg5) * t2 * t36 * t39 - 0.225e1 * (0.2e1 * t43 * t44 - 0.2e1 * t43 * t47 + cg5 * (cg5 - 0.2e1 * cg1)) * t54 * t1 / cg1 / t8 * t39);
etail_ij = 2.0*MY_PI*all[0]*all[1]*t64;
ptail_ij = 2.0*MY_PI*all[0]*all[1]*t64;
}
else {
} else {
double t2 = pow(cg3, 0.2e1);
double t3 = t2 * t2;
double t7 = 0.12e2 / cg3 * cg1;

151
src/library.cpp
View File

@ -1615,7 +1615,7 @@ lists the available options.
* \return pointer (cast to ``void *``) to the location of the
* requested data or ``NULL`` if not found. */
void *lammps_extract_compute(void *handle, char *id, int style, int type)
void *lammps_extract_compute(void *handle, const char *id, int style, int type)
{
LAMMPS *lmp = (LAMMPS *) handle;
@ -1801,7 +1801,7 @@ The following table lists the available options.
* \return pointer (cast to ``void *``) to the location of the
* requested data or ``NULL`` if not found. */
void *lammps_extract_fix(void *handle, char *id, int style, int type,
void *lammps_extract_fix(void *handle, const char *id, int style, int type,
int nrow, int ncol)
{
LAMMPS *lmp = (LAMMPS *) handle;
@ -3822,8 +3822,8 @@ X(1),Y(1),Z(1),X(2),Y(2),Z(2),...,X(N),Y(N),Z(N).
* \return number of atoms created on success;
-1 on failure (no box, no atom IDs, etc.) */
int lammps_create_atoms(void *handle, int n, tagint *id, int *type,
double *x, double *v, imageint *image,
int lammps_create_atoms(void *handle, int n, const tagint *id, const int *type,
const double *x, const double *v, const imageint *image,
int bexpand)
{
LAMMPS *lmp = (LAMMPS *) handle;
@ -3859,13 +3859,17 @@ int lammps_create_atoms(void *handle, int n, tagint *id, int *type,
int nlocal_prev = nlocal;
double xdata[3];
imageint idata, *img;
for (int i = 0; i < n; i++) {
xdata[0] = x[3*i];
xdata[1] = x[3*i+1];
xdata[2] = x[3*i+2];
imageint * img = image ? image + i : nullptr;
tagint tag = id ? id[i] : 0;
if (image) {
idata = image[i];
img = &idata;
} else img = nullptr;
const tagint tag = id ? id[i] : 0;
// create atom only on MPI rank that would own it
@ -3917,40 +3921,43 @@ int lammps_create_atoms(void *handle, int n, tagint *id, int *type,
// Library functions for accessing neighbor lists
// ----------------------------------------------------------------------
/** Find neighbor list index of pair style neighbor list
/** Find index of a neighbor list requested by a pair style
*
* Try finding pair instance that matches style. If exact is set, the pair must
* match style exactly. If exact is 0, style must only be contained. If pair is
* of style pair/hybrid, style is instead matched the nsub-th hybrid sub-style.
* This function determines which of the available neighbor lists for
* pair styles matches the given conditions. It first matches the style
* name. If exact is 1 the name must match exactly, if exact is 0, a
* regular expression or sub-string match is done. If the pair style is
* hybrid or hybrid/overlay the style is matched against the sub styles
* instead.
* If a the same pair style is used multiple times as a sub-style, the
* nsub argument must be > 0 and represents the nth instance of the sub-style
* (same as for the pair_coeff command, for example). In that case
* nsub=0 will not produce a match and this function will return -1.
*
* Once the pair instance has been identified, multiple neighbor list requests
* may be found. Every neighbor list is uniquely identified by its request
* index. Thus, providing this request index ensures that the correct neighbor
* list index is returned.
* The final condition to be checked is the request ID (reqid). This
* will normally be 0, but some pair styles request multiple neighbor
* lists and set the request ID to a value > 0.
*
* \param handle pointer to a previously created LAMMPS instance cast to ``void *``.
* \param style String used to search for pair style instance
* \param exact Flag to control whether style should match exactly or only
* must be contained in pair style name
* \param nsub match nsub-th hybrid sub-style
* \param request request index that specifies which neighbor list should be
* returned, in case there are multiple neighbor lists requests
* for the found pair style
* a regular expression / sub-string match is applied.
* \param nsub match nsub-th hybrid sub-style instance of the same style
* \param reqid request id to identify neighbor list in case there are
* multiple requests from the same pair style instance
* \return return neighbor list index if found, otherwise -1 */
int lammps_find_pair_neighlist(void* handle, char * style, int exact, int nsub, int request) {
LAMMPS * lmp = (LAMMPS *) handle;
Pair* pair = lmp->force->pair_match(style, exact, nsub);
int lammps_find_pair_neighlist(void *handle, const char *style, int exact, int nsub, int reqid) {
LAMMPS *lmp = (LAMMPS *) handle;
Pair *pair = lmp->force->pair_match(style, exact, nsub);
if (pair != nullptr) {
// find neigh list
for (int i = 0; i < lmp->neighbor->nlist; i++) {
NeighList * list = lmp->neighbor->lists[i];
if (list->requestor_type != NeighList::PAIR || pair != list->requestor) continue;
if (list->index == request) {
return i;
}
NeighList *list = lmp->neighbor->lists[i];
if ( (list->requestor_type == NeighList::PAIR)
&& (pair == list->requestor)
&& (list->id == reqid) ) return i;
}
}
return -1;
@ -3958,74 +3965,60 @@ int lammps_find_pair_neighlist(void* handle, char * style, int exact, int nsub,
/* ---------------------------------------------------------------------- */
/** Find neighbor list index of fix neighbor list
/** Find index of a neighbor list requested by a fix
*
* The neighbor list request from a fix is identified by the fix ID and
* the request ID. The request ID is typically 0, but will be > 0 in
* case a fix has multiple neighbor list requests.
*
* \param handle pointer to a previously created LAMMPS instance cast to ``void *``.
* \param id Identifier of fix instance
* \param request request index that specifies which request should be returned,
* in case there are multiple neighbor lists for this fix
* \param reqid request id to identify neighbor list in case there are
* multiple requests from the same fix
* \return return neighbor list index if found, otherwise -1 */
int lammps_find_fix_neighlist(void* handle, char *id, int request) {
LAMMPS * lmp = (LAMMPS *) handle;
Fix* fix = nullptr;
const int nfix = lmp->modify->nfix;
int lammps_find_fix_neighlist(void *handle, const char *id, int reqid) {
LAMMPS *lmp = (LAMMPS *) handle;
const int ifix = lmp->modify->find_fix(id);
if (ifix < 0) return -1;
// find fix with name
for (int ifix = 0; ifix < nfix; ifix++) {
if (strcmp(lmp->modify->fix[ifix]->id, id) == 0) {
fix = lmp->modify->fix[ifix];
break;
}
}
if (fix != nullptr) {
// find neigh list
for (int i = 0; i < lmp->neighbor->nlist; i++) {
NeighList * list = lmp->neighbor->lists[i];
if (list->requestor_type != NeighList::FIX || fix != list->requestor) continue;
if (list->index == request) {
return i;
}
}
Fix *fix = lmp->modify->fix[ifix];
// find neigh list
for (int i = 0; i < lmp->neighbor->nlist; i++) {
NeighList *list = lmp->neighbor->lists[i];
if ( (list->requestor_type == NeighList::FIX)
&& (fix == list->requestor)
&& (list->id == reqid) ) return i;
}
return -1;
}
/* ---------------------------------------------------------------------- */
/** Find neighbor list index of compute neighbor list
/** Find index of a neighbor list requested by a compute
*
* The neighbor list request from a compute is identified by the compute
* ID and the request ID. The request ID is typically 0, but will be
* > 0 in case a compute has multiple neighbor list requests.
*
* \param handle pointer to a previously created LAMMPS instance cast to ``void *``.
* \param id Identifier of fix instance
* \param request request index that specifies which request should be returned,
* in case there are multiple neighbor lists for this fix
* \param id Identifier of compute instance
* \param reqid request id to identify neighbor list in case there are
* multiple requests from the same compute
* \return return neighbor list index if found, otherwise -1 */
int lammps_find_compute_neighlist(void* handle, char *id, int request) {
LAMMPS * lmp = (LAMMPS *) handle;
Compute* compute = nullptr;
const int ncompute = lmp->modify->ncompute;
int lammps_find_compute_neighlist(void* handle, const char *id, int reqid) {
LAMMPS *lmp = (LAMMPS *) handle;
const int icompute = lmp->modify->find_compute(id);
if (icompute < 0) return -1;
// find compute with name
for (int icompute = 0; icompute < ncompute; icompute++) {
if (strcmp(lmp->modify->compute[icompute]->id, id) == 0) {
compute = lmp->modify->compute[icompute];
break;
}
}
if (compute != nullptr) {
// find neigh list
for (int i = 0; i < lmp->neighbor->nlist; i++) {
NeighList * list = lmp->neighbor->lists[i];
if (list->requestor_type != NeighList::COMPUTE || compute != list->requestor) continue;
if (list->index == request) {
return i;
}
}
Compute *compute = lmp->modify->compute[icompute];
// find neigh list
for (int i = 0; i < lmp->neighbor->nlist; i++) {
NeighList * list = lmp->neighbor->lists[i];
if ( (list->requestor_type == NeighList::COMPUTE)
&& (compute == list->requestor)
&& (list->id == reqid) ) return i;
}
return -1;
}

18
src/library.h
View File

@ -138,8 +138,8 @@ void *lammps_extract_atom(void *handle, const char *name);
* Library functions to access data from computes, fixes, variables in LAMMPS
* ---------------------------------------------------------------------- */
void *lammps_extract_compute(void *handle, char *id, int, int);
void *lammps_extract_fix(void *handle, char *, int, int, int, int);
void *lammps_extract_compute(void *handle, const char *, int, int);
void *lammps_extract_fix(void *handle, const char *, int, int, int, int);
void *lammps_extract_variable(void *handle, const char *, const char *);
int lammps_set_variable(void *, char *, char *);
@ -160,20 +160,20 @@ void lammps_scatter(void *handle, char *name, int type, int count, void *data);
void lammps_scatter_subset(void *handle, char *name, int type, int count, int ndata, int *ids, void *data);
#if !defined(LAMMPS_BIGBIG)
int lammps_create_atoms(void *handle, int n, int *id, int *type,
double *x, double *v, int *image, int bexpand);
int lammps_create_atoms(void *handle, int n, const int *id, const int *type,
const double *x, const double *v, const int *image, int bexpand);
#else
int lammps_create_atoms(void *handle, int n, int64_t *id, int *type,
double *x, double *v, int64_t* image, int bexpand);
int lammps_create_atoms(void *handle, int n, const int64_t *id, const int *type,
const double *x, const double *v, const int64_t* image, int bexpand);
#endif
/* ----------------------------------------------------------------------
* Library functions for accessing neighbor lists
* ---------------------------------------------------------------------- */
int lammps_find_pair_neighlist(void *handle, char *style, int exact, int nsub, int request);
int lammps_find_fix_neighlist(void *handle, char *id, int request);
int lammps_find_compute_neighlist(void *handle, char *id, int request);
int lammps_find_pair_neighlist(void *handle, const char *style, int exact, int nsub, int request);
int lammps_find_fix_neighlist(void *handle, const char *id, int request);
int lammps_find_compute_neighlist(void *handle, const char *id, int request);
int lammps_neighlist_num_elements(void *handle, int idx);
void lammps_neighlist_element_neighbors(void *handle, int idx, int element, int *iatom, int *numneigh, int **neighbors);

1
src/neigh_list.cpp
View File

@ -143,6 +143,7 @@ void NeighList::post_constructor(NeighRequest *nq)
respamiddle = nq->respamiddle;
respainner = nq->respainner;
copy = nq->copy;
id = nq->id;
if (nq->copy) {
listcopy = neighbor->lists[nq->copylist];

1
src/neigh_list.h
View File

@ -44,6 +44,7 @@ class NeighList : protected Pointers {
int copy; // 1 if this list is copied from another list
int kk2cpu; // 1 if this list is copied from Kokkos to CPU
int copymode; // 1 if this is a Kokkos on-device copy
int id; // copied from neighbor list request
// data structs to store neighbor pairs I,J and associated values

11
src/random_mars.cpp
View File

@ -136,13 +136,16 @@ double RanMars::gaussian(double mu, double sigma)
double RanMars::rayleigh(double sigma)
{
double first,v1;
double v1;
if (sigma <= 0) error->all(FLERR,"Invalid Rayleigh parameter");
if (sigma <= 0.0) error->all(FLERR,"Invalid Rayleigh parameter");
v1 = uniform();
first = sigma*sqrt(-2.0*log(v1));
return first;
// avoid a floating point exception due to log(0.0)
// and just return a very big number
if (v1 == 0.0) return 1.0e300;
return sigma*sqrt(-2.0*log(v1));
}
/* ----------------------------------------------------------------------

12
src/rerun.cpp
View File

@ -51,6 +51,7 @@ void Rerun::command(int narg, char **arg)
if (strcmp(arg[iarg],"start") == 0) break;
if (strcmp(arg[iarg],"stop") == 0) break;
if (strcmp(arg[iarg],"dump") == 0) break;
if (strcmp(arg[iarg],"post") == 0) break;
iarg++;
}
int nfile = iarg;
@ -65,6 +66,7 @@ void Rerun::command(int narg, char **arg)
int nskip = 1;
int startflag = 0;
int stopflag = 0;
int postflag = 0;
bigint start = -1;
bigint stop = -1;
@ -101,6 +103,14 @@ void Rerun::command(int narg, char **arg)
stop = utils::bnumeric(FLERR,arg[iarg+1],false,lmp);
if (stop < 0) error->all(FLERR,"Illegal rerun command");
iarg += 2;
} else if (strcmp(arg[iarg],"post") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal rerun command");
if (strcmp(arg[iarg+1],"yes") == 0) {
postflag = 1;
} else if (strcmp(arg[iarg+1],"no") == 0) {
postflag = 0;
} else error->all(FLERR,"Illegal rerun command");
iarg += 2;
} else if (strcmp(arg[iarg],"dump") == 0) {
break;
} else error->all(FLERR,"Illegal rerun command");
@ -182,7 +192,7 @@ void Rerun::command(int narg, char **arg)
update->nsteps = ndump;
Finish finish(lmp);
finish.end(1);
finish.end(postflag);
update->whichflag = 0;
update->firststep = update->laststep = 0;

4
src/tokenizer.h
View File

@ -41,7 +41,7 @@ public:
Tokenizer& operator=(Tokenizer&&) = default;
void reset();
void skip(int n);
void skip(int n=1);
bool has_next() const;
bool contains(const std::string &str) const;
std::string next();
@ -104,7 +104,7 @@ public:
bool has_next() const;
bool contains(const std::string &value) const;
void skip(int ntokens);
void skip(int ntokens=1);
size_t count();
};

2
src/velocity.cpp
View File

@ -108,7 +108,7 @@ void Velocity::command(int narg, char **arg)
int initcomm = 0;
if (style == ZERO && rfix >= 0 &&
utils::strmatch(modify->fix[rfix]->style,"^rigid/small")) initcomm = 1;
utils::strmatch(modify->fix[rfix]->style,"^rigid.*/small.*")) initcomm = 1;
if ((style == CREATE || style == SET) && temperature &&
strcmp(temperature->style,"temp/cs") == 0) initcomm = 1;

2
tools/singularity/ubuntu18.04_amd_rocm.def
View File

@ -94,7 +94,7 @@ From: ubuntu:18.04
###########################################################################
export PATH=$PATH:/opt/rocm/bin:/opt/rocm/profiler/bin:/opt/rocm/opencl/bin/x86_64
git clone -b rocm-3.7.x https://github.com/ROCmSoftwarePlatform/hipCUB.git
git clone -b rocm-4.1.x https://github.com/ROCmSoftwarePlatform/hipCUB.git
mkdir hipCUB/build
cd hipCUB/build
CXX=hipcc cmake -D BUILD_TEST=off ..

4
tools/singularity/ubuntu20.04_amd_rocm.def
View File

@ -2,7 +2,7 @@ BootStrap: docker
From: ubuntu:20.04
%environment
export PATH=/usr/lib/ccache:/usr/local/cuda-11.0/bin:${PATH}:/opt/rocm/bin:/opt/rocm/profiler/bin:/opt/rocm/opencl/bin/x86_64
export PATH=/usr/lib/ccache:${PATH}:/opt/rocm/bin:/opt/rocm/profiler/bin:/opt/rocm/opencl/bin/x86_64
%post
export DEBIAN_FRONTEND=noninteractive
apt-get update
@ -90,7 +90,7 @@ From: ubuntu:20.04
###########################################################################
export PATH=$PATH:/opt/rocm/bin:/opt/rocm/profiler/bin:/opt/rocm/opencl/bin/x86_64
git clone -b rocm-3.7.x https://github.com/ROCmSoftwarePlatform/hipCUB.git
git clone -b rocm-4.1.x https://github.com/ROCmSoftwarePlatform/hipCUB.git
mkdir hipCUB/build
cd hipCUB/build
CXX=hipcc cmake -D BUILD_TEST=off ..

1
unittest/c-library/CMakeLists.txt
View File

@ -11,6 +11,7 @@ add_executable(test_library_properties test_library_properties.cpp test_main.cpp
target_link_libraries(test_library_properties PRIVATE lammps GTest::GTest GTest::GMock)
target_compile_definitions(test_library_properties PRIVATE -DTEST_INPUT_FOLDER=${CMAKE_CURRENT_SOURCE_DIR})
add_test(LibraryProperties test_library_properties)
set_tests_properties(LibraryProperties PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR}")
set(TEST_CONFIG_DEFS "-DTEST_INPUT_FOLDER=${CMAKE_CURRENT_SOURCE_DIR};-DLAMMPS_${LAMMPS_SIZES}")
set(PKG_COUNT 0)

117
unittest/c-library/test_library_properties.cpp
View File

@ -160,7 +160,9 @@ TEST_F(LibraryProperties, box)
boxlo[0] = -6.1;
boxhi[1] = 7.3;
xy = 0.1;
if (!verbose) ::testing::internal::CaptureStdout();
lammps_reset_box(lmp, boxlo, boxhi, xy, yz, xz);
if (!verbose) ::testing::internal::GetCapturedStdout();
lammps_extract_box(lmp, boxlo, boxhi, &xy, &yz, &xz, pflags, &boxflag);
EXPECT_DOUBLE_EQ(boxlo[0], -6.1);
EXPECT_DOUBLE_EQ(boxlo[1], -7.692866);
@ -313,6 +315,121 @@ TEST_F(LibraryProperties, global)
EXPECT_DOUBLE_EQ((*d_ptr), 0.1);
};
TEST_F(LibraryProperties, neighlist)
{
if (!lammps_has_style(lmp, "pair", "sw")) GTEST_SKIP();
const char sysinit[] = "boundary f f f\n"
"units real\n"
"region box block -5 5 -5 5 -5 5\n"
"create_box 2 box\n"
"mass 1 1.0\n"
"mass 2 1.0\n"
"pair_style hybrid/overlay lj/cut 4.0 lj/cut 4.0 morse 4.0 sw\n"
"pair_coeff * * sw Si.sw Si NULL\n"
"pair_coeff 1 2 morse 0.2 2.0 2.0\n"
"pair_coeff 2 2 lj/cut 1 0.1 2.0\n"
"pair_coeff * * lj/cut 2 0.01 2.0\n"
"compute dist all pair/local dist\n"
"fix dist all ave/histo 1 1 1 0.0 3.0 4 c_dist mode vector\n"
"thermo_style custom f_dist[*]";
const double pos[] = {0.0, 0.0, 0.0, -1.1, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, -1.1,
0.0, 0.0, 1.0, 0.0, 0.0, 0.0, -1.1, 0.0, 0.0, 1.0};
const tagint ids[] = {1, 2, 3, 4, 5, 6, 7};
const int types[] = {1, 1, 1, 1, 2, 2, 2};
const int numatoms = sizeof(ids) / sizeof(tagint);
if (!verbose) ::testing::internal::CaptureStdout();
lammps_commands_string(lmp, sysinit);
lammps_create_atoms(lmp, numatoms, ids, types, pos, nullptr, nullptr, 0);
lammps_command(lmp, "run 0 post no");
if (!verbose) ::testing::internal::GetCapturedStdout();
int nhisto =
*(double *)lammps_extract_fix(lmp, "dist", LMP_STYLE_GLOBAL, LMP_TYPE_VECTOR, 0, 0);
int nskip = *(double *)lammps_extract_fix(lmp, "dist", LMP_STYLE_GLOBAL, LMP_TYPE_VECTOR, 1, 0);
double minval =
*(double *)lammps_extract_fix(lmp, "dist", LMP_STYLE_GLOBAL, LMP_TYPE_VECTOR, 2, 0);
double maxval =
*(double *)lammps_extract_fix(lmp, "dist", LMP_STYLE_GLOBAL, LMP_TYPE_VECTOR, 3, 0);
// 21 pair distances counted, none skipped, smallest 1.0, largest 2.1
EXPECT_EQ(nhisto, 21);
EXPECT_EQ(nskip, 0);
EXPECT_DOUBLE_EQ(minval, 1.0);
EXPECT_DOUBLE_EQ(maxval, 2.1);
const int nlocal = lammps_extract_setting(lmp, "nlocal");
EXPECT_EQ(nlocal, numatoms);
EXPECT_NE(lammps_find_pair_neighlist(lmp, "sw", 1, 0, 0), -1);
EXPECT_NE(lammps_find_pair_neighlist(lmp, "morse", 1, 0, 0), -1);
EXPECT_NE(lammps_find_pair_neighlist(lmp, "lj/cut", 1, 1, 0), -1);
EXPECT_NE(lammps_find_pair_neighlist(lmp, "lj/cut", 1, 2, 0), -1);
EXPECT_EQ(lammps_find_pair_neighlist(lmp, "lj/cut", 1, 0, 0), -1);
EXPECT_EQ(lammps_find_pair_neighlist(lmp, "hybrid/overlay", 1, 0, 0), -1);
EXPECT_NE(lammps_find_compute_neighlist(lmp, "dist", 0), -1);
EXPECT_EQ(lammps_find_fix_neighlist(lmp, "dist", 0), -1);
EXPECT_EQ(lammps_find_compute_neighlist(lmp, "xxx", 0), -1);
// full neighbor list for 4 type 1 atoms
// all have 3 type 1 atom neighbors
int idx = lammps_find_pair_neighlist(lmp, "sw", 1, 0, 0);
int num = lammps_neighlist_num_elements(lmp, idx);
EXPECT_EQ(num, 4);
int iatom, inum, *neighbors;
for (int i = 0; i < num; ++i) {
lammps_neighlist_element_neighbors(lmp, idx, i, &iatom, &inum, &neighbors);
EXPECT_EQ(iatom, i);
EXPECT_EQ(inum, 3);
EXPECT_NE(neighbors, nullptr);
}
// half neighbor list for all pairs between type 1 and type 2
// 4 type 1 atoms with 3 type 2 neighbors and 3 type 2 atoms without neighbors
idx = lammps_find_pair_neighlist(lmp, "morse", 0, 0, 0);
num = lammps_neighlist_num_elements(lmp, idx);
EXPECT_EQ(num, nlocal);
for (int i = 0; i < num; ++i) {
lammps_neighlist_element_neighbors(lmp, idx, i, &iatom, &inum, &neighbors);
if (i < 4)
EXPECT_EQ(inum, 3);
else
EXPECT_EQ(inum, 0);
EXPECT_NE(neighbors, nullptr);
}
// half neighbor list between type 2 atoms only
// 3 pairs with 2, 1, 0 neighbors
idx = lammps_find_pair_neighlist(lmp, "lj/cut", 1, 1, 0);
num = lammps_neighlist_num_elements(lmp, idx);
EXPECT_EQ(num, 3);
for (int i = 0; i < num; ++i) {
lammps_neighlist_element_neighbors(lmp, idx, i, &iatom, &inum, &neighbors);
EXPECT_EQ(inum, 2 - i);
EXPECT_NE(neighbors, nullptr);
}
// half neighbor list between all pairs. same as simple lj/cut case
idx = lammps_find_pair_neighlist(lmp, "lj/cut", 1, 2, 0);
num = lammps_neighlist_num_elements(lmp, idx);
EXPECT_EQ(num, nlocal);
for (int i = 0; i < num; ++i) {
lammps_neighlist_element_neighbors(lmp, idx, i, &iatom, &inum, &neighbors);
EXPECT_EQ(inum, nlocal - 1 - i);
EXPECT_NE(neighbors, nullptr);
}
// the compute has a half neighbor list
idx = lammps_find_compute_neighlist(lmp, "dist", 0);
num = lammps_neighlist_num_elements(lmp, idx);
EXPECT_EQ(num, nlocal);
for (int i = 0; i < num; ++i) {
lammps_neighlist_element_neighbors(lmp, idx, i, &iatom, &inum, &neighbors);
EXPECT_EQ(inum, nlocal - 1 - i);
EXPECT_NE(neighbors, nullptr);
}
};
class AtomProperties : public ::testing::Test {
protected:
void *lmp;

4
unittest/commands/test_simple_commands.cpp
View File

@ -328,11 +328,7 @@ TEST_F(SimpleCommandsTest, Units)
#if defined(LMP_PLUGIN)
TEST_F(SimpleCommandsTest, Plugin)
{
#if defined(__APPLE__)
std::string loadfmt("plugin load {}plugin.dylib");
#else
std::string loadfmt("plugin load {}plugin.so");
#endif
::testing::internal::CaptureStdout();
lmp->input->one(fmt::format(loadfmt, "hello"));
auto text = ::testing::internal::GetCapturedStdout();

44
unittest/force-styles/test_angle_style.cpp
View File

@ -38,7 +38,6 @@
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <mpi.h>
#include <map>
@ -235,41 +234,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
YamlWriter writer(outfile);
// lammps_version
writer.emit("lammps_version", lmp->version);
// date_generated
std::time_t now = time(NULL);
block = ctime(&now);
block = block.substr(0, block.find("\n") - 1);
writer.emit("date_generated", block);
// epsilon
writer.emit("epsilon", config.epsilon);
// prerequisites
block.clear();
for (auto &prerequisite : config.prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer.emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : config.pre_commands) {
block += command + "\n";
}
writer.emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : config.post_commands) {
block += command + "\n";
}
writer.emit_block("post_commands", block);
// input_file
writer.emit("input_file", config.input_file);
// write yaml header
write_yaml_header(&writer, &test_config, lmp->version);
// angle_style
writer.emit("angle_style", config.angle_style);
@ -340,6 +306,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
TEST(AngleStyle, plain)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"AngleStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -566,6 +534,8 @@ TEST(AngleStyle, plain)
TEST(AngleStyle, omp)
{
if (!LAMMPS::is_installed_pkg("USER-OMP")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"AngleStyle", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "omp", "4", "-sf", "omp"};
@ -739,6 +709,8 @@ TEST(AngleStyle, omp)
TEST(AngleStyle, single)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"AngleStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;

46
unittest/force-styles/test_bond_style.cpp
View File

@ -38,7 +38,6 @@
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <mpi.h>
#include <map>
@ -235,41 +234,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
YamlWriter writer(outfile);
// lammps_version
writer.emit("lammps_version", lmp->version);
// date_generated
std::time_t now = time(NULL);
block = ctime(&now);
block = block.substr(0, block.find("\n") - 1);
writer.emit("date_generated", block);
// epsilon
writer.emit("epsilon", config.epsilon);
// prerequisites
block.clear();
for (auto &prerequisite : config.prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer.emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : config.pre_commands) {
block += command + "\n";
}
writer.emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : config.post_commands) {
block += command + "\n";
}
writer.emit_block("post_commands", block);
// input_file
writer.emit("input_file", config.input_file);
// write yaml header
write_yaml_header(&writer, &test_config, lmp->version);
// bond_style
writer.emit("bond_style", config.bond_style);
@ -340,6 +306,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
TEST(BondStyle, plain)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"BondStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -566,6 +534,8 @@ TEST(BondStyle, plain)
TEST(BondStyle, omp)
{
if (!LAMMPS::is_installed_pkg("USER-OMP")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"BondStyle", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "omp", "4", "-sf", "omp"};
@ -738,6 +708,8 @@ TEST(BondStyle, omp)
TEST(BondStyle, single)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"BondStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -994,6 +966,8 @@ TEST(BondStyle, single)
TEST(BondStyle, extract)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"BondStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;

3
unittest/force-styles/test_config.h
View File

@ -14,6 +14,7 @@
#ifndef TEST_CONFIG_H
#define TEST_CONFIG_H
#include <set>
#include <string>
#include <utility>
#include <vector>
@ -32,6 +33,7 @@ public:
std::string date_generated;
std::string basename;
double epsilon;
std::set<std::string> skip_tests;
std::vector<std::pair<std::string, std::string>> prerequisites;
std::vector<std::string> pre_commands;
std::vector<std::string> post_commands;
@ -74,6 +76,7 @@ public:
init_vdwl(0), run_vdwl(0), init_coul(0), run_coul(0), init_stress({0, 0, 0, 0, 0, 0}),
run_stress({0, 0, 0, 0, 0, 0}), global_scalar(0)
{
skip_tests.clear();
prerequisites.clear();
pre_commands.clear();
post_commands.clear();

11
unittest/force-styles/test_config_reader.cpp
View File

@ -15,6 +15,7 @@
#include "test_config.h"
#include "yaml.h"
#include "yaml_reader.h"
#include "utils.h"
#include <cstdlib>
#include <cstring>
@ -25,11 +26,14 @@
#include <utility>
#include <vector>
using LAMMPS_NS::utils::split_words;
TestConfigReader::TestConfigReader(TestConfig &config) : YamlReader(), config(config)
{
consumers["lammps_version"] = &TestConfigReader::lammps_version;
consumers["date_generated"] = &TestConfigReader::date_generated;
consumers["epsilon"] = &TestConfigReader::epsilon;
consumers["skip_tests"] = &TestConfigReader::skip_tests;
consumers["prerequisites"] = &TestConfigReader::prerequisites;
consumers["pre_commands"] = &TestConfigReader::pre_commands;
consumers["post_commands"] = &TestConfigReader::post_commands;
@ -66,6 +70,13 @@ TestConfigReader::TestConfigReader(TestConfig &config) : YamlReader(), config(co
consumers["equilibrium"] = &TestConfigReader::equilibrium;
}
void TestConfigReader::skip_tests(const yaml_event_t &event)
{
config.skip_tests.clear();
for (auto &word : split_words((char *)event.data.scalar.value))
config.skip_tests.insert(word);
}
void TestConfigReader::prerequisites(const yaml_event_t &event)
{
config.prerequisites.clear();

1
unittest/force-styles/test_config_reader.h
View File

@ -23,6 +23,7 @@ class TestConfigReader : public YamlReader<TestConfigReader> {
public:
TestConfigReader(TestConfig &config);
void skip_tests(const yaml_event_t &event);
void prerequisites(const yaml_event_t &event);
void pre_commands(const yaml_event_t &event);
void post_commands(const yaml_event_t &event);

42
unittest/force-styles/test_dihedral_style.cpp
View File

@ -38,7 +38,6 @@
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <mpi.h>
#include <map>
@ -244,41 +243,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
YamlWriter writer(outfile);
// lammps_version
writer.emit("lammps_version", lmp->version);
// date_generated
std::time_t now = time(NULL);
block = ctime(&now);
block = block.substr(0, block.find("\n") - 1);
writer.emit("date_generated", block);
// epsilon
writer.emit("epsilon", config.epsilon);
// prerequisites
block.clear();
for (auto &prerequisite : config.prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer.emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : config.pre_commands) {
block += command + "\n";
}
writer.emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : config.post_commands) {
block += command + "\n";
}
writer.emit_block("post_commands", block);
// input_file
writer.emit("input_file", config.input_file);
// write yaml header
write_yaml_header(&writer, &test_config, lmp->version);
// dihedral_style
writer.emit("dihedral_style", config.dihedral_style);
@ -343,6 +309,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
TEST(DihedralStyle, plain)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"DihedralStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -569,6 +537,8 @@ TEST(DihedralStyle, plain)
TEST(DihedralStyle, omp)
{
if (!LAMMPS::is_installed_pkg("USER-OMP")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"DihedralStyle", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "omp", "4", "-sf", "omp"};

42
unittest/force-styles/test_fix_timestep.cpp
View File

@ -42,7 +42,6 @@
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <mpi.h>
#include <map>
@ -193,41 +192,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
std::string block("");
YamlWriter writer(outfile);
// lammps_version
writer.emit("lammps_version", lmp->version);
// date_generated
std::time_t now = time(NULL);
block = ctime(&now);
block = block.substr(0, block.find("\n") - 1);
writer.emit("date_generated", block);
// epsilon
writer.emit("epsilon", config.epsilon);
// prerequisites
block.clear();
for (auto &prerequisite : config.prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer.emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : config.pre_commands) {
block += command + "\n";
}
writer.emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : config.post_commands) {
block += command + "\n";
}
writer.emit_block("post_commands", block);
// input_file
writer.emit("input_file", config.input_file);
// write yaml header
write_yaml_header(&writer, &test_config, lmp->version);
// natoms
writer.emit("natoms", natoms);
@ -288,6 +254,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
TEST(FixTimestep, plain)
{
if (!LAMMPS::is_installed_pkg("MOLECULE")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"FixTimestep", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -730,6 +698,8 @@ TEST(FixTimestep, omp)
{
if (!LAMMPS::is_installed_pkg("USER-OMP")) GTEST_SKIP();
if (!LAMMPS::is_installed_pkg("MOLECULE")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"FixTimestep", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "omp", "4", "-sf", "omp"};

42
unittest/force-styles/test_improper_style.cpp
View File

@ -38,7 +38,6 @@
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <mpi.h>
#include <map>
@ -235,41 +234,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
YamlWriter writer(outfile);
// lammps_version
writer.emit("lammps_version", lmp->version);
// date_generated
std::time_t now = time(NULL);
block = ctime(&now);
block = block.substr(0, block.find("\n") - 1);
writer.emit("date_generated", block);
// epsilon
writer.emit("epsilon", config.epsilon);
// prerequisites
block.clear();
for (auto &prerequisite : config.prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer.emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : config.pre_commands) {
block += command + "\n";
}
writer.emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : config.post_commands) {
block += command + "\n";
}
writer.emit_block("post_commands", block);
// input_file
writer.emit("input_file", config.input_file);
// write yaml header
write_yaml_header(&writer, &test_config, lmp->version);
// improper_style
writer.emit("improper_style", config.improper_style);
@ -334,6 +300,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
TEST(ImproperStyle, plain)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"ImproperStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -560,6 +528,8 @@ TEST(ImproperStyle, plain)
TEST(ImproperStyle, omp)
{
if (!LAMMPS::is_installed_pkg("USER-OMP")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"ImproperStyle", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "omp", "4", "-sf", "omp"};

52
unittest/force-styles/test_main.cpp
View File

@ -16,16 +16,19 @@
#include "test_config.h"
#include "test_config_reader.h"
#include "utils.h"
#include "yaml_writer.h"
#include "gmock/gmock.h"
#include "gtest/gtest.h"
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <iostream>
#include <mpi.h>
#include <vector>
using LAMMPS_NS::utils::split_words;
using LAMMPS_NS::utils::trim;
// common read_yaml_file function
bool read_yaml_file(const char *infile, TestConfig &config)
@ -37,6 +40,55 @@ bool read_yaml_file(const char *infile, TestConfig &config)
return true;
}
// write out common header items for yaml files
void write_yaml_header(YamlWriter *writer,
TestConfig *cfg,
const char *version)
{
// lammps_version
writer->emit("lammps_version", version);
// date_generated
std::time_t now = time(NULL);
std::string block = trim(ctime(&now));
writer->emit("date_generated", block);
// epsilon
writer->emit("epsilon", cfg->epsilon);
// skip tests
block.clear();
for (auto &skip : cfg->skip_tests) {
if (block.empty()) block = skip;
else block += " " + skip;
}
writer->emit("skip_tests", block);
// prerequisites
block.clear();
for (auto &prerequisite : cfg->prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer->emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : cfg->pre_commands) {
block += command + "\n";
}
writer->emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : cfg->post_commands) {
block += command + "\n";
}
writer->emit_block("post_commands", block);
// input_file
writer->emit("input_file", cfg->input_file);
}
// need to be defined in unit test body
extern void generate_yaml_file(const char *, const TestConfig &);

9
unittest/force-styles/test_main.h
View File

@ -22,6 +22,10 @@ extern bool print_stats;
extern bool verbose;
extern std::string INPUT_FOLDER;
// convenience method to write out common entries
void write_yaml_header(class YamlWriter *writer, TestConfig *cfg,
const char *version);
#define EXPECT_FP_LE_WITH_EPS(val1, val2, eps) \
do { \
const double diff = fabs(val1 - val2); \
@ -31,10 +35,11 @@ extern std::string INPUT_FOLDER;
EXPECT_PRED_FORMAT2(::testing::DoubleLE, err, eps); \
} while (0);
#endif
#if defined _WIN32
static const char PATH_SEP = '\\';
#else
static const char PATH_SEP = '/';
#endif
#endif

206
unittest/force-styles/test_pair_style.cpp
View File

@ -40,7 +40,6 @@
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <mpi.h>
#include <map>
@ -125,6 +124,13 @@ LAMMPS *init_lammps(int argc, char **argv, const TestConfig &cfg, const bool new
command("pair_coeff " + pair_coeff);
}
// set this here explicitly to a setting different
// from the default, so we can spot YAML files for
// long-range interactions that do not include these
// settings. they will fail after restart or read data.
command("pair_modify table 0");
command("pair_modify table/disp 0");
for (auto &post_command : cfg.post_commands) {
command(post_command);
}
@ -238,41 +244,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
YamlWriter writer(outfile);
// lammps_version
writer.emit("lammps_version", lmp->version);
// date_generated
std::time_t now = time(NULL);
block = ctime(&now);
block = block.substr(0, block.find("\n") - 1);
writer.emit("date_generated", block);
// epsilon
writer.emit("epsilon", config.epsilon);
// prerequisites
block.clear();
for (auto &prerequisite : config.prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer.emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : config.pre_commands) {
block += command + "\n";
}
writer.emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : config.post_commands) {
block += command + "\n";
}
writer.emit_block("post_commands", block);
// input_file
writer.emit("input_file", config.input_file);
// write yaml header
write_yaml_header(&writer, &test_config, lmp->version);
// pair_style
writer.emit("pair_style", config.pair_style);
@ -343,6 +316,8 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
TEST(PairStyle, plain)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -625,9 +600,14 @@ TEST(PairStyle, plain)
TEST(PairStyle, omp)
{
if (!LAMMPS::is_installed_pkg("USER-OMP")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "omp", "4", "-sf", "omp"};
// cannot run dpd styles with more than 1 thread due to using multiple pRNGs
if (utils::strmatch(test_config.pair_style, "^dpd")) args[8] = "1";
char **argv = (char **)args;
int argc = sizeof(args) / sizeof(char *);
@ -649,14 +629,6 @@ TEST(PairStyle, omp)
EXPECT_THAT(output, StartsWith("LAMMPS ("));
EXPECT_THAT(output, HasSubstr("Loop time"));
if (utils::strmatch(test_config.pair_style, "^dpd")) {
std::cerr << "Skipping pair style " << lmp->force->pair_style << "\n";
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
GTEST_SKIP();
}
// abort if running in parallel and not all atoms are local
const int nlocal = lmp->atom->nlocal;
ASSERT_EQ(lmp->atom->natoms, nlocal);
@ -804,13 +776,136 @@ TEST(PairStyle, omp)
if (!verbose) ::testing::internal::GetCapturedStdout();
};
TEST(PairStyle, gpu)
{
if (!LAMMPS::is_installed_pkg("GPU")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args_neigh[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite", "-sf", "gpu"};
const char *args_noneigh[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite", "-sf", "gpu", "-pk", "gpu", "0", "neigh", "no"};
char **argv = (char **)args_neigh;
int argc = sizeof(args_neigh) / sizeof(char *);
// cannot use GPU neighbor list with hybrid pair style (yet)
if (test_config.pair_style.substr(0, 6) == "hybrid") {
argv = (char **)args_noneigh;
argc = sizeof(args_noneigh) / sizeof(char *);
}
::testing::internal::CaptureStdout();
LAMMPS *lmp = init_lammps(argc, argv, test_config, false);
std::string output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
if (!lmp) {
std::cerr << "One or more prerequisite styles with /gpu suffix\n"
"are not available in this LAMMPS configuration:\n";
for (auto &prerequisite : test_config.prerequisites) {
std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n";
}
GTEST_SKIP();
}
EXPECT_THAT(output, StartsWith("LAMMPS ("));
EXPECT_THAT(output, HasSubstr("Loop time"));
// abort if running in parallel and not all atoms are local
const int nlocal = lmp->atom->nlocal;
ASSERT_EQ(lmp->atom->natoms, nlocal);
// relax error for GPU package depending on precision setting
double epsilon = test_config.epsilon;
if (Info::has_accelerator_feature("GPU","precision","double"))
epsilon *= 7.5;
else if (Info::has_accelerator_feature("GPU","precision","mixed"))
epsilon *= 5.0e8;
else epsilon *= 1.0e10;
// relax test precision when using pppm and single precision FFTs, but only when also running with double precision
#if defined(FFT_SINGLE)
if (lmp->force->kspace && lmp->force->kspace->compute_flag && Info::has_accelerator_feature("GPU","precision","double"))
if (utils::strmatch(lmp->force->kspace_style, "^pppm")) epsilon *= 2.0e8;
#endif
const std::vector<coord_t> &f_ref = test_config.init_forces;
const std::vector<coord_t> &f_run = test_config.run_forces;
ErrorStats stats;
auto f = lmp->atom->f;
auto tag = lmp->atom->tag;
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_ref[tag[i]].x, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_ref[tag[i]].y, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_ref[tag[i]].z, epsilon);
}
if (print_stats) std::cerr << "init_forces stats, newton off:" << stats << std::endl;
auto pair = lmp->force->pair;
auto stress = pair->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.init_stress.xx, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.init_stress.yy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.init_stress.zz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.init_stress.xy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.init_stress.xz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.init_stress.yz, 10 * epsilon);
if (print_stats) std::cerr << "init_stress stats, newton off:" << stats << std::endl;
stats.reset();
EXPECT_FP_LE_WITH_EPS(pair->eng_vdwl, test_config.init_vdwl, epsilon);
EXPECT_FP_LE_WITH_EPS(pair->eng_coul, test_config.init_coul, epsilon);
if (print_stats) std::cerr << "init_energy stats, newton off:" << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
run_lammps(lmp);
if (!verbose) ::testing::internal::GetCapturedStdout();
f = lmp->atom->f;
tag = lmp->atom->tag;
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_run[tag[i]].x, 5 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_run[tag[i]].y, 5 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_run[tag[i]].z, 5 * epsilon);
}
if (print_stats) std::cerr << "run_forces stats, newton off:" << stats << std::endl;
stress = pair->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.run_stress.xx, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.run_stress.yy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.run_stress.zz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.run_stress.xy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.run_stress.xz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.run_stress.yz, 10 * epsilon);
if (print_stats) std::cerr << "run_stress stats, newton off:" << stats << std::endl;
stats.reset();
auto id = lmp->modify->find_compute("sum");
auto energy = lmp->modify->compute[id]->compute_scalar();
EXPECT_FP_LE_WITH_EPS(pair->eng_vdwl, test_config.run_vdwl, epsilon);
EXPECT_FP_LE_WITH_EPS(pair->eng_coul, test_config.run_coul, epsilon);
EXPECT_FP_LE_WITH_EPS((pair->eng_vdwl + pair->eng_coul), energy, epsilon);
if (print_stats) std::cerr << "run_energy stats, newton off:" << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
};
TEST(PairStyle, intel)
{
if (!LAMMPS::is_installed_pkg("USER-INTEL")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "intel", "0", "mode", "double", "omp",
"4", "lrt", "no", "-sf", "intel"};
// cannot use more than 1 thread for dpd styles due to pRNG
if (utils::strmatch(test_config.pair_style, "^dpd")) args[12] = "1";
char **argv = (char **)args;
int argc = sizeof(args) / sizeof(char *);
@ -829,23 +924,6 @@ TEST(PairStyle, intel)
GTEST_SKIP();
}
if ((test_config.pair_style == "rebo") || utils::strmatch(test_config.pair_style, "^dpd") ||
utils::strmatch(test_config.pair_style, "^tersoff.* shift ")) {
std::cerr << "Skipping pair style " << lmp->force->pair_style << "\n";
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
GTEST_SKIP();
}
if (lmp->force->kspace && utils::strmatch(lmp->force->kspace_style, "^pppm/disp/intel")) {
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
std::cerr << "Skipping kspace style pppm/disp/intel\n";
GTEST_SKIP();
}
// relax error a bit for USER-INTEL package
double epsilon = 7.5 * test_config.epsilon;
// relax test precision when using pppm and single precision FFTs
@ -946,6 +1024,8 @@ TEST(PairStyle, intel)
TEST(PairStyle, opt)
{
if (!LAMMPS::is_installed_pkg("OPT")) GTEST_SKIP();
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite", "-sf", "opt"};
char **argv = (char **)args;
@ -1051,6 +1131,8 @@ TEST(PairStyle, opt)
TEST(PairStyle, single)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
@ -1308,6 +1390,8 @@ TEST(PairStyle, single)
TEST(PairStyle, extract)
{
if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP();
const char *args[] = {"PairStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;

Some files were not shown because too many files have changed in this diff Show More