diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
index 71070aff450fa41a03a759746d2bc7187255d24c..d7fb090e4ab977cc907d3b624dcbed62ae3fbfeb 100644
--- a/cmake/CMakeLists.txt
+++ b/cmake/CMakeLists.txt
@@ -4,7 +4,7 @@
# Created by Christoph Junghans and Richard Berger
cmake_minimum_required(VERSION 3.1)
-project(lammps)
+project(lammps LANGUAGES CXX)
set(SOVERSION 0)
set(LAMMPS_SOURCE_DIR ${CMAKE_CURRENT_SOURCE_DIR}/../src)
set(LAMMPS_LIB_SOURCE_DIR ${CMAKE_CURRENT_SOURCE_DIR}/../lib)
@@ -23,14 +23,22 @@ if(NOT CMAKE_BUILD_TYPE AND NOT CMAKE_CXX_FLAGS)
set(CMAKE_BUILD_TYPE Release CACHE STRING "Choose the type of build, options are: None Debug Release RelWithDebInfo MinSizeRel." FORCE)
endif(NOT CMAKE_BUILD_TYPE AND NOT CMAKE_CXX_FLAGS)
-# remove any style headers in the src dir
-file(GLOB SRC_STYLE_FILES ${LAMMPS_SOURCE_DIR}/style_*.h)
-if(SRC_STYLE_FILES)
- file(REMOVE ${SRC_STYLE_FILES})
+file(GLOB SRC_FILES ${LAMMPS_SOURCE_DIR}/*.cpp)
+list(SORT SRC_FILES)
+# check for files installed by make-based buildsystem
+# only run this time consuming check if there are new files
+if(NOT SRC_FILES STREQUAL SRC_FILES_CACHED)
+ file(GLOB SRC_PKG_FILES ${LAMMPS_SOURCE_DIR}/*/*.cpp)
+ message(STATUS "Running check for installed package (this might take a while)")
+ foreach(_SRC SRC_PKG_FILES)
+ get_filename_component(FILENAME "${_SRC}" NAME)
+ if(EXISTS ${LAMMPS_SOURCE_DIR}/${FILENAME})
+ message(FATAL_ERROR "Found packages installed by the make-based buildsystem, please run 'make -C ${LAMMPS_SOURCE_DIR} no-all purge'")
+ endif()
+ endforeach()
+ set(SRC_FILES_CACHED "${SRC_FILES}" CACHE INTERNAL "List of file in LAMMPS_SOURCE_DIR" FORCE)
endif()
-enable_language(CXX)
-
######################################################################
# compiler tests
# these need ot be done early (before further tests).
@@ -53,16 +61,19 @@ include(GNUInstallDirs)
set(LAMMPS_LINK_LIBS)
set(LAMMPS_DEPS)
set(LAMMPS_API_DEFINES)
-option(ENABLE_MPI "Build MPI version" OFF)
-if(ENABLE_MPI)
+
+find_package(MPI QUIET)
+option(BUILD_MPI "Build MPI version" ${MPI_FOUND})
+if(BUILD_MPI)
find_package(MPI REQUIRED)
- include_directories(${MPI_C_INCLUDE_PATH})
+ include_directories(${MPI_CXX_INCLUDE_PATH})
list(APPEND LAMMPS_LINK_LIBS ${MPI_CXX_LIBRARIES})
option(LAMMPS_LONGLONG_TO_LONG "Workaround if your system or MPI version does not recognize 'long long' data types" OFF)
if(LAMMPS_LONGLONG_TO_LONG)
add_definitions(-DLAMMPS_LONGLONG_TO_LONG)
endif()
else()
+ enable_language(C)
file(GLOB MPI_SOURCES ${LAMMPS_SOURCE_DIR}/STUBS/mpi.c)
add_library(mpi_stubs STATIC ${MPI_SOURCES})
include_directories(${LAMMPS_SOURCE_DIR}/STUBS)
@@ -108,14 +119,14 @@ set(OTHER_PACKAGES KIM PYTHON MSCG MPIIO VORONOI POEMS LATTE
USER-SMTBQ USER-SPH USER-TALLY USER-UEF USER-VTK USER-QUIP USER-QMMM)
set(ACCEL_PACKAGES USER-OMP KOKKOS OPT USER-INTEL GPU)
foreach(PKG ${DEFAULT_PACKAGES})
- option(ENABLE_${PKG} "Build ${PKG} Package" ${ENABLE_ALL})
+ option(PKG_${PKG} "Build ${PKG} Package" ${ENABLE_ALL})
endforeach()
foreach(PKG ${ACCEL_PACKAGES} ${OTHER_PACKAGES})
- option(ENABLE_${PKG} "Build ${PKG} Package" OFF)
+ option(PKG_${PKG} "Build ${PKG} Package" OFF)
endforeach()
macro(pkg_depends PKG1 PKG2)
- if(ENABLE_${PKG1} AND NOT ENABLE_${PKG2})
+ if(PKG_${PKG1} AND NOT PKG_${PKG2})
message(FATAL_ERROR "${PKG1} package needs LAMMPS to be build with ${PKG2}")
endif()
endmacro()
@@ -131,31 +142,50 @@ pkg_depends(USER-PHONON KSPACE)
######################################################
# packages with special compiler needs or external libs
######################################################
-if(ENABLE_REAX OR ENABLE_MEAM OR ENABLE_USER-QUIP OR ENABLE_USER-QMMM OR ENABLE_LATTE)
+if(PKG_REAX OR PKG_MEAM OR PKG_USER-QUIP OR PKG_USER-QMMM OR PKG_LATTE)
enable_language(Fortran)
- include(CheckFortranCompilerFlag)
- check_Fortran_compiler_flag("-fno-second-underscore" FC_HAS_NO_SECOND_UNDERSCORE)
endif()
-if(ENABLE_KOKKOS OR ENABLE_MSCG)
+if(PKG_MEAM)
+ enable_language(C)
+endif()
+
+if(PKG_KOKKOS OR PKG_MSCG)
# starting with CMake 3.1 this is all you have to do to enforce C++11
set(CMAKE_CXX_STANDARD 11) # C++11...
set(CMAKE_CXX_STANDARD_REQUIRED ON) #...is required...
set(CMAKE_CXX_EXTENSIONS OFF) #...without compiler extensions like gnu++11
endif()
-if(ENABLE_USER-OMP OR ENABLE_KOKKOS OR ENABLE_USER-INTEL)
+find_package(OpenMP QUIET)
+option(BUILD_OMP "Build with OpenMP support" ${OpenMP_FOUND})
+if(BUILD_OMP OR PKG_USER-OMP OR PKG_KOKKOS OR PKG_USER-INTEL)
find_package(OpenMP REQUIRED)
set (CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS}")
set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS}")
endif()
-if(ENABLE_KSPACE)
- set(FFT "KISSFFT" CACHE STRING "FFT library for KSPACE package")
- set_property(CACHE FFT PROPERTY STRINGS KISSFFT FFTW3 MKL FFTW2)
+if(PKG_KSPACE)
+ option(FFT_SINGLE "Use single precision FFT instead of double" OFF)
+ set(FFTW "FFTW3")
+ if(FFT_SINGLE)
+ set(FFTW "FFTW3F")
+ add_definitions(-DFFT_SINGLE)
+ endif()
+ find_package(${FFTW} QUIET)
+ if(${FFTW}_FOUND)
+ set(FFT "${FFTW}" CACHE STRING "FFT library for KSPACE package")
+ else()
+ set(FFT "KISSFFT" CACHE STRING "FFT library for KSPACE package")
+ endif()
+ set_property(CACHE FFT PROPERTY STRINGS KISSFFT ${FFTW} MKL)
if(NOT FFT STREQUAL "KISSFFT")
find_package(${FFT} REQUIRED)
- add_definitions(-DFFT_${FFT})
+ if(NOT FFT STREQUAL "FFTW3F")
+ add_definitions(-DFFT_FFTW)
+ else()
+ add_definitions(-DFFT_${FFT})
+ endif()
include_directories(${${FFT}_INCLUDE_DIRS})
list(APPEND LAMMPS_LINK_LIBS ${${FFT}_LIBRARIES})
endif()
@@ -166,22 +196,17 @@ if(ENABLE_KSPACE)
endif()
endif()
-if(ENABLE_MSCG OR ENABLE_USER-ATC OR ENABLE_USER-AWPMD OR ENABLE_USER-QUIP OR ENABLE_LATTE)
+if(PKG_MSCG OR PKG_USER-ATC OR PKG_USER-AWPMD OR PKG_USER-QUIP OR PKG_LATTE)
find_package(LAPACK)
if(NOT LAPACK_FOUND)
enable_language(Fortran)
file(GLOB LAPACK_SOURCES ${LAMMPS_LIB_SOURCE_DIR}/linalg/*.f)
add_library(linalg STATIC ${LAPACK_SOURCES})
- include(CheckFortranCompilerFlag)
- check_Fortran_compiler_flag("-fno-second-underscore" FC_HAS_NO_SECOND_UNDERSCORE)
- if(FC_HAS_NO_SECOND_UNDERSCORE)
- target_compile_options(linalg PRIVATE -fno-second-underscore)
- endif()
set(LAPACK_LIBRARIES linalg)
endif()
endif()
-if(ENABLE_PYTHON)
+if(PKG_PYTHON)
find_package(PythonInterp REQUIRED)
find_package(PythonLibs REQUIRED)
add_definitions(-DLMP_PYTHON)
@@ -197,16 +222,25 @@ if(ENABLE_PYTHON)
endif()
endif()
-find_package(JPEG)
-if(JPEG_FOUND)
+find_package(JPEG QUIET)
+option(WITH_JPEG "Enable JPEG support" ${JPEG_FOUND})
+if(WITH_JPEG)
+ find_package(JPEG REQUIRED)
add_definitions(-DLAMMPS_JPEG)
include_directories(${JPEG_INCLUDE_DIR})
list(APPEND LAMMPS_LINK_LIBS ${JPEG_LIBRARIES})
endif()
-find_package(PNG)
-find_package(ZLIB)
+find_package(PNG QUIET)
+find_package(ZLIB QUIET)
if(PNG_FOUND AND ZLIB_FOUND)
+ option(WITH_PNG "Enable PNG support" ON)
+else()
+ option(WITH_PNG "Enable PNG support" OFF)
+endif()
+if(WITH_PNG)
+ find_package(PNG REQUIRED)
+ find_package(ZLIB REQUIRED)
include_directories(${PNG_INCLUDE_DIRS} ${ZLIB_INCLUDE_DIRS})
list(APPEND LAMMPS_LINK_LIBS ${PNG_LIBRARIES} ${ZLIB_LIBRARIES})
add_definitions(-DLAMMPS_PNG)
@@ -214,25 +248,50 @@ endif()
find_program(GZIP_EXECUTABLE gzip)
find_package_handle_standard_args(GZIP REQUIRED_VARS GZIP_EXECUTABLE)
-if(GZIP_FOUND)
+option(WITH_GZIP "Enable GZIP support" ${GZIP_FOUND})
+if(WITH_GZIP)
+ if(NOT GZIP_FOUND)
+ message(FATAL_ERROR "gzip executable not found")
+ endif()
add_definitions(-DLAMMPS_GZIP)
endif()
find_program(FFMPEG_EXECUTABLE ffmpeg)
find_package_handle_standard_args(FFMPEG REQUIRED_VARS FFMPEG_EXECUTABLE)
-if(FFMPEG_FOUND)
+option(WITH_FFMPEG "Enable FFMPEG support" ${FFMPEG_FOUND})
+if(WITH_FFMPEG)
+ if(NOT FFMPEG_FOUND)
+ message(FATAL_ERROR "ffmpeg executable not found")
+ endif()
add_definitions(-DLAMMPS_FFMPEG)
endif()
-if(ENABLE_VORONOI)
- find_package(VORO REQUIRED) #some distros
+if(PKG_VORONOI)
+ option(DOWNLOAD_VORO "Download voro++ (instead of using the system's one)" OFF)
+ if(DOWNLOAD_VORO)
+ include(ExternalProject)
+ ExternalProject_Add(voro_build
+ URL http://math.lbl.gov/voro++/download/dir/voro++-0.4.6.tar.gz
+ URL_MD5 2338b824c3b7b25590e18e8df5d68af9
+ CONFIGURE_COMMAND "" BUILD_IN_SOURCE 1 INSTALL_COMMAND ""
+ )
+ ExternalProject_get_property(voro_build SOURCE_DIR)
+ set(VORO_LIBRARIES ${SOURCE_DIR}/src/libvoro++.a)
+ set(VORO_INCLUDE_DIRS ${SOURCE_DIR}/src)
+ list(APPEND LAMMPS_DEPS voro_build)
+ else()
+ find_package(VORO)
+ if(NOT VORO_FOUND)
+ message(FATAL_ERROR "VORO not found, help CMake to find it by setting VORO_LIBRARY and VORO_INCLUDE_DIR, or set DOWNLOAD_VORO=ON to download it")
+ endif()
+ endif()
include_directories(${VORO_INCLUDE_DIRS})
list(APPEND LAMMPS_LINK_LIBS ${VORO_LIBRARIES})
endif()
-if(ENABLE_LATTE)
- find_package(LATTE QUIET)
- if(NOT LATTE_FOUND)
+if(PKG_LATTE)
+ option(DOWNLOAD_LATTE "Download latte (instead of using the system's one)" OFF)
+ if(DOWNLOAD_LATTE)
message(STATUS "LATTE not found - we will build our own")
include(ExternalProject)
ExternalProject_Add(latte_build
@@ -244,48 +303,73 @@ if(ENABLE_LATTE)
ExternalProject_get_property(latte_build INSTALL_DIR)
set(LATTE_LIBRARIES ${INSTALL_DIR}/${CMAKE_INSTALL_LIBDIR}/liblatte.a)
list(APPEND LAMMPS_DEPS latte_build)
+ else()
+ find_package(LATTE)
+ if(NOT LATTE_FOUND)
+ message(FATAL_ERROR "LATTE not found, help CMake to find it by setting LATTE_LIBRARY, or set DOWNLOAD_LATTE=ON to download it")
+ endif()
endif()
list(APPEND LAMMPS_LINK_LIBS ${LATTE_LIBRARIES} ${LAPACK_LIBRARIES} ${CMAKE_Fortran_IMPLICIT_LINK_LIBRARIES})
endif()
-if(ENABLE_USER-MOLFILE)
+if(PKG_USER-MOLFILE)
add_library(molfile INTERFACE)
target_include_directories(molfile INTERFACE ${LAMMPS_LIB_SOURCE_DIR}/molfile)
target_link_libraries(molfile INTERFACE ${CMAKE_DL_LIBS})
list(APPEND LAMMPS_LINK_LIBS molfile)
endif()
-if(ENABLE_USER-NETCDF)
+if(PKG_USER-NETCDF)
find_package(NetCDF REQUIRED)
include_directories(NETCDF_INCLUDE_DIR)
list(APPEND LAMMPS_LINK_LIBS ${NETCDF_LIBRARY})
add_definitions(-DLMP_HAS_NETCDF -DNC_64BIT_DATA=0x0020)
endif()
-if(ENABLE_USER-SMD)
- find_package(Eigen3 REQUIRED)
+if(PKG_USER-SMD)
+ option(DOWNLOAD_Eigen3 "Download Eigen3 (instead of using the system's one)" OFF)
+ if(DOWNLOAD_Eigen3)
+ include(ExternalProject)
+ ExternalProject_Add(Eigen3_build
+ URL http://bitbucket.org/eigen/eigen/get/3.3.4.tar.gz
+ URL_MD5 1a47e78efe365a97de0c022d127607c3
+ CMAKE_ARGS -DCMAKE_INSTALL_PREFIX=<INSTALL_DIR> -DEIGEN_TEST_NOQT=ON
+ -DCMAKE_DISABLE_FIND_PACKAGE_LAPACK=ON -DCMAKE_DISABLE_FIND_PACKAGE_Cholmod=ON -DCMAKE_DISABLE_FIND_PACKAGE_Umfpack=ON -DCMAKE_DISABLE_FIND_PACKAGE_SuperLU=ON
+ -DCMAKE_DISABLE_FIND_PACKAGE_PASTIX=ON -DCMAKE_DISABLE_FIND_PACKAGE_SPQR=ON -DCMAKE_DISABLE_FIND_PACKAGE_Boost=ON -DCMAKE_DISABLE_FIND_PACKAGE_CUDA=ON
+ -DCMAKE_DISABLE_FIND_PACKAGE_FFTW=ON -DCMAKE_DISABLE_FIND_PACKAGE_MPFR=ON -DCMAKE_DISABLE_FIND_PACKAGE_OpenGL=ON
+ )
+ ExternalProject_get_property(Eigen3_build INSTALL_DIR)
+ set(EIGEN3_INCLUDE_DIR ${INSTALL_DIR}/include/eigen3)
+ list(APPEND LAMMPS_DEPS Eigen3_build)
+ else()
+ find_package(Eigen3)
+ if(NOT Eigen3_FOUND)
+ message(FATAL_ERROR "Eigen3 not found, help CMake to find it by setting EIGEN3_INCLUDE_DIR, or set DOWNLOAD_Eigen3=ON to download it")
+ endif()
+ endif()
include_directories(${EIGEN3_INCLUDE_DIR})
endif()
-if(ENABLE_USER-QUIP)
+if(PKG_USER-QUIP)
find_package(QUIP REQUIRED)
list(APPEND LAMMPS_LINK_LIBS ${QUIP_LIBRARIES} ${LAPACK_LIBRARIES} ${CMAKE_Fortran_IMPLICIT_LINK_LIBRARIES})
endif()
-if(ENABLE_USER-QMMM)
+if(PKG_USER-QMMM)
+ message(WARNING "Building QMMM with CMake is still experimental")
find_package(QE REQUIRED)
include_directories(${QE_INCLUDE_DIRS})
list(APPEND LAMMPS_LINK_LIBS ${QE_LIBRARIES} ${CMAKE_Fortran_IMPLICIT_LINK_LIBRARIES})
endif()
-if(ENABLE_USER-VTK)
+if(PKG_USER-VTK)
find_package(VTK REQUIRED NO_MODULE)
include(${VTK_USE_FILE})
add_definitions(-DLAMMPS_VTK)
list(APPEND LAMMPS_LINK_LIBS ${VTK_LIBRARIES})
endif()
-if(ENABLE_KIM)
+if(PKG_KIM)
find_package(KIM QUIET)
if(NOT KIM_FOUND)
message(STATUS "KIM not found - we will build our own")
@@ -305,7 +389,7 @@ if(ENABLE_KIM)
include_directories(${KIM_INCLUDE_DIRS})
endif()
-if(ENABLE_MSCG)
+if(PKG_MSCG)
find_package(GSL REQUIRED)
set(LAMMPS_LIB_MSCG_BIN_DIR ${LAMMPS_LIB_BINARY_DIR}/mscg)
set(MSCG_TARBALL ${LAMMPS_LIB_MSCG_BIN_DIR}/MS-CG-master.zip)
@@ -329,12 +413,18 @@ if(ENABLE_MSCG)
target_link_libraries(mscg ${GSL_LIBRARIES} ${LAPACK_LIBRARIES})
endif()
+if(PKG_COMPRESS)
+ find_package(ZLIB REQUIRED)
+ include_directories(${ZLIB_INCLUDE_DIRS})
+ list(APPEND LAMMPS_LINK_LIBS ${ZLIB_LIBRARIES})
+endif()
+
########################################################################
# Basic system tests (standard libraries, headers, functions, types) #
########################################################################
-include(CheckIncludeFile)
+include(CheckIncludeFileCXX)
foreach(HEADER math.h)
- check_include_file(${HEADER} FOUND_${HEADER})
+ check_include_file_cxx(${HEADER} FOUND_${HEADER})
if(NOT FOUND_${HEADER})
message(FATAL_ERROR "Could not find needed header - ${HEADER}")
endif(NOT FOUND_${HEADER})
@@ -378,7 +468,7 @@ foreach(PKG ${DEFAULT_PACKAGES} ${OTHER_PACKAGES})
DetectAndRemovePackageHeader(${LAMMPS_SOURCE_DIR}/${FNAME})
endforeach()
- if(ENABLE_${PKG})
+ if(PKG_${PKG})
# detects styles in package and adds them to global list
RegisterStyles(${${PKG}_SOURCES_DIR})
@@ -392,7 +482,7 @@ endforeach()
############################################
foreach(SIMPLE_LIB REAX MEAM POEMS USER-ATC USER-AWPMD USER-COLVARS USER-H5MD
USER-QMMM)
- if(ENABLE_${SIMPLE_LIB})
+ if(PKG_${SIMPLE_LIB})
string(REGEX REPLACE "^USER-" "" PKG_LIB "${SIMPLE_LIB}")
string(TOLOWER "${PKG_LIB}" PKG_LIB)
file(GLOB_RECURSE ${PKG_LIB}_SOURCES ${LAMMPS_LIB_SOURCE_DIR}/${PKG_LIB}/*.F
@@ -413,40 +503,26 @@ foreach(SIMPLE_LIB REAX MEAM POEMS USER-ATC USER-AWPMD USER-COLVARS USER-H5MD
endif()
endforeach()
-if(ENABLE_USER-AWPMD)
+if(PKG_USER-AWPMD)
target_link_libraries(awpmd ${LAPACK_LIBRARIES})
endif()
-if(ENABLE_USER-ATC)
+if(PKG_USER-ATC)
target_link_libraries(atc ${LAPACK_LIBRARIES})
endif()
-if(ENABLE_USER-H5MD)
+if(PKG_USER-H5MD)
find_package(HDF5 REQUIRED)
target_link_libraries(h5md ${HDF5_LIBRARIES})
target_include_directories(h5md PRIVATE ${HDF5_INCLUDE_DIRS})
endif()
-if(ENABLE_MEAM AND FC_HAS_NO_SECOND_UNDERSCORE)
- foreach(FSRC ${meam_SOURCES})
- string(REGEX REPLACE "^.*\\." "" FEXT "${FSRC}")
- list(FIND CMAKE_Fortran_SOURCE_FILE_EXTENSIONS "${FEXT}" FINDEX)
- if(FINDEX GREATER -1)
- set_property(SOURCE ${FSRC} APPEND PROPERTY COMPILE_FLAGS "-fno-second-underscore")
- endif()
- endforeach()
-endif()
-
-if(ENABLE_REAX AND FC_HAS_NO_SECOND_UNDERSCORE)
- target_compile_options(reax PRIVATE -fno-second-underscore)
-endif()
-
######################################################################
# packages which selectively include variants based on enabled styles
# e.g. accelerator packages
######################################################################
-if(ENABLE_USER-OMP)
+if(PKG_USER-OMP)
set(USER-OMP_SOURCES_DIR ${LAMMPS_SOURCE_DIR}/USER-OMP)
set(USER-OMP_SOURCES ${USER-OMP_SOURCES_DIR}/thr_data.cpp
${USER-OMP_SOURCES_DIR}/thr_omp.cpp
@@ -463,7 +539,7 @@ if(ENABLE_USER-OMP)
include_directories(${USER-OMP_SOURCES_DIR})
endif()
-if(ENABLE_KOKKOS)
+if(PKG_KOKKOS)
set(LAMMPS_LIB_KOKKOS_SRC_DIR ${LAMMPS_LIB_SOURCE_DIR}/kokkos)
set(LAMMPS_LIB_KOKKOS_BIN_DIR ${LAMMPS_LIB_BINARY_DIR}/kokkos)
add_definitions(-DLMP_KOKKOS)
@@ -499,7 +575,7 @@ if(ENABLE_KOKKOS)
RegisterNBinStyle(${KOKKOS_PKG_SOURCES_DIR}/nbin_kokkos.h)
RegisterNPairStyle(${KOKKOS_PKG_SOURCES_DIR}/npair_kokkos.h)
- if(ENABLE_USER-DPD)
+ if(PKG_USER-DPD)
get_property(KOKKOS_PKG_SOURCES GLOBAL PROPERTY KOKKOS_PKG_SOURCES)
list(APPEND KOKKOS_PKG_SOURCES ${KOKKOS_PKG_SOURCES_DIR}/npair_ssa_kokkos.cpp)
RegisterNPairStyle(${KOKKOS_PKG_SOURCES_DIR}/npair_ssa_kokkos.h)
@@ -512,7 +588,7 @@ if(ENABLE_KOKKOS)
include_directories(${KOKKOS_PKG_SOURCES_DIR})
endif()
-if(ENABLE_OPT)
+if(PKG_OPT)
set(OPT_SOURCES_DIR ${LAMMPS_SOURCE_DIR}/OPT)
set(OPT_SOURCES)
set_property(GLOBAL PROPERTY "OPT_SOURCES" "${OPT_SOURCES}")
@@ -526,7 +602,7 @@ if(ENABLE_OPT)
include_directories(${OPT_SOURCES_DIR})
endif()
-if(ENABLE_USER-INTEL)
+if(PKG_USER-INTEL)
set(USER-INTEL_SOURCES_DIR ${LAMMPS_SOURCE_DIR}/USER-INTEL)
set(USER-INTEL_SOURCES ${USER-INTEL_SOURCES_DIR}/intel_preprocess.h
${USER-INTEL_SOURCES_DIR}/intel_buffers.h
@@ -550,7 +626,7 @@ if(ENABLE_USER-INTEL)
include_directories(${USER-INTEL_SOURCES_DIR})
endif()
-if(ENABLE_GPU)
+if(PKG_GPU)
set(GPU_SOURCES_DIR ${LAMMPS_SOURCE_DIR}/GPU)
set(GPU_SOURCES ${GPU_SOURCES_DIR}/gpu_extra.h
${GPU_SOURCES_DIR}/fix_gpu.h
@@ -687,6 +763,7 @@ include_directories(${LAMMPS_STYLE_HEADERS_DIR})
# Actually add executable and lib to build
############################################
add_library(lammps ${LIB_SOURCES})
+list(REMOVE_DUPLICATES LAMMPS_LINK_LIBS)
target_link_libraries(lammps ${LAMMPS_LINK_LIBS})
if(LAMMPS_DEPS)
add_dependencies(lammps ${LAMMPS_DEPS})
@@ -712,31 +789,60 @@ endif()
# Print package summary
##################################
foreach(PKG ${DEFAULT_PACKAGES} ${OTHER_PACKAGES} ${ACCEL_PACKAGES})
- if(ENABLE_${PKG})
+ if(PKG_${PKG})
message(STATUS "Building package: ${PKG}")
endif()
endforeach()
string(TOUPPER "${CMAKE_BUILD_TYPE}" BTYPE)
+get_directory_property(CPPFLAGS DIRECTORY ${CMAKE_SOURCE_DIR} COMPILE_DEFINITIONS)
+include(FeatureSummary)
+feature_summary(INCLUDE_QUIET_PACKAGES WHAT ALL)
message(STATUS "<<< Build configuration >>>
Build type ${CMAKE_BUILD_TYPE}
Install path ${CMAKE_INSTALL_PREFIX}
Compilers and Flags:
C++ Compiler ${CMAKE_CXX_COMPILER}
Type ${CMAKE_CXX_COMPILER_ID}
- C++ Flags ${CMAKE_CXX_FLAGS} ${CMAKE_CXX_FLAGS_${BTYPE}}")
+ Version ${CMAKE_CXX_COMPILER_VERSION}
+ C++ Flags ${CMAKE_CXX_FLAGS} ${CMAKE_CXX_FLAGS_${BTYPE}}
+ Defines ${CPPFLAGS}")
get_property(LANGUAGES GLOBAL PROPERTY ENABLED_LANGUAGES)
-if(LANGUAGES MATCHES ".*Fortran.*")
+list (FIND LANGUAGES "Fortran" _index)
+if (${_index} GREATER -1)
message(STATUS "Fortran Compiler ${CMAKE_Fortran_COMPILER}
Type ${CMAKE_Fortran_COMPILER_ID}
+ Version ${CMAKE_Fortran_COMPILER_VERSION}
Fortran Flags ${CMAKE_Fortran_FLAGS} ${CMAKE_Fortran_FLAGS_${BTYPE}}")
endif()
-message(STATUS "Linker flags:
+list (FIND LANGUAGES "C" _index)
+if (${_index} GREATER -1)
+ message(STATUS "C Compiler ${CMAKE_C_COMPILER}
+ Type ${CMAKE_C_COMPILER_ID}
+ Version ${CMAKE_C_COMPILER_VERSION}
+ C Flags ${CMAKE_C_FLAGS} ${CMAKE_C_FLAGS_${BTYPE}}")
+endif()
+if(CMAKE_EXE_LINKER_FLAGS)
+ message(STATUS "Linker flags:
Executable ${CMAKE_EXE_LINKER_FLAGS}")
+ endif()
if(BUILD_SHARED_LIBS)
- message(STATUS "Shared libries ${CMAKE_SHARED_LINKER_FLAGS}")
+ message(STATUS "Shared libraries ${CMAKE_SHARED_LINKER_FLAGS}")
else()
- message(STATUS "Static libries ${CMAKE_STATIC_LINKER_FLAGS}")
+ message(STATUS "Static libraries ${CMAKE_STATIC_LINKER_FLAGS}")
endif()
message(STATUS "Link libraries: ${LAMMPS_LINK_LIBS}")
-
+if(BUILD_MPI)
+ message(STATUS "Using mpi with headers in ${MPI_CXX_INCLUDE_PATH} and ${MPI_CXX_LIBRARIES}")
+endif()
+if(ENABLED_GPU)
+ message(STATUS "GPU Api: ${GPU_API}")
+ if(GPU_API STREQUAL "CUDA")
+ message(STATUS "GPU Arch: ${GPU_ARCH}")
+ elseif(GPU_API STREQUAL "OpenCL")
+ message(STATUS "OCL Tune: ${OCL_TUNE}")
+ endif()
+endif()
+if(PKG_KSPACE)
+ message(STATUS "Using ${FFT} as FFT")
+endif()
diff --git a/cmake/Modules/FindFFTW2.cmake b/cmake/Modules/FindFFTW2.cmake
deleted file mode 100644
index c77e6cf8e9d126d5ebb4223e6e7cf53ef25b60b2..0000000000000000000000000000000000000000
--- a/cmake/Modules/FindFFTW2.cmake
+++ /dev/null
@@ -1,22 +0,0 @@
-# - Find fftw2
-# Find the native FFTW2 headers and libraries.
-#
-# FFTW2_INCLUDE_DIRS - where to find fftw2.h, etc.
-# FFTW2_LIBRARIES - List of libraries when using fftw2.
-# FFTW2_FOUND - True if fftw2 found.
-#
-
-find_path(FFTW2_INCLUDE_DIR fftw.h)
-
-find_library(FFTW2_LIBRARY NAMES fftw)
-
-set(FFTW2_LIBRARIES ${FFTW2_LIBRARY})
-set(FFTW2_INCLUDE_DIRS ${FFTW2_INCLUDE_DIR})
-
-include(FindPackageHandleStandardArgs)
-# handle the QUIETLY and REQUIRED arguments and set FFTW2_FOUND to TRUE
-# if all listed variables are TRUE
-
-find_package_handle_standard_args(FFTW2 DEFAULT_MSG FFTW2_LIBRARY FFTW2_INCLUDE_DIR)
-
-mark_as_advanced(FFTW2_INCLUDE_DIR FFTW2_LIBRARY )
diff --git a/cmake/Modules/FindFFTW3F.cmake b/cmake/Modules/FindFFTW3F.cmake
new file mode 100644
index 0000000000000000000000000000000000000000..92d1e85e791ec362a15909bd09dccdc046efc2a5
--- /dev/null
+++ b/cmake/Modules/FindFFTW3F.cmake
@@ -0,0 +1,25 @@
+# - Find fftw3f
+# Find the native FFTW3F headers and libraries.
+#
+# FFTW3F_INCLUDE_DIRS - where to find fftw3f.h, etc.
+# FFTW3F_LIBRARIES - List of libraries when using fftw3f.
+# FFTW3F_FOUND - True if fftw3f found.
+#
+
+find_package(PkgConfig)
+
+pkg_check_modules(PC_FFTW3F fftw3f)
+find_path(FFTW3F_INCLUDE_DIR fftw3.h HINTS ${PC_FFTW3F_INCLUDE_DIRS})
+
+find_library(FFTW3F_LIBRARY NAMES fftw3f HINTS ${PC_FFTW3F_LIBRARY_DIRS})
+
+set(FFTW3F_LIBRARIES ${FFTW3F_LIBRARY})
+set(FFTW3F_INCLUDE_DIRS ${FFTW3F_INCLUDE_DIR})
+
+include(FindPackageHandleStandardArgs)
+# handle the QUIETLY and REQUIRED arguments and set FFTW3F_FOUND to TRUE
+# if all listed variables are TRUE
+
+find_package_handle_standard_args(FFTW3F DEFAULT_MSG FFTW3F_LIBRARY FFTW3F_INCLUDE_DIR)
+
+mark_as_advanced(FFTW3F_INCLUDE_DIR FFTW3F_LIBRARY )
diff --git a/doc/src/Eqs/dihedral_table_cut.jpg b/doc/src/Eqs/dihedral_table_cut.jpg
new file mode 100644
index 0000000000000000000000000000000000000000..c124184c1773c60fe2564d8cabce2f6c7d7883c0
Binary files /dev/null and b/doc/src/Eqs/dihedral_table_cut.jpg differ
diff --git a/doc/src/Eqs/dihedral_table_cut.tex b/doc/src/Eqs/dihedral_table_cut.tex
new file mode 100644
index 0000000000000000000000000000000000000000..3cc1d331a8978d76e668de3f83cf0d4264be3301
--- /dev/null
+++ b/doc/src/Eqs/dihedral_table_cut.tex
@@ -0,0 +1,11 @@
+\documentclass[12pt]{article}
+\pagestyle{empty}
+\begin{document}
+
+\begin{eqnarray*}
+ f(\theta) & = & K \qquad\qquad\qquad\qquad\qquad\qquad \theta < \theta_1 \\
+ f(\theta) & = & K \left(1-\frac{(\theta - \theta_1)^2}{(\theta_2 - \theta_1)^2}\right) \qquad \theta_1 < \theta < \theta_2
+\end{eqnarray*}
+
+\end{document}
+
diff --git a/doc/src/Manual.txt b/doc/src/Manual.txt
index 077927ccc7ac2a62d60b988480cde725187a2736..6b2c7daac7d682c05469e44db60b9e0aad02f5c5 100644
--- a/doc/src/Manual.txt
+++ b/doc/src/Manual.txt
@@ -1,7 +1,7 @@
<!-- HTML_ONLY -->
<HEAD>
<TITLE>LAMMPS Users Manual</TITLE>
-<META NAME="docnumber" CONTENT="20 Apr 2018 version">
+<META NAME="docnumber" CONTENT="11 May 2018 version">
<META NAME="author" CONTENT="http://lammps.sandia.gov - Sandia National Laboratories">
<META NAME="copyright" CONTENT="Copyright (2003) Sandia Corporation. This software and manual is distributed under the GNU General Public License.">
</HEAD>
@@ -19,7 +19,7 @@
:line
LAMMPS Documentation :c,h1
-20 Apr 2018 version :c,h2
+11 May 2018 version :c,h2
Version info: :h3
diff --git a/doc/src/Section_commands.txt b/doc/src/Section_commands.txt
index cdbdec9ee3ef6e7a8d13cbb28b53d9d41034c487..3dabdbeaa1e1091db35b798f19cf9f3b2b74848f 100644
--- a/doc/src/Section_commands.txt
+++ b/doc/src/Section_commands.txt
@@ -1212,7 +1212,8 @@ package"_Section_start.html#start_3.
"nharmonic (o)"_dihedral_nharmonic.html,
"quadratic (o)"_dihedral_quadratic.html,
"spherical (o)"_dihedral_spherical.html,
-"table (o)"_dihedral_table.html :tb(c=4,ea=c)
+"table (o)"_dihedral_table.html,
+"table/cut"_dihedral_table_cut.html :tb(c=4,ea=c)
:line
diff --git a/doc/src/dihedral_table_cut.txt b/doc/src/dihedral_table_cut.txt
new file mode 100644
index 0000000000000000000000000000000000000000..1c83d4ffa09f7221ac8c3859594abc754b9b28c4
--- /dev/null
+++ b/doc/src/dihedral_table_cut.txt
@@ -0,0 +1,205 @@
+"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
+
+:link(lws,http://lammps.sandia.gov)
+:link(ld,Manual.html)
+:link(lc,Section_commands.html#comm)
+
+:line
+
+dihedral_style table/cut command :h3
+
+[Syntax:]
+
+dihedral_style table/cut style Ntable :pre
+
+style = {linear} or {spline} = method of interpolation
+Ntable = size of the internal lookup table :ul
+
+[Examples:]
+
+dihedral_style table/cut spline 400
+dihedral_style table/cut linear 1000
+dihedral_coeff 1 aat 1.0 177 180 file.table DIH_TABLE1
+dihedral_coeff 2 aat 0.5 170 180 file.table DIH_TABLE2 :pre
+
+[Description:]
+
+The {table/cut} dihedral style creates interpolation tables of length
+{Ntable} from dihedral potential and derivative values listed in a
+file(s) as a function of the dihedral angle "phi". In addition, an
+analytic cutoff that is quadratic in the bond-angle (theta) is applied
+in order to regularize the dihedral interaction. The dihedral table
+files are read by the "dihedral_coeff"_dihedral_coeff.html command.
+
+The interpolation tables are created by fitting cubic splines to the
+file values and interpolating energy and derivative values at each of
+{Ntable} dihedral angles. During a simulation, these tables are used
+to interpolate energy and force values on individual atoms as
+needed. The interpolation is done in one of 2 styles: {linear} or
+{spline}.
+
+For the {linear} style, the dihedral angle (phi) is used to find 2
+surrounding table values from which an energy or its derivative is
+computed by linear interpolation.
+
+For the {spline} style, cubic spline coefficients are computed and
+stored at each of the {Ntable} evenly-spaced values in the
+interpolated table. For a given dihedral angle (phi), the appropriate
+coefficients are chosen from this list, and a cubic polynomial is used
+to compute the energy and the derivative at this angle.
+
+The following coefficients must be defined for each dihedral type via
+the "dihedral_coeff"_dihedral_coeff.html command as in the example
+above.
+
+style (aat)
+cutoff prefactor
+cutoff angle1
+cutoff angle2
+filename
+keyword :ul
+
+The cutoff dihedral style uses a tabulated dihedral interaction with a
+cutoff function:
+
+:c,image(Eqs/dihedral_table_cut.jpg)
+
+The cutoff specifies an prefactor to the cutoff function. While this value
+would ordinarily equal 1 there may be situations where the value should change.
+
+The cutoff angle1 specifies the angle (in degrees) below which the dihedral
+interaction is unmodified, i.e. the cutoff function is 1.
+
+The cutoff function is applied between angle1 and angle2, which is the angle at
+which the cutoff function drops to zero. The value of zero effectively "turns
+off" the dihedral interaction.
+
+The filename specifies a file containing tabulated energy and
+derivative values. The keyword specifies a section of the file. The
+format of this file is described below.
+
+:line
+
+The format of a tabulated file is as follows (without the
+parenthesized comments). It can begin with one or more comment
+or blank lines.
+
+# Table of the potential and its negative derivative :pre
+
+DIH_TABLE1 (keyword is the first text on line)
+N 30 DEGREES (N, NOF, DEGREES, RADIANS, CHECKU/F)
+ (blank line)
+1 -168.0 -1.40351172223 0.0423346818422
+2 -156.0 -1.70447981034 0.00811786522531
+3 -144.0 -1.62956100432 -0.0184129719987
+...
+30 180.0 -0.707106781187 0.0719306095245 :pre
+
+# Example 2: table of the potential. Forces omitted :pre
+
+DIH_TABLE2
+N 30 NOF CHECKU testU.dat CHECKF testF.dat :pre
+
+1 -168.0 -1.40351172223
+2 -156.0 -1.70447981034
+3 -144.0 -1.62956100432
+...
+30 180.0 -0.707106781187 :pre
+
+A section begins with a non-blank line whose 1st character is not a
+"#"; blank lines or lines starting with "#" can be used as comments
+between sections. The first line begins with a keyword which
+identifies the section. The line can contain additional text, but the
+initial text must match the argument specified in the
+"dihedral_coeff"_dihedral_coeff.html command. The next line lists (in
+any order) one or more parameters for the table. Each parameter is a
+keyword followed by one or more numeric values.
+
+Following a blank line, the next N lines list the tabulated values. On
+each line, the 1st value is the index from 1 to N, the 2nd value is
+the angle value, the 3rd value is the energy (in energy units), and
+the 4th is -dE/d(phi) also in energy units). The 3rd term is the
+energy of the 4-atom configuration for the specified angle. The 4th
+term (when present) is the negative derivative of the energy with
+respect to the angle (in degrees, or radians depending on whether the
+user selected DEGREES or RADIANS). Thus the units of the last term
+are still energy, not force. The dihedral angle values must increase
+from one line to the next.
+
+Dihedral table splines are cyclic. There is no discontinuity at 180
+degrees (or at any other angle). Although in the examples above, the
+angles range from -180 to 180 degrees, in general, the first angle in
+the list can have any value (positive, zero, or negative). However
+the {range} of angles represented in the table must be {strictly} less
+than 360 degrees (2pi radians) to avoid angle overlap. (You may not
+supply entries in the table for both 180 and -180, for example.) If
+the user's table covers only a narrow range of dihedral angles,
+strange numerical behavior can occur in the large remaining gap.
+
+[Parameters:]
+
+The parameter "N" is required and its value is the number of table
+entries that follow. Note that this may be different than the N
+specified in the "dihedral_style table"_dihedral_style.html command.
+Let {Ntable} is the number of table entries requested dihedral_style
+command, and let {Nfile} be the parameter following "N" in the
+tabulated file ("30" in the sparse example above). What LAMMPS does
+is a preliminary interpolation by creating splines using the {Nfile}
+tabulated values as nodal points. It uses these to interpolate as
+needed to generate energy and derivative values at {Ntable} different
+points (which are evenly spaced over a 360 degree range, even if the
+angles in the file are not). The resulting tables of length {Ntable}
+are then used as described above, when computing energy and force for
+individual dihedral angles and their atoms. This means that if you
+want the interpolation tables of length {Ntable} to match exactly what
+is in the tabulated file (with effectively nopreliminary
+interpolation), you should set {Ntable} = {Nfile}. To insure the
+nodal points in the user's file are aligned with the interpolated
+table entries, the angles in the table should be integer multiples of
+360/{Ntable} degrees, or 2*PI/{Ntable} radians (depending on your
+choice of angle units).
+
+The optional "NOF" keyword allows the user to omit the forces
+(negative energy derivatives) from the table file (normally located in
+the 4th column). In their place, forces will be calculated
+automatically by differentiating the potential energy function
+indicated by the 3rd column of the table (using either linear or
+spline interpolation).
+
+The optional "DEGREES" keyword allows the user to specify angles in
+degrees instead of radians (default).
+
+The optional "RADIANS" keyword allows the user to specify angles in
+radians instead of degrees. (Note: This changes the way the forces
+are scaled in the 4th column of the data file.)
+
+The optional "CHECKU" keyword is followed by a filename. This allows
+the user to save all of the the {Ntable} different entries in the
+interpolated energy table to a file to make sure that the interpolated
+function agrees with the user's expectations. (Note: You can
+temporarily increase the {Ntable} parameter to a high value for this
+purpose. "{Ntable}" is explained above.)
+
+The optional "CHECKF" keyword is analogous to the "CHECKU" keyword.
+It is followed by a filename, and it allows the user to check the
+interpolated force table. This option is available even if the user
+selected the "NOF" option.
+
+Note that one file can contain many sections, each with a tabulated
+potential. LAMMPS reads the file section by section until it finds one
+that matches the specified keyword.
+
+[Restrictions:]
+
+This dihedral style can only be used if LAMMPS was built with the
+USER-MISC package. See the "Making LAMMPS"_Section_start.html#start_3
+section for more info on packages.
+
+[Related commands:]
+
+"dihedral_coeff"_dihedral_coeff.html, "dihedral_style table"_dihedral_table.html
+
+[Default:] none
+
+:link(dihedralcut-Salerno)
+[(Salerno)] Salerno, Bernstein, J Chem Theory Comput, --, ---- (2018).
diff --git a/doc/src/dihedrals.txt b/doc/src/dihedrals.txt
index 500a6a52bf69b7beb641c414f6dee354d6ca6b7b..a862bf50a0a9c357ffeed1129e3b67423259447c 100644
--- a/doc/src/dihedrals.txt
+++ b/doc/src/dihedrals.txt
@@ -19,6 +19,7 @@ Dihedral Styles :h1
dihedral_quadratic
dihedral_spherical
dihedral_table
+ dihedral_table_cut
dihedral_zero
dihedral_charmm
dihedral_class2
diff --git a/doc/src/dump_modify.txt b/doc/src/dump_modify.txt
index 6de6de545e720d518316ea9d31c82fb825fd8c9c..3230507dc3eaebc8bce634acf6a7ac2e3c2c7707 100644
--- a/doc/src/dump_modify.txt
+++ b/doc/src/dump_modify.txt
@@ -15,7 +15,7 @@ dump_modify dump-ID keyword values ... :pre
dump-ID = ID of dump to modify :ulb,l
one or more keyword/value pairs may be appended :l
these keywords apply to various dump styles :l
-keyword = {append} or {at} or {buffer} or {delay} or {element} or {every} or {fileper} or {first} or {flush} or {format} or {image} or {label} or {nfile} or {pad} or {precision} or {region} or {scale} or {sort} or {thresh} or {unwrap} :l
+keyword = {append} or {at} or {buffer} or {delay} or {element} or {every} or {fileper} or {first} or {flush} or {format} or {image} or {label} or {maxfiles} or {nfile} or {pad} or {precision} or {region} or {scale} or {sort} or {thresh} or {unwrap} :l
{append} arg = {yes} or {no}
{at} arg = N
N = index of frame written upon first dump
@@ -37,6 +37,8 @@ keyword = {append} or {at} or {buffer} or {delay} or {element} or {every} or {fi
{image} arg = {yes} or {no}
{label} arg = string
string = character string (e.g. BONDS) to use in header of dump local file
+ {maxfiles} arg = Fmax
+ Fmax = keep only the most recent {Fmax} snapshots (one snapshot per file)
{nfile} arg = Nf
Nf = write this many files, one from each of Nf processors
{pad} arg = Nchar = # of characters to convert timestep to
@@ -364,6 +366,20 @@ e.g. BONDS or ANGLES.
:line
+The {maxfiles} keyword can only be used when a '*' wildcard is
+included in the dump file name, i.e. when writing a new file(s) for
+each snapshot. The specified {Fmax} is how many snapshots will be
+kept. Once this number is reached, the file(s) containing the oldest
+snapshot is deleted before a new dump file is written. If the
+specified {Fmax} <= 0, then all files are retained.
+
+This can be useful for debugging, especially if you don't know on what
+timestep something bad will happen, e.g. when LAMMPS will exit with an
+error. You can dump every timestep, and limit the number of dump
+files produced, even if you run for 1000s of steps.
+
+:line
+
The {nfile} or {fileper} keywords can be used in conjunction with the
"%" wildcard character in the specified dump file name, for all dump
styles except the {dcd}, {image}, {movie}, {xtc}, and {xyz} styles
@@ -901,6 +917,7 @@ flush = yes
format = %d and %g for each integer or floating point value
image = no
label = ENTRIES
+maxifiles = -1
nfile = 1
pad = 0
pbc = no
diff --git a/doc/src/lammps.book b/doc/src/lammps.book
index 1f4bf157f4829df2446970b8c9550f7098b765f8..0764c593f7e56cc251588d800a53d06a11e232b4 100644
--- a/doc/src/lammps.book
+++ b/doc/src/lammps.book
@@ -582,6 +582,7 @@ dihedral_opls.html
dihedral_quadratic.html
dihedral_spherical.html
dihedral_table.html
+dihedral_table_cut.html
dihedral_zero.html
lammps_commands_improper.html
diff --git a/lib/atc/Function.cpp b/lib/atc/Function.cpp
index dfc42d514836c653d7d0866aa350248d6be9057d..6004937a018ffe58fcc6f4c0ace6cf916b881f83 100644
--- a/lib/atc/Function.cpp
+++ b/lib/atc/Function.cpp
@@ -1,4 +1,6 @@
+#ifndef _WIN32
#include <alloca.h>
+#endif
#include "Function.h"
#include "ATC_Error.h"
#include "LammpsInterface.h"
@@ -59,7 +61,11 @@ namespace ATC {
{
string type = args[0];
int narg = nargs -1;
- double *dargs = alloca(sizeof(double) * narg);
+#ifdef _WIN32
+ double *dargs = (double *) _alloca(sizeof(double) * narg);
+#endif
+ double *dargs = (double *) alloca(sizeof(double) * narg);
+#endif
for (int i = 0; i < narg; ++i) dargs[i] = atof(args[i+1]);
return function(type, narg, dargs);
@@ -193,7 +199,11 @@ XT_Function_Mgr * XT_Function_Mgr::myInstance_ = NULL;
{
string type = args[0];
int narg = nargs -1;
- double *dargs = alloca(sizeof(double) * narg);
+#ifdef _WIN32
+ double *dargs = (double *) _alloca(sizeof(double) * narg);
+#else
+ double *dargs = (double *) alloca(sizeof(double) * narg);
+#endif
for (int i = 0; i < narg; ++i) dargs[i] = atof(args[i+1]);
return function(type, narg, dargs);
diff --git a/lib/atc/Makefile.mpic++ b/lib/atc/Makefile.mpic++
new file mode 120000
index 0000000000000000000000000000000000000000..1f5f55d2ad24386005ef7961cca68ac1bd7f5d66
--- /dev/null
+++ b/lib/atc/Makefile.mpic++
@@ -0,0 +1 @@
+Makefile.mpi
\ No newline at end of file
diff --git a/lib/atc/SparseMatrix-inl.h b/lib/atc/SparseMatrix-inl.h
index 8d6aec78dca94f04503caf662761ff203d4e75e7..f70dd0361d3cf52e519fef06abc07f0518c86b21 100644
--- a/lib/atc/SparseMatrix-inl.h
+++ b/lib/atc/SparseMatrix-inl.h
@@ -67,9 +67,9 @@ void SparseMatrix<T>::_create(INDEX size, INDEX nrows)
// assign memory to hold matrix
try
{
- _val = (_size*nrows) ? new T [_size] : NULL;
- _ia = (_size*nrows) ? new INDEX [_nRowsCRS+1] : NULL;
- _ja = (_size*nrows) ? new INDEX [_size] : NULL;
+ _val = (_size && nrows) ? new T [_size] : NULL;
+ _ia = (_size && nrows) ? new INDEX [_nRowsCRS+1] : NULL;
+ _ja = (_size && nrows) ? new INDEX [_size] : NULL;
}
catch (std::exception &e)
{
diff --git a/src/.gitignore b/src/.gitignore
index 3e06a33580122b07bf88f6e5a0fb8bf992752073..3ae05079d211d2c043615cf93be2013449090506 100644
--- a/src/.gitignore
+++ b/src/.gitignore
@@ -324,6 +324,8 @@
/dihedral_spherical.h
/dihedral_table.cpp
/dihedral_table.h
+/dihedral_table_cut.cpp
+/dihedral_table_cut.h
/dump_atom_gz.cpp
/dump_atom_gz.h
/dump_xyz_gz.cpp
diff --git a/src/COLLOID/pair_lubricate.cpp b/src/COLLOID/pair_lubricate.cpp
index 5410a4c7731b493bd456f1a6e8658362a9ad6adf..f0f7ce8c143f259c298e053ae1b509f1531837ce 100644
--- a/src/COLLOID/pair_lubricate.cpp
+++ b/src/COLLOID/pair_lubricate.cpp
@@ -43,10 +43,6 @@
using namespace LAMMPS_NS;
using namespace MathConst;
-// same as fix_deform.cpp
-
-enum{NO_REMAP,X_REMAP,V_REMAP};
-
// same as fix_wall.cpp
enum{EDGE,CONSTANT,VARIABLE};
@@ -570,7 +566,7 @@ void PairLubricate::init_style()
for (int i = 0; i < modify->nfix; i++){
if (strcmp(modify->fix[i]->style,"deform") == 0) {
shearing = flagdeform = 1;
- if (((FixDeform *) modify->fix[i])->remapflag != V_REMAP)
+ if (((FixDeform *) modify->fix[i])->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using pair lubricate with inconsistent "
"fix deform remap option");
}
diff --git a/src/COLLOID/pair_lubricate_poly.cpp b/src/COLLOID/pair_lubricate_poly.cpp
index 0ed87c37be8834d051163d8da2666d5a860fdc19..5e5293336491ba2d78f15bea004eb6712ab76442 100644
--- a/src/COLLOID/pair_lubricate_poly.cpp
+++ b/src/COLLOID/pair_lubricate_poly.cpp
@@ -44,11 +44,6 @@
using namespace LAMMPS_NS;
using namespace MathConst;
-// same as fix_deform.cpp
-
-enum{NO_REMAP,X_REMAP,V_REMAP};
-
-
// same as fix_wall.cpp
enum{EDGE,CONSTANT,VARIABLE};
@@ -474,7 +469,7 @@ void PairLubricatePoly::init_style()
for (int i = 0; i < modify->nfix; i++){
if (strcmp(modify->fix[i]->style,"deform") == 0) {
shearing = flagdeform = 1;
- if (((FixDeform *) modify->fix[i])->remapflag != V_REMAP)
+ if (((FixDeform *) modify->fix[i])->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using pair lubricate with inconsistent "
"fix deform remap option");
}
@@ -550,7 +545,7 @@ void PairLubricatePoly::init_style()
for (int i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
shearing = 1;
- if (((FixDeform *) modify->fix[i])->remapflag != V_REMAP)
+ if (((FixDeform *) modify->fix[i])->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using pair lubricate/poly with inconsistent "
"fix deform remap option");
}
diff --git a/src/COMPRESS/dump_atom_gz.cpp b/src/COMPRESS/dump_atom_gz.cpp
index fb605b74cee56bae62e53ab4f0251f7e612d3aa1..ef7e6583be5e28d1af258aa636fa6c509816bc47 100644
--- a/src/COMPRESS/dump_atom_gz.cpp
+++ b/src/COMPRESS/dump_atom_gz.cpp
@@ -71,6 +71,19 @@ void DumpAtomGZ::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
// each proc with filewriter = 1 opens a file
diff --git a/src/COMPRESS/dump_cfg_gz.cpp b/src/COMPRESS/dump_cfg_gz.cpp
index a1ffa92ec859958c616a310c7705c088f35114c0..aaeb878584b0cca77110b1a7bbef58d8de6f6cab 100644
--- a/src/COMPRESS/dump_cfg_gz.cpp
+++ b/src/COMPRESS/dump_cfg_gz.cpp
@@ -75,6 +75,19 @@ void DumpCFGGZ::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
// each proc with filewriter = 1 opens a file
diff --git a/src/COMPRESS/dump_custom_gz.cpp b/src/COMPRESS/dump_custom_gz.cpp
index 52d67445be0ca940fdafc3b08ef8ba9bb93c51ab..9c30f4742f92de2dcf92c2c4c2e22f28b9c1ca02 100644
--- a/src/COMPRESS/dump_custom_gz.cpp
+++ b/src/COMPRESS/dump_custom_gz.cpp
@@ -73,6 +73,19 @@ void DumpCustomGZ::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
// each proc with filewriter = 1 opens a file
diff --git a/src/COMPRESS/dump_xyz_gz.cpp b/src/COMPRESS/dump_xyz_gz.cpp
index 047ad8652f9587bdba8d306be030a357711ccc78..7be1a10fe2b15b19daa5f3b93d0659d75959922f 100644
--- a/src/COMPRESS/dump_xyz_gz.cpp
+++ b/src/COMPRESS/dump_xyz_gz.cpp
@@ -73,6 +73,19 @@ void DumpXYZGZ::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
// each proc with filewriter = 1 opens a file
diff --git a/src/GRANULAR/fix_pour.cpp b/src/GRANULAR/fix_pour.cpp
index 835d6c42ddd49aff7d3d24c07b9729ac3b6060ac..318c54c99a8844bea5f3f86d92c9ebd67f23abd3 100644
--- a/src/GRANULAR/fix_pour.cpp
+++ b/src/GRANULAR/fix_pour.cpp
@@ -39,7 +39,6 @@ using namespace MathConst;
enum{ATOM,MOLECULE};
enum{ONE,RANGE,POLY};
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
#define EPSILON 0.001
#define SMALL 1.0e-10
@@ -609,7 +608,7 @@ void FixPour::pre_exchange()
newcoord[1] >= sublo[1] && newcoord[1] < subhi[1] &&
newcoord[2] >= sublo[2] && newcoord[2] < subhi[2]) flag = 1;
else if (dimension == 3 && newcoord[2] >= domain->boxhi[2]) {
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (comm->myloc[2] == comm->procgrid[2]-1 &&
newcoord[0] >= sublo[0] && newcoord[0] < subhi[0] &&
newcoord[1] >= sublo[1] && newcoord[1] < subhi[1]) flag = 1;
@@ -619,7 +618,7 @@ void FixPour::pre_exchange()
newcoord[1] >= sublo[1] && newcoord[1] < subhi[1]) flag = 1;
}
} else if (dimension == 2 && newcoord[1] >= domain->boxhi[1]) {
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (comm->myloc[1] == comm->procgrid[1]-1 &&
newcoord[0] >= sublo[0] && newcoord[0] < subhi[0]) flag = 1;
} else {
diff --git a/src/KOKKOS/comm_kokkos.cpp b/src/KOKKOS/comm_kokkos.cpp
index d24fa8a536d56cbb550d5b7e547b32aef57bef0b..96d1c64d0d1aa389caca8c7ef8a635dd406fc2e9 100644
--- a/src/KOKKOS/comm_kokkos.cpp
+++ b/src/KOKKOS/comm_kokkos.cpp
@@ -36,8 +36,6 @@ using namespace LAMMPS_NS;
#define BUFMIN 10000
#define BUFEXTRA 1000
-enum{SINGLE,MULTI};
-
/* ----------------------------------------------------------------------
setup MPI and allocate buffer space
------------------------------------------------------------------------- */
@@ -709,7 +707,7 @@ void CommKokkos::borders()
if (!exchange_comm_classic) {
static int print = 1;
- if (style != SINGLE || bordergroup || ghost_velocity) {
+ if (style != Comm::SINGLE || bordergroup || ghost_velocity) {
if (print && comm->me==0) {
error->warning(FLERR,"Required border comm not yet implemented in Kokkos communication, "
"switching to classic communication");
@@ -817,7 +815,7 @@ void CommKokkos::borders_device() {
// store sent atom indices in list for use in future timesteps
x = atom->x;
- if (style == SINGLE) {
+ if (style == Comm::SINGLE) {
lo = slablo[iswap];
hi = slabhi[iswap];
} else {
@@ -846,7 +844,7 @@ void CommKokkos::borders_device() {
if (sendflag) {
if (!bordergroup || ineed >= 2) {
- if (style == SINGLE) {
+ if (style == Comm::SINGLE) {
k_total_send.h_view() = 0;
k_total_send.template modify<LMPHostType>();
k_total_send.template sync<LMPDeviceType>();
@@ -894,7 +892,7 @@ void CommKokkos::borders_device() {
} else {
error->all(FLERR,"Required border comm not yet "
"implemented with Kokkos");
- if (style == SINGLE) {
+ if (style == Comm::SINGLE) {
ngroup = atom->nfirst;
for (i = 0; i < ngroup; i++)
if (x[i][dim] >= lo && x[i][dim] <= hi) {
@@ -1099,7 +1097,7 @@ void CommKokkos::grow_swap(int n)
{
free_swap();
allocate_swap(n);
- if (style == MULTI) {
+ if (style == Comm::MULTI) {
free_multi();
allocate_multi(n);
}
diff --git a/src/KOKKOS/comm_tiled_kokkos.cpp b/src/KOKKOS/comm_tiled_kokkos.cpp
index 4809f4cff2eebae63bcd002d0199d8e2fba1c0bc..81cf1f0563ac799f801a2937afa481227c7d63c4 100644
--- a/src/KOKKOS/comm_tiled_kokkos.cpp
+++ b/src/KOKKOS/comm_tiled_kokkos.cpp
@@ -39,9 +39,6 @@ using namespace LAMMPS_NS;
#define DELTA_PROCS 16
-enum{SINGLE,MULTI}; // same as in Comm
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
/* ---------------------------------------------------------------------- */
CommTiledKokkos::CommTiledKokkos(LAMMPS *lmp) : CommTiled(lmp)
diff --git a/src/KOKKOS/fix_deform_kokkos.cpp b/src/KOKKOS/fix_deform_kokkos.cpp
index 050364979550f25b74be0961ce5930ba4aa62f1d..7b9336ed3b03822e78cfa91b4e1e920be550a02d 100644
--- a/src/KOKKOS/fix_deform_kokkos.cpp
+++ b/src/KOKKOS/fix_deform_kokkos.cpp
@@ -41,10 +41,6 @@ using namespace MathConst;
enum{NONE,FINAL,DELTA,SCALE,VEL,ERATE,TRATE,VOLUME,WIGGLE,VARIABLE};
enum{ONE_FROM_ONE,ONE_FROM_TWO,TWO_FROM_ONE};
-// same as domain.cpp, fix_nvt_sllod.cpp, compute_temp_deform.cpp
-
-enum{NO_REMAP,X_REMAP,V_REMAP};
-
/* ---------------------------------------------------------------------- */
FixDeformKokkos::FixDeformKokkos(LAMMPS *lmp, int narg, char **arg) : FixDeform(lmp, narg, arg)
@@ -315,7 +311,7 @@ void FixDeformKokkos::end_of_step()
// convert atoms and rigid bodies to lamda coords
- if (remapflag == X_REMAP) {
+ if (remapflag == Domain::X_REMAP) {
int nlocal = atom->nlocal;
domainKK->x2lamda(nlocal);
@@ -355,7 +351,7 @@ void FixDeformKokkos::end_of_step()
// convert atoms and rigid bodies back to box coords
- if (remapflag == X_REMAP) {
+ if (remapflag == Domain::X_REMAP) {
int nlocal = atom->nlocal;
domainKK->lamda2x(nlocal);
diff --git a/src/KOKKOS/nbin_kokkos.cpp b/src/KOKKOS/nbin_kokkos.cpp
index 2ca2a83b5ae5a7284ec5c65af956fdcf9fc131c3..c35c34967517291b67407a3f78f81773d7e9fe68 100644
--- a/src/KOKKOS/nbin_kokkos.cpp
+++ b/src/KOKKOS/nbin_kokkos.cpp
@@ -23,8 +23,6 @@
using namespace LAMMPS_NS;
-enum{NSQ,BIN,MULTI}; // also in Neighbor
-
#define SMALL 1.0e-6
#define CUT2BIN_RATIO 100
diff --git a/src/KOKKOS/neigh_bond_kokkos.cpp b/src/KOKKOS/neigh_bond_kokkos.cpp
index 4dd39a1d383cb5da9063a5c439ca725d511d7c48..615dc94ee2faf6d2c690557909757b9a0525181d 100644
--- a/src/KOKKOS/neigh_bond_kokkos.cpp
+++ b/src/KOKKOS/neigh_bond_kokkos.cpp
@@ -37,8 +37,6 @@ using namespace LAMMPS_NS;
#define BONDDELTA 10000
#define LB_FACTOR 1.5
-enum{IGNORE,WARN,ERROR}; // same as thermo.cpp
-
/* ---------------------------------------------------------------------- */
template<class DeviceType>
@@ -288,7 +286,7 @@ void NeighBondKokkos<DeviceType>::bond_all()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Bond atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -296,7 +294,7 @@ void NeighBondKokkos<DeviceType>::bond_all()
}
if (neighbor->cluster_check) bond_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -317,7 +315,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondBondAll, const int &i,
int atom1 = map_array(bond_atom(i,m));
if (atom1 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
@@ -383,7 +381,7 @@ void NeighBondKokkos<DeviceType>::bond_partial()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Bond atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -391,7 +389,7 @@ void NeighBondKokkos<DeviceType>::bond_partial()
}
if (neighbor->cluster_check) bond_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -413,7 +411,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondBondPartial, const int
int atom1 = map_array(bond_atom(i,m));
if (atom1 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
@@ -505,7 +503,7 @@ void NeighBondKokkos<DeviceType>::angle_all()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Angle atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -513,7 +511,7 @@ void NeighBondKokkos<DeviceType>::angle_all()
}
if (neighbor->cluster_check) angle_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -536,7 +534,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondAngleAll, const int &i,
int atom3 = map_array(angle_atom3(i,m));
if (atom1 == -1 || atom2 == -1 || atom3 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
@@ -607,7 +605,7 @@ void NeighBondKokkos<DeviceType>::angle_partial()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Angle atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -615,7 +613,7 @@ void NeighBondKokkos<DeviceType>::angle_partial()
}
if (neighbor->cluster_check) angle_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -639,7 +637,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondAnglePartial, const int
int atom3 = map_array(angle_atom3(i,m));
if (atom1 == -1 || atom2 == -1 || atom3 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
@@ -749,7 +747,7 @@ void NeighBondKokkos<DeviceType>::dihedral_all()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Dihedral atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -757,7 +755,7 @@ void NeighBondKokkos<DeviceType>::dihedral_all()
}
if (neighbor->cluster_check) dihedral_check(neighbor->ndihedrallist,v_dihedrallist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -781,7 +779,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondDihedralAll, const int
int atom4 = map_array(dihedral_atom4(i,m));
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
@@ -856,7 +854,7 @@ void NeighBondKokkos<DeviceType>::dihedral_partial()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Dihedral atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -864,7 +862,7 @@ void NeighBondKokkos<DeviceType>::dihedral_partial()
}
if (neighbor->cluster_check) dihedral_check(neighbor->ndihedrallist,v_dihedrallist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -889,7 +887,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondDihedralPartial, const
int atom4 = map_array(dihedral_atom4(i,m));
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
@@ -1020,7 +1018,7 @@ void NeighBondKokkos<DeviceType>::improper_all()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Improper atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -1028,7 +1026,7 @@ void NeighBondKokkos<DeviceType>::improper_all()
}
if (neighbor->cluster_check) dihedral_check(neighbor->nimproperlist,v_improperlist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -1052,7 +1050,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondImproperAll, const int
int atom4 = map_array(improper_atom4(i,m));
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
@@ -1127,7 +1125,7 @@ void NeighBondKokkos<DeviceType>::improper_partial()
}
} while (h_fail_flag());
- if (nmissing && lostbond == ERROR) {
+ if (nmissing && lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Improper atoms missing on proc %d at step " BIGINT_FORMAT,
me,update->ntimestep);
@@ -1135,7 +1133,7 @@ void NeighBondKokkos<DeviceType>::improper_partial()
}
if (neighbor->cluster_check) dihedral_check(neighbor->nimproperlist,v_improperlist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
@@ -1160,7 +1158,7 @@ void NeighBondKokkos<DeviceType>::operator()(TagNeighBondImproperPartial, const
int atom4 = map_array(improper_atom4(i,m));
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) return;
+ if (lostbond == Thermo::ERROR) return;
continue;
}
atom1 = closest_image(i,atom1);
diff --git a/src/KOKKOS/neigh_list_kokkos.cpp b/src/KOKKOS/neigh_list_kokkos.cpp
index 98294a802aaa88e16a2c1ca8269cf235d149ff3a..93cf0af937a5e695e8f19f1d6bb18a12032deab7 100644
--- a/src/KOKKOS/neigh_list_kokkos.cpp
+++ b/src/KOKKOS/neigh_list_kokkos.cpp
@@ -17,8 +17,6 @@
using namespace LAMMPS_NS;
-enum{NSQ,BIN,MULTI};
-
/* ---------------------------------------------------------------------- */
template<class Device>
diff --git a/src/KOKKOS/neighbor_kokkos.cpp b/src/KOKKOS/neighbor_kokkos.cpp
index 72092aef0b17bfffaea8865c9156f010f9f829c0..6cd5c4e6e597c93c306156550bc7bd73787d6c28 100644
--- a/src/KOKKOS/neighbor_kokkos.cpp
+++ b/src/KOKKOS/neighbor_kokkos.cpp
@@ -33,8 +33,6 @@
using namespace LAMMPS_NS;
-enum{NSQ,BIN,MULTI}; // also in neigh_list.cpp
-
/* ---------------------------------------------------------------------- */
NeighborKokkos::NeighborKokkos(LAMMPS *lmp) : Neighbor(lmp),
@@ -98,7 +96,7 @@ void NeighborKokkos::init_cutneighsq_kokkos(int n)
void NeighborKokkos::create_kokkos_list(int i)
{
- if (style != BIN)
+ if (style != Neighbor::BIN)
error->all(FLERR,"KOKKOS package only supports 'bin' neighbor lists");
if (requests[i]->kokkos_device) {
@@ -300,7 +298,7 @@ void NeighborKokkos::build_kokkos(int topoflag)
// if bin then, atoms may have moved outside of proc domain & bin extent,
// leading to errors or even a crash
- if (style != NSQ) {
+ if (style != Neighbor::NSQ) {
for (int i = 0; i < nbin; i++) {
neigh_bin[i]->bin_atoms_setup(nall);
neigh_bin[i]->bin_atoms();
diff --git a/src/KSPACE/ewald_disp.cpp b/src/KSPACE/ewald_disp.cpp
index fb48edc117ae48ce986ba9f794bcc58bfe0b062f..8cdaa044aa737917ebfa9c04ff26391cbb2eb3eb 100644
--- a/src/KSPACE/ewald_disp.cpp
+++ b/src/KSPACE/ewald_disp.cpp
@@ -39,8 +39,6 @@ using namespace MathSpecial;
#define SMALL 0.00001
-enum{GEOMETRIC,ARITHMETIC,SIXTHPOWER}; // same as in pair.h
-
//#define DEBUG
/* ---------------------------------------------------------------------- */
@@ -116,7 +114,7 @@ void EwaldDisp::init()
if (!(ptr||cutoff))
error->all(FLERR,"KSpace style is incompatible with Pair style");
int ewald_order = ptr ? *((int *) ptr) : 1<<1;
- int ewald_mix = ptr ? *((int *) pair->extract("ewald_mix",tmp)) : GEOMETRIC;
+ int ewald_mix = ptr ? *((int *) pair->extract("ewald_mix",tmp)) : Pair::GEOMETRIC;
memset(function, 0, EWALD_NFUNCS*sizeof(int));
for (int i=0; i<=EWALD_NORDER; ++i) // transcribe order
if (ewald_order&(1<<i)) { // from pair_style
@@ -127,8 +125,8 @@ void EwaldDisp::init()
case 3:
k = 3; break;
case 6:
- if (ewald_mix==GEOMETRIC) { k = 1; break; }
- else if (ewald_mix==ARITHMETIC) { k = 2; break; }
+ if (ewald_mix==Pair::GEOMETRIC) { k = 1; break; }
+ else if (ewald_mix==Pair::ARITHMETIC) { k = 2; break; }
error->all(FLERR,
"Unsupported mixing rule in kspace_style ewald/disp");
default:
diff --git a/src/KSPACE/pppm_disp.cpp b/src/KSPACE/pppm_disp.cpp
index 4fbc20c0fc80d8d338308c22724a49a5ff3fbdb9..aea57da326e7f75b15a23ad3945a552f3eb56be9 100644
--- a/src/KSPACE/pppm_disp.cpp
+++ b/src/KSPACE/pppm_disp.cpp
@@ -46,7 +46,6 @@ using namespace MathConst;
#define LARGE 10000.0
#define EPS_HOC 1.0e-7
-enum{GEOMETRIC,ARITHMETIC,SIXTHPOWER};
enum{REVERSE_RHO, REVERSE_RHO_G, REVERSE_RHO_A, REVERSE_RHO_NONE};
enum{FORWARD_IK, FORWARD_AD, FORWARD_IK_PERATOM, FORWARD_AD_PERATOM,
FORWARD_IK_G, FORWARD_AD_G, FORWARD_IK_PERATOM_G, FORWARD_AD_PERATOM_G,
@@ -304,7 +303,7 @@ void PPPMDisp::init()
// check out which types of potentials will have to be calculated
int ewald_order = ptr ? *((int *) ptr) : 1<<1;
- int ewald_mix = ptr ? *((int *) pair->extract("ewald_mix",tmp)) : GEOMETRIC;
+ int ewald_mix = ptr ? *((int *) pair->extract("ewald_mix",tmp)) : Pair::GEOMETRIC;
memset(function, 0, EWALD_FUNCS*sizeof(int));
for (int i=0; i<=EWALD_MAXORDER; ++i) // transcribe order
if (ewald_order&(1<<i)) { // from pair_style
@@ -314,9 +313,9 @@ void PPPMDisp::init()
case 1:
k = 0; break;
case 6:
- if ((ewald_mix==GEOMETRIC || ewald_mix==SIXTHPOWER ||
+ if ((ewald_mix==Pair::GEOMETRIC || ewald_mix==Pair::SIXTHPOWER ||
mixflag == 1) && mixflag!= 2) { k = 1; break; }
- else if (ewald_mix==ARITHMETIC && mixflag!=2) { k = 2; break; }
+ else if (ewald_mix==Pair::ARITHMETIC && mixflag!=2) { k = 2; break; }
else if (mixflag == 2) { k = 3; break; }
default:
sprintf(str, "Unsupported order in kspace_style "
diff --git a/src/MISC/fix_deposit.cpp b/src/MISC/fix_deposit.cpp
index 4215df07a62c3ef7d1f8e276ef9a8015bf581d41..1428299eb7d8366446cff2323f7c5f7496311f5e 100644
--- a/src/MISC/fix_deposit.cpp
+++ b/src/MISC/fix_deposit.cpp
@@ -37,7 +37,6 @@ using namespace FixConst;
using namespace MathConst;
enum{ATOM,MOLECULE};
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
enum{DIST_UNIFORM,DIST_GAUSSIAN};
#define EPSILON 1.0e6
@@ -506,7 +505,7 @@ void FixDeposit::pre_exchange()
newcoord[1] >= sublo[1] && newcoord[1] < subhi[1] &&
newcoord[2] >= sublo[2] && newcoord[2] < subhi[2]) flag = 1;
else if (dimension == 3 && newcoord[2] >= domain->boxhi[2]) {
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (comm->myloc[2] == comm->procgrid[2]-1 &&
newcoord[0] >= sublo[0] && newcoord[0] < subhi[0] &&
newcoord[1] >= sublo[1] && newcoord[1] < subhi[1]) flag = 1;
@@ -516,7 +515,7 @@ void FixDeposit::pre_exchange()
newcoord[1] >= sublo[1] && newcoord[1] < subhi[1]) flag = 1;
}
} else if (dimension == 2 && newcoord[1] >= domain->boxhi[1]) {
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (comm->myloc[1] == comm->procgrid[1]-1 &&
newcoord[0] >= sublo[0] && newcoord[0] < subhi[0]) flag = 1;
} else {
diff --git a/src/MPIIO/dump_atom_mpiio.cpp b/src/MPIIO/dump_atom_mpiio.cpp
index 4422ef1cdfb5c16df75cdd8ff745d67a212d889d..f2f299144261d9d6035e49db8e1f30dedeea38ac 100644
--- a/src/MPIIO/dump_atom_mpiio.cpp
+++ b/src/MPIIO/dump_atom_mpiio.cpp
@@ -77,6 +77,19 @@ void DumpAtomMPIIO::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
if (append_flag) { // append open
diff --git a/src/MPIIO/dump_cfg_mpiio.cpp b/src/MPIIO/dump_cfg_mpiio.cpp
index 5086d1dc6cdf70cd93ad2f44760ca3c8795fc64b..c580d2e7af71cdb1b7ca037392e0e0656b6ac5b0 100644
--- a/src/MPIIO/dump_cfg_mpiio.cpp
+++ b/src/MPIIO/dump_cfg_mpiio.cpp
@@ -47,8 +47,6 @@ using namespace LAMMPS_NS;
#define DUMP_BUF_CHUNK_SIZE 16384
#define DUMP_BUF_INCREMENT_SIZE 4096
-enum{INT,DOUBLE,STRING,BIGINT}; // same as in DumpCustom
-
#define UNWRAPEXPAND 10.0
#define ONEFIELD 32
#define DELTA 1048576
@@ -96,6 +94,19 @@ void DumpCFGMPIIO::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
if (append_flag) { // append open
@@ -377,13 +388,13 @@ int DumpCFGMPIIO::convert_string_omp(int n, double *mybuf)
} else if (j == 1) {
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),"%s \n",typenames[(int) mybuf[bufOffset[tid]+m]]);
} else if (j >= 2) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],static_cast<int> (mybuf[bufOffset[tid]+m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],mybuf[bufOffset[tid]+m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],typenames[(int) mybuf[bufOffset[tid]+m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],static_cast<bigint> (mybuf[bufOffset[tid]+m]));
}
m++;
@@ -410,18 +421,18 @@ int DumpCFGMPIIO::convert_string_omp(int n, double *mybuf)
//offset += sprintf(&sbuf[offset],vformat[j],unwrap_coord);
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],unwrap_coord);
} else if (j >= 5 ) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
//offset +=
// sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],static_cast<int> (mybuf[bufOffset[tid]+m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
// offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],mybuf[bufOffset[tid]+m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
// offset +=
// sprintf(&sbuf[offset],vformat[j],typenames[(int) mybuf[m]]);
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],typenames[(int) mybuf[bufOffset[tid]+m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
// offset +=
// sprintf(&sbuf[offset],vformat[j],static_cast<bigint> (mybuf[m]));
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],static_cast<bigint> (mybuf[bufOffset[tid]+m]));
diff --git a/src/MPIIO/dump_custom_mpiio.cpp b/src/MPIIO/dump_custom_mpiio.cpp
index 2da6dd700d1348a15a96f68c14d6363b0a04ce14..3650ca994e4d519d319ae0d182cb9b2011afd764 100644
--- a/src/MPIIO/dump_custom_mpiio.cpp
+++ b/src/MPIIO/dump_custom_mpiio.cpp
@@ -53,7 +53,6 @@ enum{ID,MOL,TYPE,ELEMENT,MASS,
TQX,TQY,TQZ,SPIN,ERADIUS,ERVEL,ERFORCE,
COMPUTE,FIX,VARIABLE};
enum{LT,LE,GT,GE,EQ,NEQ};
-enum{INT,DOUBLE,STRING,BIGINT}; // same as in DumpCustom
/* ---------------------------------------------------------------------- */
@@ -98,6 +97,19 @@ void DumpCustomMPIIO::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
if (append_flag) { // append open
@@ -246,13 +258,13 @@ void DumpCustomMPIIO::init_style()
if (format_column_user[i]) {
vformat[i] = new char[strlen(format_column_user[i]) + 2];
strcpy(vformat[i],format_column_user[i]);
- } else if (vtype[i] == INT && format_int_user) {
+ } else if (vtype[i] == Dump::INT && format_int_user) {
vformat[i] = new char[strlen(format_int_user) + 2];
strcpy(vformat[i],format_int_user);
- } else if (vtype[i] == DOUBLE && format_float_user) {
+ } else if (vtype[i] == Dump::DOUBLE && format_float_user) {
vformat[i] = new char[strlen(format_float_user) + 2];
strcpy(vformat[i],format_float_user);
- } else if (vtype[i] == BIGINT && format_bigint_user) {
+ } else if (vtype[i] == Dump::BIGINT && format_bigint_user) {
vformat[i] = new char[strlen(format_bigint_user) + 2];
strcpy(vformat[i],format_bigint_user);
} else {
@@ -618,11 +630,11 @@ int DumpCustomMPIIO::convert_string_omp(int n, double *mybuf)
}
for (int j = 0; j < size_one; j++) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],static_cast<int> (mybuf[bufOffset[tid]+m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],mybuf[bufOffset[tid]+m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
mpifhStringCountPerThread[tid] += sprintf(&(mpifh_buffer_line_per_thread[tid][mpifhStringCountPerThread[tid]]),vformat[j],typenames[(int) mybuf[bufOffset[tid]+m]]);
m ++;
}
diff --git a/src/MPIIO/dump_xyz_mpiio.cpp b/src/MPIIO/dump_xyz_mpiio.cpp
index 4725efccd14b990d28ea7afe6991c40117f5d688..f15d340cd14f3f470ac981e80d29c7428c2f6e09 100644
--- a/src/MPIIO/dump_xyz_mpiio.cpp
+++ b/src/MPIIO/dump_xyz_mpiio.cpp
@@ -53,7 +53,6 @@ enum{ID,MOL,TYPE,ELEMENT,MASS,
TQX,TQY,TQZ,SPIN,ERADIUS,ERVEL,ERFORCE,
COMPUTE,FIX,VARIABLE};
enum{LT,LE,GT,GE,EQ,NEQ};
-enum{INT,DOUBLE,STRING}; // same as in DumpCFG
/* ---------------------------------------------------------------------- */
@@ -98,6 +97,19 @@ void DumpXYZMPIIO::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
if (append_flag) { // append open
diff --git a/src/POEMS/fix_poems.cpp b/src/POEMS/fix_poems.cpp
index edc72ea70ce953ba2fc28d21efbedcb3fd49f00c..2c2fefb250b14252fa31e88e489ea9f307e918a5 100644
--- a/src/POEMS/fix_poems.cpp
+++ b/src/POEMS/fix_poems.cpp
@@ -369,7 +369,7 @@ void FixPOEMS::init()
if (pflag && (modify->fmask[i] & POST_FORCE) &&
!modify->fix[i]->rigid_flag) {
char str[128];
- sprintf(str,"Fix %d alters forces after fix poems",modify->fix[i]->id);
+ sprintf(str,"Fix %s alters forces after fix poems",modify->fix[i]->id);
error->warning(FLERR,str);
}
}
diff --git a/src/RIGID/fix_rigid.cpp b/src/RIGID/fix_rigid.cpp
index 07ebf884c7450cf53c4d6e028f758ace69499100..7104d3084994d2c7d9c4b4ff5d89d8aa387d6c00 100644
--- a/src/RIGID/fix_rigid.cpp
+++ b/src/RIGID/fix_rigid.cpp
@@ -730,7 +730,7 @@ void FixRigid::init()
if (rflag && (modify->fmask[i] & POST_FORCE) &&
!modify->fix[i]->rigid_flag) {
char str[128];
- sprintf(str,"Fix %d alters forces after fix rigid",modify->fix[i]->id);
+ sprintf(str,"Fix %s alters forces after fix rigid",modify->fix[i]->id);
error->warning(FLERR,str);
}
}
diff --git a/src/RIGID/fix_rigid_small.cpp b/src/RIGID/fix_rigid_small.cpp
index d0cb596c82c8736af9a0e92457224e021c0c7f41..c1de89efc70da5221fd8c0ee6c75a0ec5312d450 100644
--- a/src/RIGID/fix_rigid_small.cpp
+++ b/src/RIGID/fix_rigid_small.cpp
@@ -569,7 +569,7 @@ void FixRigidSmall::init()
if (rflag && (modify->fmask[i] & POST_FORCE) &&
!modify->fix[i]->rigid_flag) {
char str[128];
- sprintf(str,"Fix %d alters forces after fix rigid",modify->fix[i]->id);
+ sprintf(str,"Fix %s alters forces after fix rigid",modify->fix[i]->id);
error->warning(FLERR,str);
}
}
diff --git a/src/SHOCK/fix_append_atoms.cpp b/src/SHOCK/fix_append_atoms.cpp
index 27a7f5514abcfe821ddc88a416dbc6d2037f2c90..d898d2a790d62a03c8451eaf9c564fe52334f3f1 100644
--- a/src/SHOCK/fix_append_atoms.cpp
+++ b/src/SHOCK/fix_append_atoms.cpp
@@ -29,8 +29,6 @@
using namespace LAMMPS_NS;
using namespace FixConst;
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
#define BIG 1.0e30
#define EPSILON 1.0e-6
@@ -417,7 +415,7 @@ void FixAppendAtoms::pre_exchange()
if (spatflag==1) if (get_spatial()==0) return;
int addflag = 0;
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (comm->myloc[2] == comm->procgrid[2]-1) addflag = 1;
} else {
if (comm->mysplit[2][1] == 1.0) addflag = 1;
diff --git a/src/SRD/fix_srd.cpp b/src/SRD/fix_srd.cpp
index 6bc72dca4bd72974cd39a054263364b51bab9dc8..564c1e98d7362a5947bc281a8b766ad5b6f4e479 100644
--- a/src/SRD/fix_srd.cpp
+++ b/src/SRD/fix_srd.cpp
@@ -52,8 +52,6 @@ enum{BIG_MOVE,SRD_MOVE,SRD_ROTATE};
enum{CUBIC_ERROR,CUBIC_WARN};
enum{SHIFT_NO,SHIFT_YES,SHIFT_POSSIBLE};
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-
#define EINERTIA 0.2 // moment of inertia prefactor for ellipsoid
#define ATOMPERBIN 30
@@ -384,7 +382,7 @@ void FixSRD::init()
if (strcmp(modify->fix[i]->style,"deform") == 0) {
deformflag = 1;
FixDeform *deform = (FixDeform *) modify->fix[i];
- if (deform->box_change_shape && deform->remapflag != V_REMAP)
+ if (deform->box_change_shape && deform->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using fix srd with inconsistent "
"fix deform remap option");
}
diff --git a/src/USER-EFF/compute_temp_deform_eff.cpp b/src/USER-EFF/compute_temp_deform_eff.cpp
index 53fb7683a76bcde71ecf5dcfd32270c48c877844..d3856e431d0d427e1f90c3435f99e63c37097c60 100644
--- a/src/USER-EFF/compute_temp_deform_eff.cpp
+++ b/src/USER-EFF/compute_temp_deform_eff.cpp
@@ -35,8 +35,6 @@
using namespace LAMMPS_NS;
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-
/* ---------------------------------------------------------------------- */
ComputeTempDeformEff::ComputeTempDeformEff(LAMMPS *lmp, int narg, char **arg) :
@@ -76,7 +74,7 @@ void ComputeTempDeformEff::init()
int i;
for (i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
- if (((FixDeform *) modify->fix[i])->remapflag == X_REMAP &&
+ if (((FixDeform *) modify->fix[i])->remapflag == Domain::X_REMAP &&
comm->me == 0)
error->warning(FLERR,"Using compute temp/deform/eff with inconsistent "
"fix deform remap option");
diff --git a/src/USER-EFF/fix_nvt_sllod_eff.cpp b/src/USER-EFF/fix_nvt_sllod_eff.cpp
index f7b1a4ab2b162925999df9797050993902491772..13cb34048c105d0e5a7f5b753bafbdcdab044eba 100644
--- a/src/USER-EFF/fix_nvt_sllod_eff.cpp
+++ b/src/USER-EFF/fix_nvt_sllod_eff.cpp
@@ -27,8 +27,6 @@
using namespace LAMMPS_NS;
using namespace FixConst;
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-
/* ---------------------------------------------------------------------- */
FixNVTSllodEff::FixNVTSllodEff(LAMMPS *lmp, int narg, char **arg) :
@@ -78,7 +76,7 @@ void FixNVTSllodEff::init()
int i;
for (i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
- if (((FixDeform *) modify->fix[i])->remapflag != V_REMAP)
+ if (((FixDeform *) modify->fix[i])->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using fix nvt/sllod/eff with inconsistent fix deform "
"remap option");
break;
diff --git a/src/USER-INTEL/fix_intel.cpp b/src/USER-INTEL/fix_intel.cpp
index 124c7a9c76f403ea29dcee43a368e1ad97e3774f..7a206f472a92eed49e814e9029d69a936b196daa 100644
--- a/src/USER-INTEL/fix_intel.cpp
+++ b/src/USER-INTEL/fix_intel.cpp
@@ -55,8 +55,6 @@ using namespace FixConst;
#endif
#endif
-enum{NSQ,BIN,MULTI};
-
/* ---------------------------------------------------------------------- */
FixIntel::FixIntel(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
@@ -353,7 +351,7 @@ void FixIntel::init()
void FixIntel::setup(int vflag)
{
- if (neighbor->style != BIN)
+ if (neighbor->style != Neighbor::BIN)
error->all(FLERR,
"Currently, neighbor style BIN must be used with Intel package.");
if (neighbor->exclude_setting() != 0)
diff --git a/src/USER-INTEL/fix_nvt_sllod_intel.cpp b/src/USER-INTEL/fix_nvt_sllod_intel.cpp
index 5684f469038d023d71cbc4dd95f79cd1b9907601..ad7e1384f06bd14d90dc241beb6ed1d25f465f49 100644
--- a/src/USER-INTEL/fix_nvt_sllod_intel.cpp
+++ b/src/USER-INTEL/fix_nvt_sllod_intel.cpp
@@ -27,8 +27,6 @@
using namespace LAMMPS_NS;
using namespace FixConst;
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-
/* ---------------------------------------------------------------------- */
FixNVTSllodIntel::FixNVTSllodIntel(LAMMPS *lmp, int narg, char **arg) :
@@ -79,7 +77,7 @@ void FixNVTSllodIntel::init()
int i;
for (i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
- if (((FixDeform *) modify->fix[i])->remapflag != V_REMAP)
+ if (((FixDeform *) modify->fix[i])->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using fix nvt/sllod with inconsistent fix deform "
"remap option");
break;
diff --git a/src/USER-MISC/README b/src/USER-MISC/README
index 7fc931b9625ed9fa6884bb5b1564692cf77564dd..222c06b4c426ff1adf271dcccf7438921eb5e5ff 100644
--- a/src/USER-MISC/README
+++ b/src/USER-MISC/README
@@ -37,6 +37,7 @@ dihedral_style nharmonic, Loukas Peristeras, loukas.peristeras at scienomics.com
dihedral_style quadratic, Loukas Peristeras, loukas.peristeras at scienomics.com, 27 Oct 12
dihedral_style spherical, Andrew Jewett, jewett.aij@gmail.com, 15 Jul 16
dihedral_style table, Andrew Jewett, jewett.aij@gmail.com, 10 Jan 12
+dihedral_style table/cut, Mike Salerno, ksalerno@pha.jhu.edu, 11 May 18
fix addtorque, Laurent Joly (U Lyon), ljoly.ulyon at gmail.com, 8 Aug 11
fix ave/correlate/long, Jorge Ramirez (UPM Madrid), jorge.ramirez at upm.es, 21 Oct 2015
fix bond/react, Jacob Gissinger (CU Boulder), info at disarmmd.org, 24 Feb 2018
diff --git a/src/USER-MISC/dihedral_table_cut.cpp b/src/USER-MISC/dihedral_table_cut.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..7345bcdad164c2df0cd598c8606b8f26cf1fc511
--- /dev/null
+++ b/src/USER-MISC/dihedral_table_cut.cpp
@@ -0,0 +1,1500 @@
+/* ----------------------------------------------------------------------
+ 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.
+------------------------------------------------------------------------- */
+
+/* ----------------------------------------------------------------------
+ Contributing author: K. Michael Salerno (NRL)
+ Based on tabulated dihedral (dihedral_table.cpp) by Andrew Jewett
+------------------------------------------------------------------------- */
+
+#include <cmath>
+#include <cstdlib>
+#include <cstring>
+#include <cassert>
+#include <string>
+#include <fstream>
+#include <iostream>
+#include <sstream>
+
+#include "dihedral_table_cut.h"
+#include "atom.h"
+#include "neighbor.h"
+#include "update.h"
+#include "domain.h"
+#include "comm.h"
+#include "force.h"
+#include "citeme.h"
+#include "math_const.h"
+#include "math_extra.h"
+#include "memory.h"
+#include "error.h"
+
+using namespace LAMMPS_NS;
+using namespace MathConst;
+using namespace std;
+using namespace MathExtra;
+
+
+static const char cite_dihedral_tablecut[] =
+ "dihedral_style tablecut command:\n\n"
+ "@Article{Salerno17,\n"
+ " author = {K. M. Salerno and N. Bernstein},\n"
+ " title = {Persistence Length, End-to-End Distance, and Structure of Coarse-Grained Polymers},\n"
+ " journal = {J.~Chem.~Theory Comput.},\n"
+ " year = 2018,\n"
+ " DOI = 10.1021/acs.jctc.7b01229"
+ "}\n\n";
+
+/* ---------------------------------------------------------------------- */
+
+#define TOLERANCE 0.05
+#define SMALL 0.0000001
+
+// ------------------------------------------------------------------------
+// The following auxiliary functions were left out of the
+// DihedralTable class either because they require template parameters,
+// or because they have nothing to do with dihedral angles.
+// ------------------------------------------------------------------------
+
+// -------------------------------------------------------------------
+// --------- The function was stolen verbatim from the ---------
+// --------- GNU Scientific Library (GSL, version 1.15) ---------
+// -------------------------------------------------------------------
+
+/* Author: Gerard Jungman */
+/* for description of method see [Engeln-Mullges + Uhlig, p. 96]
+ *
+ * diag[0] offdiag[0] 0 ..... offdiag[N-1]
+ * offdiag[0] diag[1] offdiag[1] .....
+ * 0 offdiag[1] diag[2]
+ * 0 0 offdiag[2] .....
+ * ... ...
+ * offdiag[N-1] ...
+ *
+ */
+// -- (A non-symmetric version of this function is also available.) --
+
+enum { //GSL status return codes.
+ GSL_FAILURE = -1,
+ GSL_SUCCESS = 0,
+ GSL_ENOMEM = 8,
+ GSL_EZERODIV = 12,
+ GSL_EBADLEN = 19
+};
+
+
+static int solve_cyc_tridiag( const double diag[], size_t d_stride,
+ const double offdiag[], size_t o_stride,
+ const double b[], size_t b_stride,
+ double x[], size_t x_stride,
+ size_t N, bool warn)
+{
+ int status = GSL_SUCCESS;
+ double * delta = (double *) malloc (N * sizeof (double));
+ double * gamma = (double *) malloc (N * sizeof (double));
+ double * alpha = (double *) malloc (N * sizeof (double));
+ double * c = (double *) malloc (N * sizeof (double));
+ double * z = (double *) malloc (N * sizeof (double));
+
+ if (delta == 0 || gamma == 0 || alpha == 0 || c == 0 || z == 0) {
+ if (warn)
+ fprintf(stderr,"Internal Cyclic Spline Error: failed to allocate working space\n");
+
+ if (delta) free(delta);
+ if (gamma) free(gamma);
+ if (alpha) free(alpha);
+ if (c) free(c);
+ if (z) free(z);
+ return GSL_ENOMEM;
+ }
+ else
+ {
+ size_t i, j;
+ double sum = 0.0;
+
+ /* factor */
+
+ if (N == 1)
+ {
+ x[0] = b[0] / diag[0];
+ free(delta);
+ free(gamma);
+ free(alpha);
+ free(c);
+ free(z);
+ return GSL_SUCCESS;
+ }
+
+ alpha[0] = diag[0];
+ gamma[0] = offdiag[0] / alpha[0];
+ delta[0] = offdiag[o_stride * (N-1)] / alpha[0];
+
+ if (alpha[0] == 0) {
+ status = GSL_EZERODIV;
+ }
+
+ for (i = 1; i < N - 2; i++)
+ {
+ alpha[i] = diag[d_stride * i] - offdiag[o_stride * (i-1)] * gamma[i - 1];
+ gamma[i] = offdiag[o_stride * i] / alpha[i];
+ delta[i] = -delta[i - 1] * offdiag[o_stride * (i-1)] / alpha[i];
+ if (alpha[i] == 0) {
+ status = GSL_EZERODIV;
+ }
+ }
+
+ for (i = 0; i < N - 2; i++)
+ {
+ sum += alpha[i] * delta[i] * delta[i];
+ }
+
+ alpha[N - 2] = diag[d_stride * (N - 2)] - offdiag[o_stride * (N - 3)] * gamma[N - 3];
+
+ gamma[N - 2] = (offdiag[o_stride * (N - 2)] - offdiag[o_stride * (N - 3)] * delta[N - 3]) / alpha[N - 2];
+
+ alpha[N - 1] = diag[d_stride * (N - 1)] - sum - alpha[(N - 2)] * gamma[N - 2] * gamma[N - 2];
+
+ /* update */
+ z[0] = b[0];
+ for (i = 1; i < N - 1; i++)
+ {
+ z[i] = b[b_stride * i] - z[i - 1] * gamma[i - 1];
+ }
+ sum = 0.0;
+ for (i = 0; i < N - 2; i++)
+ {
+ sum += delta[i] * z[i];
+ }
+ z[N - 1] = b[b_stride * (N - 1)] - sum - gamma[N - 2] * z[N - 2];
+ for (i = 0; i < N; i++)
+ {
+ c[i] = z[i] / alpha[i];
+ }
+
+ /* backsubstitution */
+ x[x_stride * (N - 1)] = c[N - 1];
+ x[x_stride * (N - 2)] = c[N - 2] - gamma[N - 2] * x[x_stride * (N - 1)];
+ if (N >= 3)
+ {
+ for (i = N - 3, j = 0; j <= N - 3; j++, i--)
+ {
+ x[x_stride * i] = c[i] - gamma[i] * x[x_stride * (i + 1)] - delta[i] * x[x_stride * (N - 1)];
+ }
+ }
+ }
+
+ free (z);
+ free (c);
+ free (alpha);
+ free (gamma);
+ free (delta);
+
+ if ((status == GSL_EZERODIV) && warn)
+ fprintf(stderr, "Internal Cyclic Spline Error: Matrix must be positive definite.\n");
+
+ return status;
+} //solve_cyc_tridiag()
+
+/* ----------------------------------------------------------------------
+ spline and splint routines modified from Numerical Recipes
+------------------------------------------------------------------------- */
+
+static int cyc_spline(double const *xa,
+ double const *ya,
+ int n,
+ double period,
+ double *y2a, bool warn)
+{
+
+ double *diag = new double[n];
+ double *offdiag = new double[n];
+ double *rhs = new double[n];
+ double xa_im1, xa_ip1;
+
+ // In the cyclic case, there are n equations with n unknows.
+ // The for loop sets up the equations we need to solve.
+ // Later we invoke the GSL tridiagonal matrix solver to solve them.
+
+ for(int i=0; i < n; i++) {
+
+ // I have to lookup xa[i+1] and xa[i-1]. This gets tricky because of
+ // periodic boundary conditions. We handle that now.
+ int im1 = i-1;
+ if (im1<0) {
+ im1 += n;
+ xa_im1 = xa[im1] - period;
+ }
+ else
+ xa_im1 = xa[im1];
+
+ int ip1 = i+1;
+ if (ip1>=n) {
+ ip1 -= n;
+ xa_ip1 = xa[ip1] + period;
+ }
+ else
+ xa_ip1 = xa[ip1];
+
+ // Recall that we want to find the y2a[] parameters (there are n of them).
+ // To solve for them, we have a linear equation with n unknowns
+ // (in the cyclic case that is). For details, the non-cyclic case is
+ // explained in equation 3.3.7 in Numerical Recipes in C, p. 115.
+ diag[i] = (xa_ip1 - xa_im1) / 3.0;
+ offdiag[i] = (xa_ip1 - xa[i]) / 6.0;
+ rhs[i] = ((ya[ip1] - ya[i]) / (xa_ip1 - xa[i])) -
+ ((ya[i] - ya[im1]) / (xa[i] - xa_im1));
+ }
+
+ // Because this matrix is tridiagonal (and cyclic), we can use the following
+ // cheap method to invert it.
+ if (solve_cyc_tridiag(diag, 1,
+ offdiag, 1,
+ rhs, 1,
+ y2a, 1,
+ n, warn) != GSL_SUCCESS) {
+ if (warn)
+ fprintf(stderr,"Error in inverting matrix for splines.\n");
+
+ delete [] diag;
+ delete [] offdiag;
+ delete [] rhs;
+ return 1;
+ }
+ delete [] diag;
+ delete [] offdiag;
+ delete [] rhs;
+ return 0;
+} // cyc_spline()
+
+/* ---------------------------------------------------------------------- */
+
+// cyc_splint(): Evaluates a spline at position x, with n control
+// points located at xa[], ya[], and with parameters y2a[]
+// The xa[] must be monotonically increasing and their
+// range should not exceed period (ie xa[n-1] < xa[0] + period).
+// x must lie in the range: [(xa[n-1]-period), (xa[0]+period)]
+// "period" is typically 2*PI.
+static double cyc_splint(double const *xa,
+ double const *ya,
+ double const *y2a,
+ int n,
+ double period,
+ double x)
+{
+ int klo = -1;
+ int khi = n;
+ int k;
+ double xlo = xa[n-1] - period;
+ double xhi = xa[0] + period;
+ while (khi-klo > 1) {
+ k = (khi+klo) >> 1; //(k=(khi+klo)/2)
+ if (xa[k] > x) {
+ khi = k;
+ xhi = xa[k];
+ }
+ else {
+ klo = k;
+ xlo = xa[k];
+ }
+ }
+
+ if (khi == n) khi = 0;
+ if (klo ==-1) klo = n-1;
+
+ double h = xhi-xlo;
+ double a = (xhi-x) / h;
+ double b = (x-xlo) / h;
+ double y = a*ya[klo] + b*ya[khi] +
+ ((a*a*a-a)*y2a[klo] + (b*b*b-b)*y2a[khi]) * (h*h)/6.0;
+
+ return y;
+
+} // cyc_splint()
+
+
+static double cyc_lin(double const *xa,
+ double const *ya,
+ int n,
+ double period,
+ double x)
+{
+ int klo = -1;
+ int khi = n;
+ int k;
+ double xlo = xa[n-1] - period;
+ double xhi = xa[0] + period;
+ while (khi-klo > 1) {
+ k = (khi+klo) >> 1; //(k=(khi+klo)/2)
+ if (xa[k] > x) {
+ khi = k;
+ xhi = xa[k];
+ }
+ else {
+ klo = k;
+ xlo = xa[k];
+ }
+ }
+
+ if (khi == n) khi = 0;
+ if (klo ==-1) klo = n-1;
+
+ double h = xhi-xlo;
+ double a = (xhi-x) / h;
+ double b = (x-xlo) / h;
+ double y = a*ya[klo] + b*ya[khi];
+
+ return y;
+
+} // cyc_lin()
+
+
+
+
+// cyc_splintD(): Evaluate the deriviative of a cyclic spline at position x,
+// with n control points at xa[], ya[], with parameters y2a[].
+// The xa[] must be monotonically increasing and their
+// range should not exceed period (ie xa[n-1] < xa[0] + period).
+// x must lie in the range: [(xa[n-1]-period), (xa[0]+period)]
+// "period" is typically 2*PI.
+static double cyc_splintD(double const *xa,
+ double const *ya,
+ double const *y2a,
+ int n,
+ double period,
+ double x)
+{
+ int klo = -1;
+ int khi = n; // (not n-1)
+ int k;
+ double xlo = xa[n-1] - period;
+ double xhi = xa[0] + period;
+ while (khi-klo > 1) {
+ k = (khi+klo) >> 1; //(k=(khi+klo)/2)
+ if (xa[k] > x) {
+ khi = k;
+ xhi = xa[k];
+ }
+ else {
+ klo = k;
+ xlo = xa[k];
+ }
+ }
+
+ if (khi == n) khi = 0;
+ if (klo ==-1) klo = n-1;
+
+ double yhi = ya[khi];
+ double ylo = ya[klo];
+ double h = xhi-xlo;
+ double g = yhi-ylo;
+ double a = (xhi-x) / h;
+ double b = (x-xlo) / h;
+ // Formula below taken from equation 3.3.5 of "numerical recipes in c"
+ // "yD" = the derivative of y
+ double yD = g/h - ( (3.0*a*a-1.0)*y2a[klo] - (3.0*b*b-1.0)*y2a[khi] ) * h/6.0;
+ // For rerefence: y = a*ylo + b*yhi +
+ // ((a*a*a-a)*y2a[klo] + (b*b*b-b)*y2a[khi]) * (h*h)/6.0;
+
+ return yD;
+
+} // cyc_splintD()
+
+
+/* ---------------------------------------------------------------------- */
+
+DihedralTableCut::DihedralTableCut(LAMMPS *lmp) : Dihedral(lmp)
+{
+ if (lmp->citeme) lmp->citeme->add(cite_dihedral_tablecut);
+ ntables = 0;
+ tables = NULL;
+ checkU_fname = checkF_fname = NULL;
+}
+
+/* ---------------------------------------------------------------------- */
+
+DihedralTableCut::~DihedralTableCut()
+{
+ if (allocated) {
+ memory->destroy(setflag);
+ memory->destroy(setflag_d);
+ memory->destroy(setflag_aat);
+
+ memory->destroy(k1);
+ memory->destroy(k2);
+ memory->destroy(k3);
+ memory->destroy(phi1);
+ memory->destroy(phi2);
+ memory->destroy(phi3);
+
+ memory->destroy(aat_k);
+ memory->destroy(aat_theta0_1);
+ memory->destroy(aat_theta0_2);
+
+ for (int m = 0; m < ntables; m++) free_table(&tables[m]);
+ memory->sfree(tables);
+ memory->sfree(checkU_fname);
+ memory->sfree(checkF_fname);
+
+ memory->destroy(setflag);
+ memory->destroy(tabindex);
+
+ }
+}
+
+/* ---------------------------------------------------------------------- */
+
+void DihedralTableCut::compute(int eflag, int vflag)
+{
+
+ int i1,i2,i3,i4,i,j,k,n,type;
+ double edihedral;
+ double vb1x,vb1y,vb1z,vb2x,vb2y,vb2z,vb3x,vb3y,vb3z,vb2xm,vb2ym,vb2zm;
+ double fphi,fpphi;
+ double r1mag2,r1,r2mag2,r2,r3mag2,r3;
+ double sb1,rb1,sb2,rb2,sb3,rb3,c0,r12c1;
+ double r12c2,costh12,costh13,costh23,sc1,sc2,s1,s2,c;
+ double phi,sinphi,a11,a22,a33,a12,a13,a23,sx1,sx2;
+ double sx12,sy1,sy2,sy12,sz1,sz2,sz12;
+ double t1,t2,t3,t4;
+ double da1,da2;
+ double s12,sin2;
+ double dcosphidr[4][3],dphidr[4][3],dthetadr[2][4][3];
+ double fabcd[4][3];
+
+ edihedral = 0.0;
+ if (eflag || vflag) ev_setup(eflag,vflag);
+ else evflag = 0;
+
+ double **x = atom->x;
+ double **f = atom->f;
+ int **dihedrallist = neighbor->dihedrallist;
+ int ndihedrallist = neighbor->ndihedrallist;
+ int nlocal = atom->nlocal;
+ int newton_bond = force->newton_bond;
+
+ for (n = 0; n < ndihedrallist; n++) {
+ i1 = dihedrallist[n][0];
+ i2 = dihedrallist[n][1];
+ i3 = dihedrallist[n][2];
+ i4 = dihedrallist[n][3];
+ type = dihedrallist[n][4];
+
+ // 1st bond
+
+ vb1x = x[i1][0] - x[i2][0];
+ vb1y = x[i1][1] - x[i2][1];
+ vb1z = x[i1][2] - x[i2][2];
+
+ // 2nd bond
+
+ vb2x = x[i3][0] - x[i2][0];
+ vb2y = x[i3][1] - x[i2][1];
+ vb2z = x[i3][2] - x[i2][2];
+
+ vb2xm = -vb2x;
+ vb2ym = -vb2y;
+ vb2zm = -vb2z;
+
+ // 3rd bond
+
+ vb3x = x[i4][0] - x[i3][0];
+ vb3y = x[i4][1] - x[i3][1];
+ vb3z = x[i4][2] - x[i3][2];
+
+ // distances
+
+ r1mag2 = vb1x*vb1x + vb1y*vb1y + vb1z*vb1z;
+ r1 = sqrt(r1mag2);
+ r2mag2 = vb2x*vb2x + vb2y*vb2y + vb2z*vb2z;
+ r2 = sqrt(r2mag2);
+ r3mag2 = vb3x*vb3x + vb3y*vb3y + vb3z*vb3z;
+ r3 = sqrt(r3mag2);
+
+ sb1 = 1.0/r1mag2;
+ rb1 = 1.0/r1;
+ sb2 = 1.0/r2mag2;
+ rb2 = 1.0/r2;
+ sb3 = 1.0/r3mag2;
+ rb3 = 1.0/r3;
+
+ c0 = (vb1x*vb3x + vb1y*vb3y + vb1z*vb3z) * rb1*rb3;
+
+ // angles
+
+ r12c1 = rb1*rb2;
+ r12c2 = rb2*rb3;
+ costh12 = (vb1x*vb2x + vb1y*vb2y + vb1z*vb2z) * r12c1;
+ costh13 = c0;
+ costh23 = (vb2xm*vb3x + vb2ym*vb3y + vb2zm*vb3z) * r12c2;
+
+ // cos and sin of 2 angles and final c
+
+ sin2 = MAX(1.0 - costh12*costh12,0.0);
+ sc1 = sqrt(sin2);
+ if (sc1 < SMALL) sc1 = SMALL;
+ sc1 = 1.0/sc1;
+
+ sin2 = MAX(1.0 - costh23*costh23,0.0);
+ sc2 = sqrt(sin2);
+ if (sc2 < SMALL) sc2 = SMALL;
+ sc2 = 1.0/sc2;
+
+ s1 = sc1 * sc1;
+ s2 = sc2 * sc2;
+ s12 = sc1 * sc2;
+ c = (c0 + costh12*costh23) * s12;
+
+ // error check
+
+ if (c > 1.0 + TOLERANCE || c < (-1.0 - TOLERANCE)) {
+ int me;
+ MPI_Comm_rank(world,&me);
+ if (screen) {
+ char str[128];
+ sprintf(str,"Dihedral problem: %d " BIGINT_FORMAT " "
+ TAGINT_FORMAT " " TAGINT_FORMAT " "
+ TAGINT_FORMAT " " TAGINT_FORMAT,
+ me,update->ntimestep,
+ atom->tag[i1],atom->tag[i2],atom->tag[i3],atom->tag[i4]);
+ error->warning(FLERR,str,0);
+ fprintf(screen," 1st atom: %d %g %g %g\n",
+ me,x[i1][0],x[i1][1],x[i1][2]);
+ fprintf(screen," 2nd atom: %d %g %g %g\n",
+ me,x[i2][0],x[i2][1],x[i2][2]);
+ fprintf(screen," 3rd atom: %d %g %g %g\n",
+ me,x[i3][0],x[i3][1],x[i3][2]);
+ fprintf(screen," 4th atom: %d %g %g %g\n",
+ me,x[i4][0],x[i4][1],x[i4][2]);
+ }
+ }
+
+ if (c > 1.0) c = 1.0;
+ if (c < -1.0) c = -1.0;
+ double phil = acos(c);
+ phi = acos(c);
+
+ sinphi = sqrt(1.0 - c*c);
+ sinphi = MAX(sinphi,SMALL);
+
+ // n123 = vb1 x vb2
+
+ double n123x = vb1y*vb2z - vb1z*vb2y;
+ double n123y = vb1z*vb2x - vb1x*vb2z;
+ double n123z = vb1x*vb2y - vb1y*vb2x;
+ double n123_dot_vb3 = n123x*vb3x + n123y*vb3y + n123z*vb3z;
+ if (n123_dot_vb3 > 0.0) {
+ phil = -phil;
+ phi = -phi;
+ sinphi = -sinphi;
+ }
+
+ a11 = -c*sb1*s1;
+ a22 = sb2 * (2.0*costh13*s12 - c*(s1+s2));
+ a33 = -c*sb3*s2;
+ a12 = r12c1 * (costh12*c*s1 + costh23*s12);
+ a13 = rb1*rb3*s12;
+ a23 = r12c2 * (-costh23*c*s2 - costh12*s12);
+
+ sx1 = a11*vb1x + a12*vb2x + a13*vb3x;
+ sx2 = a12*vb1x + a22*vb2x + a23*vb3x;
+ sx12 = a13*vb1x + a23*vb2x + a33*vb3x;
+ sy1 = a11*vb1y + a12*vb2y + a13*vb3y;
+ sy2 = a12*vb1y + a22*vb2y + a23*vb3y;
+ sy12 = a13*vb1y + a23*vb2y + a33*vb3y;
+ sz1 = a11*vb1z + a12*vb2z + a13*vb3z;
+ sz2 = a12*vb1z + a22*vb2z + a23*vb3z;
+ sz12 = a13*vb1z + a23*vb2z + a33*vb3z;
+
+ // set up d(cos(phi))/d(r) and dphi/dr arrays
+
+ dcosphidr[0][0] = -sx1;
+ dcosphidr[0][1] = -sy1;
+ dcosphidr[0][2] = -sz1;
+ dcosphidr[1][0] = sx2 + sx1;
+ dcosphidr[1][1] = sy2 + sy1;
+ dcosphidr[1][2] = sz2 + sz1;
+ dcosphidr[2][0] = sx12 - sx2;
+ dcosphidr[2][1] = sy12 - sy2;
+ dcosphidr[2][2] = sz12 - sz2;
+ dcosphidr[3][0] = -sx12;
+ dcosphidr[3][1] = -sy12;
+ dcosphidr[3][2] = -sz12;
+
+ for (i = 0; i < 4; i++)
+ for (j = 0; j < 3; j++)
+ dphidr[i][j] = -dcosphidr[i][j] / sinphi;
+
+
+ for (i = 0; i < 4; i++)
+ for (j = 0; j < 3; j++)
+ fabcd[i][j] = 0;
+ edihedral = 0;
+
+
+ // set up d(theta)/d(r) array
+ // dthetadr(i,j,k) = angle i, atom j, coordinate k
+
+ for (i = 0; i < 2; i++)
+ for (j = 0; j < 4; j++)
+ for (k = 0; k < 3; k++)
+ dthetadr[i][j][k] = 0.0;
+
+ t1 = costh12 / r1mag2;
+ t2 = costh23 / r2mag2;
+ t3 = costh12 / r2mag2;
+ t4 = costh23 / r3mag2;
+
+ // angle12
+
+ dthetadr[0][0][0] = sc1 * ((t1 * vb1x) - (vb2x * r12c1));
+ dthetadr[0][0][1] = sc1 * ((t1 * vb1y) - (vb2y * r12c1));
+ dthetadr[0][0][2] = sc1 * ((t1 * vb1z) - (vb2z * r12c1));
+
+ dthetadr[0][1][0] = sc1 * ((-t1 * vb1x) + (vb2x * r12c1) +
+ (-t3 * vb2x) + (vb1x * r12c1));
+ dthetadr[0][1][1] = sc1 * ((-t1 * vb1y) + (vb2y * r12c1) +
+ (-t3 * vb2y) + (vb1y * r12c1));
+ dthetadr[0][1][2] = sc1 * ((-t1 * vb1z) + (vb2z * r12c1) +
+ (-t3 * vb2z) + (vb1z * r12c1));
+
+ dthetadr[0][2][0] = sc1 * ((t3 * vb2x) - (vb1x * r12c1));
+ dthetadr[0][2][1] = sc1 * ((t3 * vb2y) - (vb1y * r12c1));
+ dthetadr[0][2][2] = sc1 * ((t3 * vb2z) - (vb1z * r12c1));
+
+ // angle23
+
+ dthetadr[1][1][0] = sc2 * ((t2 * vb2x) + (vb3x * r12c2));
+ dthetadr[1][1][1] = sc2 * ((t2 * vb2y) + (vb3y * r12c2));
+ dthetadr[1][1][2] = sc2 * ((t2 * vb2z) + (vb3z * r12c2));
+
+ dthetadr[1][2][0] = sc2 * ((-t2 * vb2x) - (vb3x * r12c2) +
+ (t4 * vb3x) + (vb2x * r12c2));
+ dthetadr[1][2][1] = sc2 * ((-t2 * vb2y) - (vb3y * r12c2) +
+ (t4 * vb3y) + (vb2y * r12c2));
+ dthetadr[1][2][2] = sc2 * ((-t2 * vb2z) - (vb3z * r12c2) +
+ (t4 * vb3z) + (vb2z * r12c2));
+
+ dthetadr[1][3][0] = -sc2 * ((t4 * vb3x) + (vb2x * r12c2));
+ dthetadr[1][3][1] = -sc2 * ((t4 * vb3y) + (vb2y * r12c2));
+ dthetadr[1][3][2] = -sc2 * ((t4 * vb3z) + (vb2z * r12c2));
+
+ // angle/angle/torsion cutoff
+
+ da1 = acos(costh12) - aat_theta0_1[type] ;
+ da2 = acos(costh23) - aat_theta0_1[type] ;
+ double dtheta = aat_theta0_2[type]-aat_theta0_1[type];
+
+ fphi = 0.0;
+ fpphi = 0.0;
+ if (phil < 0) phil +=MY_2PI;
+ uf_lookup(type, phil, fphi, fpphi);
+
+ double gt = aat_k[type];
+ double gtt = aat_k[type];
+ double gpt = 0;
+ double gptt = 0;
+
+ if ( acos(costh12) > aat_theta0_1[type]) {
+ gt *= 1-da1*da1/dtheta/dtheta;
+ gpt = -aat_k[type]*2*da1/dtheta/dtheta;
+ }
+
+ if ( acos(costh23) > aat_theta0_1[type]) {
+ gtt *= 1-da2*da2/dtheta/dtheta;
+ gptt = -aat_k[type]*2*da2/dtheta/dtheta;
+ }
+
+ if (eflag) edihedral = gt*gtt*fphi;
+
+ for (i = 0; i < 4; i++)
+ for (j = 0; j < 3; j++)
+ fabcd[i][j] -= - gt*gtt*fpphi*dphidr[i][j]
+ - gt*gptt*fphi*dthetadr[1][i][j] + gpt*gtt*fphi*dthetadr[0][i][j];
+
+ // apply force to each of 4 atoms
+
+ if (newton_bond || i1 < nlocal) {
+ f[i1][0] += fabcd[0][0];
+ f[i1][1] += fabcd[0][1];
+ f[i1][2] += fabcd[0][2];
+ }
+
+ if (newton_bond || i2 < nlocal) {
+ f[i2][0] += fabcd[1][0];
+ f[i2][1] += fabcd[1][1];
+ f[i2][2] += fabcd[1][2];
+ }
+
+ if (newton_bond || i3 < nlocal) {
+ f[i3][0] += fabcd[2][0];
+ f[i3][1] += fabcd[2][1];
+ f[i3][2] += fabcd[2][2];
+ }
+
+ if (newton_bond || i4 < nlocal) {
+ f[i4][0] += fabcd[3][0];
+ f[i4][1] += fabcd[3][1];
+ f[i4][2] += fabcd[3][2];
+ }
+
+ if (evflag)
+ ev_tally(i1,i2,i3,i4,nlocal,newton_bond,edihedral,
+ fabcd[0],fabcd[2],fabcd[3],
+ vb1x,vb1y,vb1z,vb2x,vb2y,vb2z,vb3x,vb3y,vb3z);
+ }
+}
+
+/* ---------------------------------------------------------------------- */
+
+void DihedralTableCut::allocate()
+{
+ allocated = 1;
+ int n = atom->ndihedraltypes;
+
+ memory->create(k1,n+1,"dihedral:k1");
+ memory->create(k2,n+1,"dihedral:k2");
+ memory->create(k3,n+1,"dihedral:k3");
+ memory->create(phi1,n+1,"dihedral:phi1");
+ memory->create(phi2,n+1,"dihedral:phi2");
+ memory->create(phi3,n+1,"dihedral:phi3");
+
+ memory->create(aat_k,n+1,"dihedral:aat_k");
+ memory->create(aat_theta0_1,n+1,"dihedral:aat_theta0_1");
+ memory->create(aat_theta0_2,n+1,"dihedral:aat_theta0_2");
+
+ memory->create(setflag,n+1,"dihedral:setflag");
+ memory->create(setflag_d,n+1,"dihedral:setflag_d");
+ memory->create(setflag_aat,n+1,"dihedral:setflag_aat");
+
+ memory->create(tabindex,n+1,"dihedral:tabindex");
+ //memory->create(phi0,n+1,"dihedral:phi0"); <-equilibrium angles not supported
+ memory->create(setflag,n+1,"dihedral:setflag");
+
+ for (int i = 1; i <= n; i++)
+ setflag[i] = setflag_d[i] = setflag_aat[i] = 0;
+}
+
+void DihedralTableCut::settings(int narg, char **arg)
+{
+ if (narg != 2) error->all(FLERR,"Illegal dihedral_style command");
+
+ if (strcmp(arg[0],"linear") == 0) tabstyle = LINEAR;
+ else if (strcmp(arg[0],"spline") == 0) tabstyle = SPLINE;
+ else error->all(FLERR,"Unknown table style in dihedral style table_cut");
+
+ tablength = force->inumeric(FLERR,arg[1]);
+ if (tablength < 3)
+ error->all(FLERR,"Illegal number of dihedral table entries");
+ // delete old tables, since cannot just change settings
+
+ for (int m = 0; m < ntables; m++) free_table(&tables[m]);
+ memory->sfree(tables);
+
+ if (allocated) {
+ memory->destroy(setflag);
+ memory->destroy(tabindex);
+ }
+ allocated = 0;
+
+ ntables = 0;
+ tables = NULL;
+}
+
+/* ----------------------------------------------------------------------
+ set coeffs for one or more types
+ arg1 = "aat" -> AngleAngleTorsion coeffs
+ arg1 -> Dihedral coeffs
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::coeff(int narg, char **arg)
+{
+
+ if (narg != 7) error->all(FLERR,"Incorrect args for dihedral coefficients");
+ if (!allocated) allocate();
+ int ilo,ihi;
+ force->bounds(FLERR,arg[0],atom->ndihedraltypes,ilo,ihi);
+
+ int count = 0;
+
+
+ double k_one = force->numeric(FLERR,arg[2]);
+ double theta0_1_one = force->numeric(FLERR,arg[3]);
+ double theta0_2_one = force->numeric(FLERR,arg[4]);
+
+ // convert theta0's from degrees to radians
+
+ for (int i = ilo; i <= ihi; i++) {
+ aat_k[i] = k_one;
+ aat_theta0_1[i] = theta0_1_one/180.0 * MY_PI;
+ aat_theta0_2[i] = theta0_2_one/180.0 * MY_PI;
+ setflag_aat[i] = 1;
+ count++;
+ }
+
+ int me;
+ MPI_Comm_rank(world,&me);
+ tables = (Table *)
+ memory->srealloc(tables,(ntables+1)*sizeof(Table), "dihedral:tables");
+ Table *tb = &tables[ntables];
+ null_table(tb);
+ if (me == 0) read_table(tb,arg[5],arg[6]);
+ bcast_table(tb);
+
+ // --- check the angle data for range errors ---
+ // --- and resolve issues with periodicity ---
+
+ if (tb->ninput < 2) {
+ string err_msg;
+ err_msg = string("Invalid dihedral table length (")
+ + string(arg[5]) + string(").");
+ error->one(FLERR,err_msg.c_str());
+ }
+ else if ((tb->ninput == 2) && (tabstyle == SPLINE)) {
+ string err_msg;
+ err_msg = string("Invalid dihedral spline table length. (Try linear)\n (")
+ + string(arg[5]) + string(").");
+ error->one(FLERR,err_msg.c_str());
+ }
+
+ // check for monotonicity
+ for (int i=0; i < tb->ninput-1; i++) {
+ if (tb->phifile[i] >= tb->phifile[i+1]) {
+ stringstream i_str;
+ i_str << i+1;
+ string err_msg =
+ string("Dihedral table values are not increasing (") +
+ string(arg[5]) + string(", ")+i_str.str()+string("th entry)");
+ if (i==0)
+ err_msg += string("\n(This is probably a mistake with your table format.)\n");
+ error->all(FLERR,err_msg.c_str());
+ }
+ }
+
+ // check the range of angles
+ double philo = tb->phifile[0];
+ double phihi = tb->phifile[tb->ninput-1];
+ if (tb->use_degrees) {
+ if ((phihi - philo) >= 360) {
+ string err_msg;
+ err_msg = string("Dihedral table angle range must be < 360 degrees (")
+ +string(arg[5]) + string(").");
+ error->all(FLERR,err_msg.c_str());
+ }
+ }
+ else {
+ if ((phihi - philo) >= MY_2PI) {
+ string err_msg;
+ err_msg = string("Dihedral table angle range must be < 2*PI radians (")
+ + string(arg[5]) + string(").");
+ error->all(FLERR,err_msg.c_str());
+ }
+ }
+
+ // convert phi from degrees to radians
+ if (tb->use_degrees) {
+ for (int i=0; i < tb->ninput; i++) {
+ tb->phifile[i] *= MY_PI/180.0;
+ // I assume that if angles are in degrees, then the forces (f=dU/dphi)
+ // are specified with "phi" in degrees as well.
+ tb->ffile[i] *= 180.0/MY_PI;
+ }
+ }
+
+ // We want all the phi dihedral angles to lie in the range from 0 to 2*PI.
+ // But I don't want to restrict users to input their data in this range.
+ // We also want the angles to be sorted in increasing order.
+ // This messy code fixes these problems with the user's data:
+ {
+ double *phifile_tmp = new double [tb->ninput]; //temporary arrays
+ double *ffile_tmp = new double [tb->ninput]; //used for sorting
+ double *efile_tmp = new double [tb->ninput];
+
+ // After re-imaging, does the range of angles cross the 0 or 2*PI boundary?
+ // If so, find the discontinuity:
+ int i_discontinuity = tb->ninput;
+ for (int i=0; i < tb->ninput; i++) {
+ double phi = tb->phifile[i];
+ // Add a multiple of 2*PI to phi until it lies in the range [0, 2*PI).
+ phi -= MY_2PI * floor(phi/MY_2PI);
+ phifile_tmp[i] = phi;
+ efile_tmp[i] = tb->efile[i];
+ ffile_tmp[i] = tb->ffile[i];
+ if ((i>0) && (phifile_tmp[i] < phifile_tmp[i-1])) {
+ //There should only be at most one discontinuity, because we have
+ //insured that the data was sorted before imaging, and because the
+ //range of angle values does not exceed 2*PI.
+ i_discontinuity = i;
+ }
+ }
+
+ int I = 0;
+ for (int i = i_discontinuity; i < tb->ninput; i++) {
+ tb->phifile[I] = phifile_tmp[i];
+ tb->efile[I] = efile_tmp[i];
+ tb->ffile[I] = ffile_tmp[i];
+ I++;
+ }
+ for (int i = 0; i < i_discontinuity; i++) {
+ tb->phifile[I] = phifile_tmp[i];
+ tb->efile[I] = efile_tmp[i];
+ tb->ffile[I] = ffile_tmp[i];
+ I++;
+ }
+
+ // clean up temporary storage
+ delete[] phifile_tmp;
+ delete[] ffile_tmp;
+ delete[] efile_tmp;
+ }
+
+ // spline read-in and compute r,e,f vectors within table
+
+ spline_table(tb);
+ compute_table(tb);
+
+ // Optional: allow the user to print out the interpolated spline tables
+
+ if (me == 0) {
+ if (checkU_fname && (strlen(checkU_fname) != 0))
+ {
+ ofstream checkU_file;
+ checkU_file.open(checkU_fname, ios::out);
+ for (int i=0; i < tablength; i++) {
+ double phi = i*MY_2PI/tablength;
+ double u = tb->e[i];
+ if (tb->use_degrees)
+ phi *= 180.0/MY_PI;
+ checkU_file << phi << " " << u << "\n";
+ }
+ checkU_file.close();
+ }
+ if (checkF_fname && (strlen(checkF_fname) != 0))
+ {
+ ofstream checkF_file;
+ checkF_file.open(checkF_fname, ios::out);
+ for (int i=0; i < tablength; i++)
+ {
+ double phi = i*MY_2PI/tablength;
+ double f;
+ if ((tabstyle == SPLINE) && (tb->f_unspecified)) {
+ double dU_dphi =
+ // (If the user did not specify the forces now, AND the user
+ // selected the "spline" option, (as opposed to "linear")
+ // THEN the tb->f array is uninitialized, so there's
+ // no point to print out the contents of the tb->f[] array.
+ // Instead, later on, we will calculate the force using the
+ // -cyc_splintD() routine to calculate the derivative of the
+ // energy spline, using the energy data (tb->e[]).
+ // To be nice and report something, I do the same thing here.)
+ cyc_splintD(tb->phi, tb->e, tb->e2, tablength, MY_2PI,phi);
+ f = -dU_dphi;
+ }
+ else
+ // Otherwise we calculated the tb->f[] array. Report its contents.
+ f = tb->f[i];
+ if (tb->use_degrees) {
+ phi *= 180.0/MY_PI;
+ // If the user wants degree angle units, we should convert our
+ // internal force tables (in energy/radians) to (energy/degrees)
+ f *= MY_PI/180.0;
+ }
+ checkF_file << phi << " " << f << "\n";
+ }
+ checkF_file.close();
+ } // if (checkF_fname && (strlen(checkF_fname) != 0))
+ } // if (me == 0)
+
+ // store ptr to table in tabindex
+ count = 0;
+ for (int i = ilo; i <= ihi; i++)
+ {
+ tabindex[i] = ntables;
+ //phi0[i] = tb->phi0; <- equilibrium dihedral angles not supported
+ setflag[i] = 1;
+ count++;
+ }
+ ntables++;
+
+
+ if (count == 0) error->all(FLERR,"Incorrect args for dihedral coefficients");
+
+ for (int i = ilo; i <= ihi; i++)
+ if (setflag_d[i] == 1 && setflag_aat[i] == 1 )
+ setflag[i] = 1;
+}
+
+/* ----------------------------------------------------------------------
+ proc 0 writes out coeffs to restart file
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::write_restart(FILE *fp)
+{
+ fwrite(&tabstyle,sizeof(int),1,fp);
+ fwrite(&tablength,sizeof(int),1,fp);
+ fwrite(&k1[1],sizeof(double),atom->ndihedraltypes,fp);
+ fwrite(&k2[1],sizeof(double),atom->ndihedraltypes,fp);
+ fwrite(&k3[1],sizeof(double),atom->ndihedraltypes,fp);
+ fwrite(&phi1[1],sizeof(double),atom->ndihedraltypes,fp);
+ fwrite(&phi2[1],sizeof(double),atom->ndihedraltypes,fp);
+ fwrite(&phi3[1],sizeof(double),atom->ndihedraltypes,fp);
+
+ fwrite(&aat_k[1],sizeof(double),atom->ndihedraltypes,fp);
+ fwrite(&aat_theta0_1[1],sizeof(double),atom->ndihedraltypes,fp);
+ fwrite(&aat_theta0_2[1],sizeof(double),atom->ndihedraltypes,fp);
+}
+
+/* ----------------------------------------------------------------------
+ proc 0 reads coeffs from restart file, bcasts them
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::read_restart(FILE *fp)
+{
+ allocate();
+
+ if (comm->me == 0) {
+ fread(&tabstyle,sizeof(int),1,fp);
+ fread(&tablength,sizeof(int),1,fp);
+ fread(&k1[1],sizeof(double),atom->ndihedraltypes,fp);
+ fread(&k2[1],sizeof(double),atom->ndihedraltypes,fp);
+ fread(&k3[1],sizeof(double),atom->ndihedraltypes,fp);
+ fread(&phi1[1],sizeof(double),atom->ndihedraltypes,fp);
+ fread(&phi2[1],sizeof(double),atom->ndihedraltypes,fp);
+ fread(&phi3[1],sizeof(double),atom->ndihedraltypes,fp);
+
+ fread(&aat_k[1],sizeof(double),atom->ndihedraltypes,fp);
+ fread(&aat_theta0_1[1],sizeof(double),atom->ndihedraltypes,fp);
+ fread(&aat_theta0_2[1],sizeof(double),atom->ndihedraltypes,fp);
+ }
+
+ MPI_Bcast(&k1[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&k2[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&k3[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&phi1[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&phi2[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&phi3[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+
+ MPI_Bcast(&aat_k[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&aat_theta0_1[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&aat_theta0_2[1],atom->ndihedraltypes,MPI_DOUBLE,0,world);
+ MPI_Bcast(&tabstyle,1,MPI_INT,0,world);
+ MPI_Bcast(&tablength,1,MPI_INT,0,world);
+
+ allocate();
+ for (int i = 1; i <= atom->ndihedraltypes; i++) setflag[i] = 1;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void DihedralTableCut::null_table(Table *tb)
+{
+ tb->phifile = tb->efile = tb->ffile = NULL;
+ tb->e2file = tb->f2file = NULL;
+ tb->phi = tb->e = tb->de = NULL;
+ tb->f = tb->df = tb->e2 = tb->f2 = NULL;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void DihedralTableCut::free_table(Table *tb)
+{
+ memory->destroy(tb->phifile);
+ memory->destroy(tb->efile);
+ memory->destroy(tb->ffile);
+ memory->destroy(tb->e2file);
+ memory->destroy(tb->f2file);
+
+ memory->destroy(tb->phi);
+ memory->destroy(tb->e);
+ memory->destroy(tb->de);
+ memory->destroy(tb->f);
+ memory->destroy(tb->df);
+ memory->destroy(tb->e2);
+ memory->destroy(tb->f2);
+}
+
+/* ----------------------------------------------------------------------
+ read table file, only called by proc 0
+------------------------------------------------------------------------- */
+
+static const int MAXLINE=2048;
+
+void DihedralTableCut::read_table(Table *tb, char *file, char *keyword)
+{
+ char line[MAXLINE];
+
+ // open file
+
+ FILE *fp = force->open_potential(file);
+ if (fp == NULL) {
+ string err_msg = string("Cannot open file ") + string(file);
+ error->one(FLERR,err_msg.c_str());
+ }
+
+ // loop until section found with matching keyword
+
+ while (1) {
+ if (fgets(line,MAXLINE,fp) == NULL) {
+ string err_msg=string("Did not find keyword \"")
+ +string(keyword)+string("\" in dihedral table file.");
+ error->one(FLERR, err_msg.c_str());
+ }
+ if (strspn(line," \t\n\r") == strlen(line)) continue; // blank line
+ if (line[0] == '#') continue; // comment
+ char *word = strtok(line," \t\n\r");
+ if (strcmp(word,keyword) == 0) break; // matching keyword
+ fgets(line,MAXLINE,fp); // no match, skip section
+ param_extract(tb,line);
+ fgets(line,MAXLINE,fp);
+ for (int i = 0; i < tb->ninput; i++)
+ fgets(line,MAXLINE,fp);
+ }
+
+ // read args on 2nd line of section
+ // allocate table arrays for file values
+
+ fgets(line,MAXLINE,fp);
+ param_extract(tb,line);
+ memory->create(tb->phifile,tb->ninput,"dihedral:phifile");
+ memory->create(tb->efile,tb->ninput,"dihedral:efile");
+ memory->create(tb->ffile,tb->ninput,"dihedral:ffile");
+
+ // read a,e,f table values from file
+
+ int itmp;
+ for (int i = 0; i < tb->ninput; i++) {
+ // Read the next line. Make sure the file is long enough.
+ if (! fgets(line,MAXLINE,fp))
+ error->one(FLERR, "Dihedral table does not contain enough entries.");
+ // Skip blank lines and delete text following a '#' character
+ char *pe = strchr(line, '#');
+ if (pe != NULL) *pe = '\0'; //terminate string at '#' character
+ char *pc = line;
+ while ((*pc != '\0') && isspace(*pc))
+ pc++;
+ if (*pc != '\0') { //If line is not a blank line
+ stringstream line_ss(line);
+ if (tb->f_unspecified) {
+ //sscanf(line,"%d %lg %lg",
+ // &itmp,&tb->phifile[i],&tb->efile[i]);
+ line_ss >> itmp;
+ line_ss >> tb->phifile[i];
+ line_ss >> tb->efile[i];
+ }
+ else {
+ //sscanf(line,"%d %lg %lg %lg",
+ // &itmp,&tb->phifile[i],&tb->efile[i],&tb->ffile[i]);
+ line_ss >> itmp;
+ line_ss >> tb->phifile[i];
+ line_ss >> tb->efile[i];
+ line_ss >> tb->ffile[i];
+ }
+ if (! line_ss) {
+ stringstream err_msg;
+ err_msg << "Read error in table "<< keyword<<", near line "<<i+1<<"\n"
+ << " (Check to make sure the number of colums is correct.)";
+ if ((! tb->f_unspecified) && (i==0))
+ err_msg << "\n (This sometimes occurs if users forget to specify the \"NOF\" option.)\n";
+ error->one(FLERR, err_msg.str().c_str());
+ }
+ }
+ else //if it is a blank line, then skip it.
+ i--;
+ } //for (int i = 0; (i < tb->ninput) && fp; i++) {
+
+ fclose(fp);
+}
+
+/* ----------------------------------------------------------------------
+ build spline representation of e,f over entire range of read-in table
+ this function sets these values in e2file,f2file.
+ I also perform a crude check for force & energy consistency.
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::spline_table(Table *tb)
+{
+ memory->create(tb->e2file,tb->ninput,"dihedral:e2file");
+ memory->create(tb->f2file,tb->ninput,"dihedral:f2file");
+
+ if (cyc_spline(tb->phifile, tb->efile, tb->ninput,
+ MY_2PI,tb->e2file,comm->me == 0))
+ error->one(FLERR,"Error computing dihedral spline tables");
+
+ if (! tb->f_unspecified) {
+ if (cyc_spline(tb->phifile, tb->ffile, tb->ninput,
+ MY_2PI, tb->f2file, comm->me == 0))
+ error->one(FLERR,"Error computing dihedral spline tables");
+ }
+
+ // CHECK to help make sure the user calculated forces in a way
+ // which is grossly numerically consistent with the energy table.
+ if (! tb->f_unspecified) {
+ int num_disagreements = 0;
+ for (int i=0; i<tb->ninput; i++) {
+
+ // Calculate what the force should be at the control points
+ // by using linear interpolation of the derivatives of the energy:
+
+ double phi_i = tb->phifile[i];
+
+ // First deal with periodicity
+ double phi_im1, phi_ip1;
+ int im1 = i-1;
+ if (im1 < 0) {
+ im1 += tb->ninput;
+ phi_im1 = tb->phifile[im1] - MY_2PI;
+ }
+ else
+ phi_im1 = tb->phifile[im1];
+ int ip1 = i+1;
+ if (ip1 >= tb->ninput) {
+ ip1 -= tb->ninput;
+ phi_ip1 = tb->phifile[ip1] + MY_2PI;
+ }
+ else
+ phi_ip1 = tb->phifile[ip1];
+
+ // Now calculate the midpoints above and below phi_i = tb->phifile[i]
+ double phi_lo= 0.5*(phi_im1 + phi_i); //midpoint between phi_im1 and phi_i
+ double phi_hi= 0.5*(phi_i + phi_ip1); //midpoint between phi_i and phi_ip1
+
+ // Use a linear approximation to the derivative at these two midpoints
+ double dU_dphi_lo =
+ (tb->efile[i] - tb->efile[im1])
+ /
+ (phi_i - phi_im1);
+ double dU_dphi_hi =
+ (tb->efile[ip1] - tb->efile[i])
+ /
+ (phi_ip1 - phi_i);
+
+ // Now calculate the derivative at position
+ // phi_i (=tb->phifile[i]) using linear interpolation
+
+ double a = (phi_i - phi_lo) / (phi_hi - phi_lo);
+ double b = (phi_hi - phi_i) / (phi_hi - phi_lo);
+ double dU_dphi = a*dU_dphi_lo + b*dU_dphi_hi;
+ double f = -dU_dphi;
+ // Alternately, we could use spline interpolation instead:
+ // double f = - splintD(tb->phifile, tb->efile, tb->e2file,
+ // tb->ninput, MY_2PI, tb->phifile[i]);
+ // This is the way I originally did it, but I trust
+ // the ugly simple linear way above better.
+ // Recall this entire block of code doess not calculate
+ // anything important. It does not have to be perfect.
+ // We are only checking for stupid user errors here.
+
+ if ((f != 0.0) &&
+ (tb->ffile[i] != 0.0) &&
+ ((f/tb->ffile[i] < 0.5) || (f/tb->ffile[i] > 2.0))) {
+ num_disagreements++;
+ }
+ } // for (int i=0; i<tb->ninput; i++)
+
+ if ((num_disagreements > tb->ninput/2) && (num_disagreements > 2)) {
+ string msg("Dihedral table has inconsistent forces and energies. (Try \"NOF\".)\n");
+ error->all(FLERR,msg.c_str());
+ }
+
+ } // check for consistency if (! tb->f_unspecified)
+
+} // DihedralTable::spline_table()
+
+
+/* ----------------------------------------------------------------------
+ compute a,e,f vectors from splined values
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::compute_table(Table *tb)
+{
+ //delta = table spacing in dihedral angle for tablength (cyclic) bins
+ tb->delta = MY_2PI / tablength;
+ tb->invdelta = 1.0/tb->delta;
+ tb->deltasq6 = tb->delta*tb->delta / 6.0;
+
+ // N evenly spaced bins in dihedral angle from 0 to 2*PI
+ // phi,e,f = value at lower edge of bin
+ // de,df values = delta values of e,f (cyclic, in this case)
+ // phi,e,f,de,df are arrays containing "tablength" number of entries
+
+ memory->create(tb->phi,tablength,"dihedral:phi");
+ memory->create(tb->e,tablength,"dihedral:e");
+ memory->create(tb->de,tablength,"dihedral:de");
+ memory->create(tb->f,tablength,"dihedral:f");
+ memory->create(tb->df,tablength,"dihedral:df");
+ memory->create(tb->e2,tablength,"dihedral:e2");
+ memory->create(tb->f2,tablength,"dihedral:f2");
+
+ if (tabstyle == SPLINE) {
+ // Use cubic spline interpolation to calculate the entries in the
+ // internal table. (This is true regardless...even if tabstyle!=SPLINE.)
+ for (int i = 0; i < tablength; i++) {
+ double phi = i*tb->delta;
+ tb->phi[i] = phi;
+ tb->e[i]= cyc_splint(tb->phifile,tb->efile,tb->e2file,tb->ninput,MY_2PI,phi);
+ if (! tb->f_unspecified)
+ tb->f[i] = cyc_splint(tb->phifile,tb->ffile,tb->f2file,tb->ninput,MY_2PI,phi);
+ }
+ } // if (tabstyle == SPLINE)
+ else if (tabstyle == LINEAR) {
+ if (! tb->f_unspecified) {
+ for (int i = 0; i < tablength; i++) {
+ double phi = i*tb->delta;
+ tb->phi[i] = phi;
+ tb->e[i]= cyc_lin(tb->phifile,tb->efile,tb->ninput,MY_2PI,phi);
+ tb->f[i]= cyc_lin(tb->phifile,tb->ffile,tb->ninput,MY_2PI,phi);
+ }
+ }
+ else {
+ for (int i = 0; i < tablength; i++) {
+ double phi = i*tb->delta;
+ tb->phi[i] = phi;
+ tb->e[i]= cyc_lin(tb->phifile,tb->efile,tb->ninput,MY_2PI,phi);
+ }
+ // In the linear case, if the user did not specify the forces, then we
+ // must generate the "f" array. Do this using linear interpolation
+ // of the e array (which itself was generated above)
+ for (int i = 0; i < tablength; i++) {
+ int im1 = i-1; if (im1 < 0) im1 += tablength;
+ int ip1 = i+1; if (ip1 >= tablength) ip1 -= tablength;
+ double dedx = (tb->e[ip1] - tb->e[im1]) / (2.0 * tb->delta);
+ // (This is the average of the linear slopes on either side of the node.
+ // Note that the nodes in the internal table are evenly spaced.)
+ tb->f[i] = -dedx;
+ }
+ }
+
+
+ // Fill the linear interpolation tables (de, df)
+ for (int i = 0; i < tablength; i++) {
+ int ip1 = i+1; if (ip1 >= tablength) ip1 -= tablength;
+ tb->de[i] = tb->e[ip1] - tb->e[i];
+ tb->df[i] = tb->f[ip1] - tb->f[i];
+ }
+ } // else if (tabstyle == LINEAR)
+
+
+
+ cyc_spline(tb->phi, tb->e, tablength, MY_2PI, tb->e2, comm->me == 0);
+ if (! tb->f_unspecified)
+ cyc_spline(tb->phi, tb->f, tablength, MY_2PI, tb->f2, comm->me == 0);
+}
+
+
+/* ----------------------------------------------------------------------
+ extract attributes from parameter line in table section
+ format of line: N value NOF DEGREES RADIANS
+ N is required, other params are optional
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::param_extract(Table *tb, char *line)
+{
+ //tb->theta0 = 180.0; <- equilibrium angles not supported
+ tb->ninput = 0;
+ tb->f_unspecified = false; //default
+ tb->use_degrees = true; //default
+
+ char *word = strtok(line," \t\n\r\f");
+ while (word) {
+ if (strcmp(word,"N") == 0) {
+ word = strtok(NULL," \t\n\r\f");
+ tb->ninput = atoi(word);
+ }
+ else if (strcmp(word,"NOF") == 0) {
+ tb->f_unspecified = true;
+ }
+ else if ((strcmp(word,"DEGREES") == 0) || (strcmp(word,"degrees") == 0)) {
+ tb->use_degrees = true;
+ }
+ else if ((strcmp(word,"RADIANS") == 0) || (strcmp(word,"radians") == 0)) {
+ tb->use_degrees = false;
+ }
+ else if (strcmp(word,"CHECKU") == 0) {
+ word = strtok(NULL," \t\n\r\f");
+ memory->sfree(checkU_fname);
+ memory->create(checkU_fname,strlen(word)+1,"dihedral_table:checkU");
+ strcpy(checkU_fname, word);
+ }
+ else if (strcmp(word,"CHECKF") == 0) {
+ word = strtok(NULL," \t\n\r\f");
+ memory->sfree(checkF_fname);
+ memory->create(checkF_fname,strlen(word)+1,"dihedral_table:checkF");
+ strcpy(checkF_fname, word);
+ }
+ // COMMENTING OUT: equilibrium angles are not supported
+ //else if (strcmp(word,"EQ") == 0) {
+ // word = strtok(NULL," \t\n\r\f");
+ // tb->theta0 = atof(word);
+ //}
+ else {
+ string err_msg("Invalid keyword in dihedral angle table parameters");
+ err_msg += string(" (") + string(word) + string(")");
+ error->one(FLERR,err_msg.c_str());
+ }
+ word = strtok(NULL," \t\n\r\f");
+ }
+
+ if (tb->ninput == 0)
+ error->one(FLERR,"Dihedral table parameters did not set N");
+
+} // DihedralTable::param_extract()
+
+
+/* ----------------------------------------------------------------------
+ broadcast read-in table info from proc 0 to other procs
+ this function communicates these values in Table:
+ ninput,phifile,efile,ffile,
+ f_unspecified,use_degrees
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::bcast_table(Table *tb)
+{
+ MPI_Bcast(&tb->ninput,1,MPI_INT,0,world);
+
+ int me;
+ MPI_Comm_rank(world,&me);
+ if (me > 0) {
+ memory->create(tb->phifile,tb->ninput,"dihedral:phifile");
+ memory->create(tb->efile,tb->ninput,"dihedral:efile");
+ memory->create(tb->ffile,tb->ninput,"dihedral:ffile");
+ }
+
+ MPI_Bcast(tb->phifile,tb->ninput,MPI_DOUBLE,0,world);
+ MPI_Bcast(tb->efile,tb->ninput,MPI_DOUBLE,0,world);
+ MPI_Bcast(tb->ffile,tb->ninput,MPI_DOUBLE,0,world);
+
+ MPI_Bcast(&tb->f_unspecified,1,MPI_INT,0,world);
+ MPI_Bcast(&tb->use_degrees,1,MPI_INT,0,world);
+
+ // COMMENTING OUT: equilibrium angles are not supported
+ //MPI_Bcast(&tb->theta0,1,MPI_DOUBLE,0,world);
+}
+
+
+
+/* ----------------------------------------------------------------------
+ proc 0 writes to data file
+------------------------------------------------------------------------- */
+
+void DihedralTableCut::write_data(FILE *fp)
+{
+ for (int i = 1; i <= atom->ndihedraltypes; i++)
+ fprintf(fp,"%d %g %g %g %g %g %g\n",i,
+ k1[i],phi1[i]*180.0/MY_PI,
+ k2[i],phi2[i]*180.0/MY_PI,
+ k3[i],phi3[i]*180.0/MY_PI);
+
+ fprintf(fp,"\nAngleAngleTorsion Coeffs\n\n");
+ for (int i = 1; i <= atom->ndihedraltypes; i++)
+ fprintf(fp,"%d %g %g %g\n",i,aat_k[i],
+ aat_theta0_1[i]*180.0/MY_PI,aat_theta0_2[i]*180.0/MY_PI);
+
+}
diff --git a/src/USER-MISC/dihedral_table_cut.h b/src/USER-MISC/dihedral_table_cut.h
new file mode 100644
index 0000000000000000000000000000000000000000..d01903c88b75d722aa2b5861a6cba4e00817ed87
--- /dev/null
+++ b/src/USER-MISC/dihedral_table_cut.h
@@ -0,0 +1,173 @@
+/* -*- 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.
+------------------------------------------------------------------------- */
+
+#ifdef DIHEDRAL_CLASS
+
+DihedralStyle(table/cut,DihedralTableCut)
+
+#else
+
+#ifndef LMP_DIHEDRAL_TABLE_CUT_H
+#define LMP_DIHEDRAL_TABLE_CUT_H
+
+#include <stdio.h>
+#include "dihedral.h"
+
+namespace LAMMPS_NS {
+
+class DihedralTableCut : public Dihedral {
+ public:
+ DihedralTableCut(class LAMMPS *);
+ virtual ~DihedralTableCut();
+ virtual void compute(int, int);
+ void settings(int, char **);
+ void coeff(int, char **);
+ void write_restart(FILE *);
+ void read_restart(FILE *);
+ void write_data(FILE *);
+ double single(int type, int i1, int i2, int i3, int i4);
+
+ protected:
+ double *k1,*k2,*k3;
+ double *phi1,*phi2,*phi3;
+ double *aat_k,*aat_theta0_1,*aat_theta0_2;
+ int *setflag_d;
+ int *setflag_aat;
+
+ void allocate();
+
+ int tabstyle,tablength;
+ // double *phi0; <- equilibrium angles not supported
+ char *checkU_fname;
+ char *checkF_fname;
+
+ struct Table {
+ int ninput;
+ //double phi0; <-equilibrium angles not supported
+ int f_unspecified; // boolean (but MPI does not like type "bool")
+ int use_degrees; // boolean (but MPI does not like type "bool")
+ double *phifile,*efile,*ffile;
+ double *e2file,*f2file;
+ double delta,invdelta,deltasq6;
+ double *phi,*e,*de,*f,*df,*e2,*f2;
+ };
+
+ int ntables;
+ Table *tables;
+ int *tabindex;
+
+ void null_table(Table *);
+ void free_table(Table *);
+ void read_table(Table *, char *, char *);
+ void bcast_table(Table *);
+ void spline_table(Table *);
+ void compute_table(Table *);
+
+ void param_extract(Table *, char *);
+
+ // --------------------------------------------
+ // ------------ inline functions --------------
+ // --------------------------------------------
+
+ // -----------------------------------------------------------
+ // uf_lookup()
+ // quickly calculate the potential u and force f at angle x,
+ // using the internal tables tb->e and tb->f (evenly spaced)
+ // -----------------------------------------------------------
+ enum{LINEAR,SPLINE};
+
+ inline void uf_lookup(int type, double x, double &u, double &f)
+ {
+ Table *tb = &tables[tabindex[type]];
+ double x_over_delta = x*tb->invdelta;
+ int i = static_cast<int> (x_over_delta);
+ double a;
+ double b = x_over_delta - i;
+ // Apply periodic boundary conditions to indices i and i+1
+ if (i >= tablength) i -= tablength;
+ int ip1 = i+1; if (ip1 >= tablength) ip1 -= tablength;
+
+ switch(tabstyle) {
+ case LINEAR:
+ u = tb->e[i] + b * tb->de[i];
+ f = -(tb->f[i] + b * tb->df[i]); //<--works even if tb->f_unspecified==true
+ break;
+ case SPLINE:
+ a = 1.0 - b;
+ u = a * tb->e[i] + b * tb->e[ip1] +
+ ((a*a*a-a)*tb->e2[i] + (b*b*b-b)*tb->e2[ip1]) *
+ tb->deltasq6;
+ if (tb->f_unspecified)
+ //Formula below taken from equation3.3.5 of "numerical recipes in c"
+ //"f"=-derivative of e with respect to x (or "phi" in this case)
+ f = -((tb->e[i]-tb->e[ip1])*tb->invdelta +
+ ((3.0*a*a-1.0)*tb->e2[i]+(1.0-3.0*b*b)*tb->e2[ip1])*tb->delta/6.0);
+ else
+ f = -(a * tb->f[i] + b * tb->f[ip1] +
+ ((a*a*a-a)*tb->f2[i] + (b*b*b-b)*tb->f2[ip1]) *
+ tb->deltasq6);
+ break;
+ } // switch(tabstyle)
+ } // uf_lookup()
+
+
+ // ----------------------------------------------------------
+ // u_lookup()
+ // quickly calculate the potential u at angle x using tb->e
+ //-----------------------------------------------------------
+
+ inline void u_lookup(int type, double x, double &u)
+ {
+ Table *tb = &tables[tabindex[type]];
+ int N = tablength;
+
+ // i = static_cast<int> ((x - tb->lo) * tb->invdelta); <-general version
+ double x_over_delta = x*tb->invdelta;
+ int i = static_cast<int> (x_over_delta);
+ double b = x_over_delta - i;
+
+ // Apply periodic boundary conditions to indices i and i+1
+ if (i >= N) i -= N;
+ int ip1 = i+1; if (ip1 >= N) ip1 -= N;
+
+ if (tabstyle == LINEAR) {
+ u = tb->e[i] + b * tb->de[i];
+ }
+ else if (tabstyle == SPLINE) {
+ double a = 1.0 - b;
+ u = a * tb->e[i] + b * tb->e[ip1] +
+ ((a*a*a-a)*tb->e2[i] + (b*b*b-b)*tb->e2[ip1]) *
+ tb->deltasq6;
+ }
+ } // u_lookup()
+
+
+};
+
+}
+
+#endif
+#endif
+
+/* ERROR/WARNING messages:
+
+W: Dihedral problem: %d %ld %d %d %d %d
+
+Conformation of the 4 listed dihedral atoms is extreme; you may want
+to check your simulation geometry.
+
+E: Incorrect args for dihedral coefficients
+
+Self-explanatory. Check the input script or data file.
+
+*/
diff --git a/src/USER-MISC/fix_bond_react.cpp b/src/USER-MISC/fix_bond_react.cpp
index 86da8a0c4fb306297f2cbdac9575e3cd32ec1dc9..4120ce0df5f8eb715008478af71d77320d440c7c 100644
--- a/src/USER-MISC/fix_bond_react.cpp
+++ b/src/USER-MISC/fix_bond_react.cpp
@@ -410,14 +410,14 @@ it will have the name 'i_limit_tags' and will be intitialized to 0 (not in group
void FixBondReact::post_constructor()
{
// let's add the limit_tags per-atom property fix
- int len = strlen("per_atom_props") + 1;
+ int len = strlen("bond_react_props_internal") + 1;
id_fix2 = new char[len];
- strcpy(id_fix2,"per_atom_props");
+ strcpy(id_fix2,"bond_react_props_internal");
int ifix = modify->find_fix(id_fix2);
if (ifix == -1) {
char **newarg = new char*[8];
- newarg[0] = (char *) "per_atom_props";
+ newarg[0] = (char *) "bond_react_props_internal";
newarg[1] = (char *) "all"; // group ID is ignored
newarg[2] = (char *) "property/atom";
newarg[3] = (char *) "i_limit_tags";
@@ -1155,6 +1155,8 @@ void FixBondReact::make_a_guess()
for (int i = 0; i < atom->ntypes; i++) {
if (mol_ntypes[i] != lcl_ntypes[i]) {
status = GUESSFAIL;
+ delete [] mol_ntypes;
+ delete [] lcl_ntypes;
return;
}
}
@@ -2022,6 +2024,7 @@ void FixBondReact::update_everything()
int nn = equivalences[n][1][rxnID]-1;
if (n!=j && bond_atom[atom->map(update_mega_glove[jj+1][i])][p] == update_mega_glove[nn+1][i] && landlocked_atoms[n][rxnID] == 1) {
for (int m = p; m < num_bond[atom->map(update_mega_glove[jj+1][i])]-1; m++) {
+ bond_type[atom->map(update_mega_glove[jj+1][i])][m] = bond_type[atom->map(update_mega_glove[jj+1][i])][m+1];
bond_atom[atom->map(update_mega_glove[jj+1][i])][m] = bond_atom[atom->map(update_mega_glove[jj+1][i])][m+1];
}
num_bond[atom->map(update_mega_glove[jj+1][i])]--;
@@ -2059,6 +2062,7 @@ void FixBondReact::update_everything()
}
}
+ // Angles! First let's delete all angle info:
if (force->angle && twomol->angleflag) {
int *num_angle = atom->num_angle;
int **angle_type = atom->angle_type;
@@ -2069,7 +2073,6 @@ void FixBondReact::update_everything()
for (int i = 0; i < update_num_mega; i++) {
rxnID = update_mega_glove[0][i];
twomol = atom->molecules[reacted_mol[rxnID]];
- // Angles! First let's delete all angle info:
for (int j = 0; j < twomol->natoms; j++) {
int jj = equivalences[j][1][rxnID]-1;
if (atom->map(update_mega_glove[jj+1][i]) < nlocal && atom->map(update_mega_glove[jj+1][i]) >= 0) {
@@ -2086,6 +2089,7 @@ void FixBondReact::update_everything()
angle_atom2[atom->map(update_mega_glove[jj+1][i])][p] == update_mega_glove[nn+1][i] ||
angle_atom3[atom->map(update_mega_glove[jj+1][i])][p] == update_mega_glove[nn+1][i])) {
for (int m = p; m < num_angle[atom->map(update_mega_glove[jj+1][i])]-1; m++) {
+ angle_type[atom->map(update_mega_glove[jj+1][i])][m] = angle_type[atom->map(update_mega_glove[jj+1][i])][m+1];
angle_atom1[atom->map(update_mega_glove[jj+1][i])][m] = angle_atom1[atom->map(update_mega_glove[jj+1][i])][m+1];
angle_atom2[atom->map(update_mega_glove[jj+1][i])][m] = angle_atom2[atom->map(update_mega_glove[jj+1][i])][m+1];
angle_atom3[atom->map(update_mega_glove[jj+1][i])][m] = angle_atom3[atom->map(update_mega_glove[jj+1][i])][m+1];
@@ -2162,6 +2166,7 @@ void FixBondReact::update_everything()
dihedral_atom3[atom->map(update_mega_glove[jj+1][i])][p] == update_mega_glove[nn+1][i] ||
dihedral_atom4[atom->map(update_mega_glove[jj+1][i])][p] == update_mega_glove[nn+1][i])) {
for (int m = p; m < num_dihedral[atom->map(update_mega_glove[jj+1][i])]-1; m++) {
+ dihedral_type[atom->map(update_mega_glove[jj+1][i])][m] = dihedral_type[atom->map(update_mega_glove[jj+1][i])][m+1];
dihedral_atom1[atom->map(update_mega_glove[jj+1][i])][m] = dihedral_atom1[atom->map(update_mega_glove[jj+1][i])][m+1];
dihedral_atom2[atom->map(update_mega_glove[jj+1][i])][m] = dihedral_atom2[atom->map(update_mega_glove[jj+1][i])][m+1];
dihedral_atom3[atom->map(update_mega_glove[jj+1][i])][m] = dihedral_atom3[atom->map(update_mega_glove[jj+1][i])][m+1];
@@ -2242,6 +2247,7 @@ void FixBondReact::update_everything()
improper_atom3[atom->map(update_mega_glove[jj+1][i])][p] == update_mega_glove[nn+1][i] ||
improper_atom4[atom->map(update_mega_glove[jj+1][i])][p] == update_mega_glove[nn+1][i])) {
for (int m = p; m < num_improper[atom->map(update_mega_glove[jj+1][i])]-1; m++) {
+ improper_type[atom->map(update_mega_glove[jj+1][i])][m] = improper_type[atom->map(update_mega_glove[jj+1][i])][m+1];
improper_atom1[atom->map(update_mega_glove[jj+1][i])][m] = improper_atom1[atom->map(update_mega_glove[jj+1][i])][m+1];
improper_atom2[atom->map(update_mega_glove[jj+1][i])][m] = improper_atom2[atom->map(update_mega_glove[jj+1][i])][m+1];
improper_atom3[atom->map(update_mega_glove[jj+1][i])][m] = improper_atom3[atom->map(update_mega_glove[jj+1][i])][m+1];
diff --git a/src/USER-OMP/fix_nvt_sllod_omp.cpp b/src/USER-OMP/fix_nvt_sllod_omp.cpp
index 33445c68af1e8000b669ba2df38db8b26a6a7ee2..6ef1188d53c7fc650377432d11d5b320e3f27352 100644
--- a/src/USER-OMP/fix_nvt_sllod_omp.cpp
+++ b/src/USER-OMP/fix_nvt_sllod_omp.cpp
@@ -31,8 +31,6 @@
using namespace LAMMPS_NS;
using namespace FixConst;
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-
typedef struct { double x,y,z; } dbl3_t;
/* ---------------------------------------------------------------------- */
@@ -85,7 +83,7 @@ void FixNVTSllodOMP::init()
int i;
for (i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
- if (((FixDeform *) modify->fix[i])->remapflag != V_REMAP)
+ if (((FixDeform *) modify->fix[i])->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using fix nvt/sllod/omp with inconsistent fix "
"deform remap option");
break;
diff --git a/src/USER-OMP/pair_lubricate_omp.cpp b/src/USER-OMP/pair_lubricate_omp.cpp
index fbfc24769b73dfddd3853ec77a42530d328a5384..a7f9fc65b59589dfad8df2362184b31727b70f3a 100644
--- a/src/USER-OMP/pair_lubricate_omp.cpp
+++ b/src/USER-OMP/pair_lubricate_omp.cpp
@@ -32,10 +32,6 @@
using namespace LAMMPS_NS;
using namespace MathConst;
-// same as fix_deform.cpp
-
-enum{NO_REMAP,X_REMAP,V_REMAP};
-
// same as fix_wall.cpp
enum{EDGE,CONSTANT,VARIABLE};
diff --git a/src/USER-OMP/pair_lubricate_poly_omp.cpp b/src/USER-OMP/pair_lubricate_poly_omp.cpp
index 825b3e8d705286ed254078ceca8dcd4fddbea070..3f0ef1dba6ca6b51d651691031a5a4137577ebf7 100644
--- a/src/USER-OMP/pair_lubricate_poly_omp.cpp
+++ b/src/USER-OMP/pair_lubricate_poly_omp.cpp
@@ -32,11 +32,6 @@
using namespace LAMMPS_NS;
using namespace MathConst;
-// same as fix_deform.cpp
-
-enum{NO_REMAP,X_REMAP,V_REMAP};
-
-
// same as fix_wall.cpp
enum{EDGE,CONSTANT,VARIABLE};
diff --git a/src/USER-SMD/fix_smd_wall_surface.cpp b/src/USER-SMD/fix_smd_wall_surface.cpp
index 4bd5078cb86558a8699d953a2c8cf79afd24af87..97c2ead5fc79d8b913173d09d709a34294888120 100644
--- a/src/USER-SMD/fix_smd_wall_surface.cpp
+++ b/src/USER-SMD/fix_smd_wall_surface.cpp
@@ -40,10 +40,6 @@ using namespace Eigen;
using namespace std;
#define DELTA 16384
#define EPSILON 1.0e-6
-enum {
- LAYOUT_UNIFORM, LAYOUT_NONUNIFORM, LAYOUT_TILED
-};
-// several files
/* ---------------------------------------------------------------------- */
@@ -151,7 +147,7 @@ void FixSMDWallSurface::setup(int vflag) {
subhi[2] = domain->subhi_lamda[2];
}
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (domain->xperiodic) {
if (comm->myloc[0] == 0)
sublo[0] -= epsilon[0];
diff --git a/src/USER-UEF/dump_cfg_uef.cpp b/src/USER-UEF/dump_cfg_uef.cpp
index f5f78297c2910877fcc5564f8db778426183a553..44af823332c156d6ed520a460f277ae1a25e2bbe 100644
--- a/src/USER-UEF/dump_cfg_uef.cpp
+++ b/src/USER-UEF/dump_cfg_uef.cpp
@@ -30,8 +30,6 @@
using namespace LAMMPS_NS;
-enum{INT,DOUBLE,STRING,BIGINT}; // same as in DumpCustom
-
#define UNWRAPEXPAND 10.0
#define ONEFIELD 32
#define DELTA 1048576
diff --git a/src/USER-VTK/dump_vtk.cpp b/src/USER-VTK/dump_vtk.cpp
index adad070f51596c9d4b28ce0a904fdc476a94ab02..5f2e3d62ec433e4cd77e80852aa1e3f5e88ad913 100644
--- a/src/USER-VTK/dump_vtk.cpp
+++ b/src/USER-VTK/dump_vtk.cpp
@@ -87,7 +87,6 @@ enum{X,Y,Z, // required for vtk, must come first
VARIABLE,COMPUTE,FIX,INAME,DNAME,
ATTRIBUTES}; // must come last
enum{LT,LE,GT,GE,EQ,NEQ};
-enum{INT,DOUBLE,STRING,BIGINT}; // same as in DumpCFG
enum{VTK,VTP,VTU,PVTP,PVTU}; // file formats
#define INVOKED_PERATOM 8
@@ -1091,7 +1090,7 @@ void DumpVTK::buf2arrays(int n, double *mybuf)
vtkAbstractArray *paa = it->second;
if (it->second->GetNumberOfComponents() == 3) {
switch (vtype[it->first]) {
- case INT:
+ case Dump::INT:
{
int iv3[3] = { static_cast<int>(mybuf[iatom*size_one+j ]),
static_cast<int>(mybuf[iatom*size_one+j+1]),
@@ -1100,7 +1099,7 @@ void DumpVTK::buf2arrays(int n, double *mybuf)
pia->InsertNextTupleValue(iv3);
break;
}
- case DOUBLE:
+ case Dump::DOUBLE:
{
vtkDoubleArray *pda = static_cast<vtkDoubleArray*>(paa);
pda->InsertNextTupleValue(&mybuf[iatom*size_one+j]);
@@ -1110,19 +1109,19 @@ void DumpVTK::buf2arrays(int n, double *mybuf)
j+=3;
} else {
switch (vtype[it->first]) {
- case INT:
+ case Dump::INT:
{
vtkIntArray *pia = static_cast<vtkIntArray*>(paa);
pia->InsertNextValue(mybuf[iatom*size_one+j]);
break;
}
- case DOUBLE:
+ case Dump::DOUBLE:
{
vtkDoubleArray *pda = static_cast<vtkDoubleArray*>(paa);
pda->InsertNextValue(mybuf[iatom*size_one+j]);
break;
}
- case STRING:
+ case Dump::STRING:
{
vtkStringArray *psa = static_cast<vtkStringArray*>(paa);
psa->InsertNextValue(typenames[static_cast<int>(mybuf[iatom*size_one+j])]);
@@ -1482,13 +1481,13 @@ void DumpVTK::reset_vtk_data_containers()
++it; ++it; ++it;
for (; it!=vtype.end(); ++it) {
switch(vtype[it->first]) {
- case INT:
+ case Dump::INT:
myarrays[it->first] = vtkSmartPointer<vtkIntArray>::New();
break;
- case DOUBLE:
+ case Dump::DOUBLE:
myarrays[it->first] = vtkSmartPointer<vtkDoubleArray>::New();
break;
- case STRING:
+ case Dump::STRING:
myarrays[it->first] = vtkSmartPointer<vtkStringArray>::New();
break;
}
@@ -1509,13 +1508,13 @@ int DumpVTK::parse_fields(int narg, char **arg)
{
pack_choice[X] = &DumpVTK::pack_x;
- vtype[X] = DOUBLE;
+ vtype[X] = Dump::DOUBLE;
name[X] = "x";
pack_choice[Y] = &DumpVTK::pack_y;
- vtype[Y] = DOUBLE;
+ vtype[Y] = Dump::DOUBLE;
name[Y] = "y";
pack_choice[Z] = &DumpVTK::pack_z;
- vtype[Z] = DOUBLE;
+ vtype[Z] = Dump::DOUBLE;
name[Z] = "z";
// customize by adding to if statement
@@ -1525,33 +1524,33 @@ int DumpVTK::parse_fields(int narg, char **arg)
if (strcmp(arg[iarg],"id") == 0) {
pack_choice[ID] = &DumpVTK::pack_id;
- vtype[ID] = INT;
+ vtype[ID] = Dump::INT;
name[ID] = arg[iarg];
} else if (strcmp(arg[iarg],"mol") == 0) {
if (!atom->molecule_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[MOL] = &DumpVTK::pack_molecule;
- vtype[MOL] = INT;
+ vtype[MOL] = Dump::INT;
name[MOL] = arg[iarg];
} else if (strcmp(arg[iarg],"proc") == 0) {
pack_choice[PROC] = &DumpVTK::pack_proc;
- vtype[PROC] = INT;
+ vtype[PROC] = Dump::INT;
name[PROC] = arg[iarg];
} else if (strcmp(arg[iarg],"procp1") == 0) {
pack_choice[PROCP1] = &DumpVTK::pack_procp1;
- vtype[PROCP1] = INT;
+ vtype[PROCP1] = Dump::INT;
name[PROCP1] = arg[iarg];
} else if (strcmp(arg[iarg],"type") == 0) {
pack_choice[TYPE] = &DumpVTK::pack_type;
- vtype[TYPE] = INT;
+ vtype[TYPE] = Dump::INT;
name[TYPE] =arg[iarg];
} else if (strcmp(arg[iarg],"element") == 0) {
pack_choice[ELEMENT] = &DumpVTK::pack_type;
- vtype[ELEMENT] = STRING;
+ vtype[ELEMENT] = Dump::STRING;
name[ELEMENT] = arg[iarg];
} else if (strcmp(arg[iarg],"mass") == 0) {
pack_choice[MASS] = &DumpVTK::pack_mass;
- vtype[MASS] = DOUBLE;
+ vtype[MASS] = Dump::DOUBLE;
name[MASS] = arg[iarg];
} else if (strcmp(arg[iarg],"x") == 0) {
@@ -1563,181 +1562,181 @@ int DumpVTK::parse_fields(int narg, char **arg)
} else if (strcmp(arg[iarg],"xs") == 0) {
if (domain->triclinic) pack_choice[XS] = &DumpVTK::pack_xs_triclinic;
else pack_choice[XS] = &DumpVTK::pack_xs;
- vtype[XS] = DOUBLE;
+ vtype[XS] = Dump::DOUBLE;
name[XS] = arg[iarg];
} else if (strcmp(arg[iarg],"ys") == 0) {
if (domain->triclinic) pack_choice[YS] = &DumpVTK::pack_ys_triclinic;
else pack_choice[YS] = &DumpVTK::pack_ys;
- vtype[YS] = DOUBLE;
+ vtype[YS] = Dump::DOUBLE;
name[YS] = arg[iarg];
} else if (strcmp(arg[iarg],"zs") == 0) {
if (domain->triclinic) pack_choice[ZS] = &DumpVTK::pack_zs_triclinic;
else pack_choice[ZS] = &DumpVTK::pack_zs;
- vtype[ZS] = DOUBLE;
+ vtype[ZS] = Dump::DOUBLE;
name[ZS] = arg[iarg];
} else if (strcmp(arg[iarg],"xu") == 0) {
if (domain->triclinic) pack_choice[XU] = &DumpVTK::pack_xu_triclinic;
else pack_choice[XU] = &DumpVTK::pack_xu;
- vtype[XU] = DOUBLE;
+ vtype[XU] = Dump::DOUBLE;
name[XU] = arg[iarg];
} else if (strcmp(arg[iarg],"yu") == 0) {
if (domain->triclinic) pack_choice[YU] = &DumpVTK::pack_yu_triclinic;
else pack_choice[YU] = &DumpVTK::pack_yu;
- vtype[YU] = DOUBLE;
+ vtype[YU] = Dump::DOUBLE;
name[YU] = arg[iarg];
} else if (strcmp(arg[iarg],"zu") == 0) {
if (domain->triclinic) pack_choice[ZU] = &DumpVTK::pack_zu_triclinic;
else pack_choice[ZU] = &DumpVTK::pack_zu;
- vtype[ZU] = DOUBLE;
+ vtype[ZU] = Dump::DOUBLE;
name[ZU] = arg[iarg];
} else if (strcmp(arg[iarg],"xsu") == 0) {
if (domain->triclinic) pack_choice[XSU] = &DumpVTK::pack_xsu_triclinic;
else pack_choice[XSU] = &DumpVTK::pack_xsu;
- vtype[XSU] = DOUBLE;
+ vtype[XSU] = Dump::DOUBLE;
name[XSU] = arg[iarg];
} else if (strcmp(arg[iarg],"ysu") == 0) {
if (domain->triclinic) pack_choice[YSU] = &DumpVTK::pack_ysu_triclinic;
else pack_choice[YSU] = &DumpVTK::pack_ysu;
- vtype[YSU] = DOUBLE;
+ vtype[YSU] = Dump::DOUBLE;
name[YSU] = arg[iarg];
} else if (strcmp(arg[iarg],"zsu") == 0) {
if (domain->triclinic) pack_choice[ZSU] = &DumpVTK::pack_zsu_triclinic;
else pack_choice[ZSU] = &DumpVTK::pack_zsu;
- vtype[ZSU] = DOUBLE;
+ vtype[ZSU] = Dump::DOUBLE;
name[ZSU] = arg[iarg];
} else if (strcmp(arg[iarg],"ix") == 0) {
pack_choice[IX] = &DumpVTK::pack_ix;
- vtype[IX] = INT;
+ vtype[IX] = Dump::INT;
name[IX] = arg[iarg];
} else if (strcmp(arg[iarg],"iy") == 0) {
pack_choice[IY] = &DumpVTK::pack_iy;
- vtype[IY] = INT;
+ vtype[IY] = Dump::INT;
name[IY] = arg[iarg];
} else if (strcmp(arg[iarg],"iz") == 0) {
pack_choice[IZ] = &DumpVTK::pack_iz;
- vtype[IZ] = INT;
+ vtype[IZ] = Dump::INT;
name[IZ] = arg[iarg];
} else if (strcmp(arg[iarg],"vx") == 0) {
pack_choice[VX] = &DumpVTK::pack_vx;
- vtype[VX] = DOUBLE;
+ vtype[VX] = Dump::DOUBLE;
name[VX] = arg[iarg];
} else if (strcmp(arg[iarg],"vy") == 0) {
pack_choice[VY] = &DumpVTK::pack_vy;
- vtype[VY] = DOUBLE;
+ vtype[VY] = Dump::DOUBLE;
name[VY] = arg[iarg];
} else if (strcmp(arg[iarg],"vz") == 0) {
pack_choice[VZ] = &DumpVTK::pack_vz;
- vtype[VZ] = DOUBLE;
+ vtype[VZ] = Dump::DOUBLE;
name[VZ] = arg[iarg];
} else if (strcmp(arg[iarg],"fx") == 0) {
pack_choice[FX] = &DumpVTK::pack_fx;
- vtype[FX] = DOUBLE;
+ vtype[FX] = Dump::DOUBLE;
name[FX] = arg[iarg];
} else if (strcmp(arg[iarg],"fy") == 0) {
pack_choice[FY] = &DumpVTK::pack_fy;
- vtype[FY] = DOUBLE;
+ vtype[FY] = Dump::DOUBLE;
name[FY] = arg[iarg];
} else if (strcmp(arg[iarg],"fz") == 0) {
pack_choice[FZ] = &DumpVTK::pack_fz;
- vtype[FZ] = DOUBLE;
+ vtype[FZ] = Dump::DOUBLE;
name[FZ] = arg[iarg];
} else if (strcmp(arg[iarg],"q") == 0) {
if (!atom->q_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[Q] = &DumpVTK::pack_q;
- vtype[Q] = DOUBLE;
+ vtype[Q] = Dump::DOUBLE;
name[Q] = arg[iarg];
} else if (strcmp(arg[iarg],"mux") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[MUX] = &DumpVTK::pack_mux;
- vtype[MUX] = DOUBLE;
+ vtype[MUX] = Dump::DOUBLE;
name[MUX] = arg[iarg];
} else if (strcmp(arg[iarg],"muy") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[MUY] = &DumpVTK::pack_muy;
- vtype[MUY] = DOUBLE;
+ vtype[MUY] = Dump::DOUBLE;
name[MUY] = arg[iarg];
} else if (strcmp(arg[iarg],"muz") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[MUZ] = &DumpVTK::pack_muz;
- vtype[MUZ] = DOUBLE;
+ vtype[MUZ] = Dump::DOUBLE;
name[MUZ] = arg[iarg];
} else if (strcmp(arg[iarg],"mu") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[MU] = &DumpVTK::pack_mu;
- vtype[MU] = DOUBLE;
+ vtype[MU] = Dump::DOUBLE;
name[MU] = arg[iarg];
} else if (strcmp(arg[iarg],"radius") == 0) {
if (!atom->radius_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[RADIUS] = &DumpVTK::pack_radius;
- vtype[RADIUS] = DOUBLE;
+ vtype[RADIUS] = Dump::DOUBLE;
name[RADIUS] = arg[iarg];
} else if (strcmp(arg[iarg],"diameter") == 0) {
if (!atom->radius_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[DIAMETER] = &DumpVTK::pack_diameter;
- vtype[DIAMETER] = DOUBLE;
+ vtype[DIAMETER] = Dump::DOUBLE;
name[DIAMETER] = arg[iarg];
} else if (strcmp(arg[iarg],"omegax") == 0) {
if (!atom->omega_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[OMEGAX] = &DumpVTK::pack_omegax;
- vtype[OMEGAX] = DOUBLE;
+ vtype[OMEGAX] = Dump::DOUBLE;
name[OMEGAX] = arg[iarg];
} else if (strcmp(arg[iarg],"omegay") == 0) {
if (!atom->omega_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[OMEGAY] = &DumpVTK::pack_omegay;
- vtype[OMEGAY] = DOUBLE;
+ vtype[OMEGAY] = Dump::DOUBLE;
name[OMEGAY] = arg[iarg];
} else if (strcmp(arg[iarg],"omegaz") == 0) {
if (!atom->omega_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[OMEGAZ] = &DumpVTK::pack_omegaz;
- vtype[OMEGAZ] = DOUBLE;
+ vtype[OMEGAZ] = Dump::DOUBLE;
name[OMEGAZ] = arg[iarg];
} else if (strcmp(arg[iarg],"angmomx") == 0) {
if (!atom->angmom_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[ANGMOMX] = &DumpVTK::pack_angmomx;
- vtype[ANGMOMX] = DOUBLE;
+ vtype[ANGMOMX] = Dump::DOUBLE;
name[ANGMOMX] = arg[iarg];
} else if (strcmp(arg[iarg],"angmomy") == 0) {
if (!atom->angmom_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[ANGMOMY] = &DumpVTK::pack_angmomy;
- vtype[ANGMOMY] = DOUBLE;
+ vtype[ANGMOMY] = Dump::DOUBLE;
name[ANGMOMY] = arg[iarg];
} else if (strcmp(arg[iarg],"angmomz") == 0) {
if (!atom->angmom_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[ANGMOMZ] = &DumpVTK::pack_angmomz;
- vtype[ANGMOMZ] = DOUBLE;
+ vtype[ANGMOMZ] = Dump::DOUBLE;
name[ANGMOMZ] = arg[iarg];
} else if (strcmp(arg[iarg],"tqx") == 0) {
if (!atom->torque_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[TQX] = &DumpVTK::pack_tqx;
- vtype[TQX] = DOUBLE;
+ vtype[TQX] = Dump::DOUBLE;
name[TQX] = arg[iarg];
} else if (strcmp(arg[iarg],"tqy") == 0) {
if (!atom->torque_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[TQY] = &DumpVTK::pack_tqy;
- vtype[TQY] = DOUBLE;
+ vtype[TQY] = Dump::DOUBLE;
name[TQY] = arg[iarg];
} else if (strcmp(arg[iarg],"tqz") == 0) {
if (!atom->torque_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[TQZ] = &DumpVTK::pack_tqz;
- vtype[TQZ] = DOUBLE;
+ vtype[TQZ] = Dump::DOUBLE;
name[TQZ] = arg[iarg];
// compute value = c_ID
@@ -1745,7 +1744,7 @@ int DumpVTK::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"c_",2) == 0) {
pack_choice[ATTRIBUTES+i] = &DumpVTK::pack_compute;
- vtype[ATTRIBUTES+i] = DOUBLE;
+ vtype[ATTRIBUTES+i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1782,7 +1781,7 @@ int DumpVTK::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"f_",2) == 0) {
pack_choice[ATTRIBUTES+i] = &DumpVTK::pack_fix;
- vtype[ATTRIBUTES+i] = DOUBLE;
+ vtype[ATTRIBUTES+i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1816,7 +1815,7 @@ int DumpVTK::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"v_",2) == 0) {
pack_choice[ATTRIBUTES+i] = &DumpVTK::pack_variable;
- vtype[ATTRIBUTES+i] = DOUBLE;
+ vtype[ATTRIBUTES+i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1837,7 +1836,7 @@ int DumpVTK::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"d_",2) == 0) {
pack_choice[ATTRIBUTES+i] = &DumpVTK::pack_custom;
- vtype[ATTRIBUTES+i] = DOUBLE;
+ vtype[ATTRIBUTES+i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1860,7 +1859,7 @@ int DumpVTK::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"i_",2) == 0) {
pack_choice[ATTRIBUTES+i] = &DumpVTK::pack_custom;
- vtype[ATTRIBUTES+i] = INT;
+ vtype[ATTRIBUTES+i] = Dump::INT;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
diff --git a/src/atom.cpp b/src/atom.cpp
index d5773a2fe1f945ccf96877b9d56d2de857c4dfe1..56bfd04de7195afbdf1eb370ebb36f61b7aa8f52 100644
--- a/src/atom.cpp
+++ b/src/atom.cpp
@@ -51,8 +51,6 @@ using namespace MathConst;
#define DELTA_MEMSTR 1024
#define EPSILON 1.0e-6
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
/* ---------------------------------------------------------------------- */
Atom::Atom(LAMMPS *lmp) : Pointers(lmp)
@@ -866,7 +864,7 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset,
sublo[2] = domain->sublo_lamda[2]; subhi[2] = domain->subhi_lamda[2];
}
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (domain->xperiodic) {
if (comm->myloc[0] == 0) sublo[0] -= epsilon[0];
if (comm->myloc[0] == comm->procgrid[0]-1) subhi[0] += epsilon[0];
diff --git a/src/balance.cpp b/src/balance.cpp
index 9be7b2a7b3267295b49d35c46af1ce007da6dd66..86deb55b47ff716f6e05f1ce0d4adcc3f7d4d793 100644
--- a/src/balance.cpp
+++ b/src/balance.cpp
@@ -48,8 +48,6 @@ using namespace LAMMPS_NS;
enum{XYZ,SHIFT,BISECTION};
enum{NONE,UNIFORM,USER};
enum{X,Y,Z};
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
/* ---------------------------------------------------------------------- */
Balance::Balance(LAMMPS *lmp) : Pointers(lmp)
@@ -281,7 +279,7 @@ void Balance::command(int narg, char **arg)
// no load-balance if imbalance doesn't exceed threshold
// unless switching from tiled to non tiled layout, then force rebalance
- if (comm->layout == LAYOUT_TILED && style != BISECTION) {
+ if (comm->layout == Comm::LAYOUT_TILED && style != BISECTION) {
} else if (imbinit < thresh) return;
// debug output of initial state
@@ -296,16 +294,16 @@ void Balance::command(int narg, char **arg)
// style XYZ = explicit setting of cutting planes of logical 3d grid
if (style == XYZ) {
- if (comm->layout == LAYOUT_UNIFORM) {
+ if (comm->layout == Comm::LAYOUT_UNIFORM) {
if (xflag == USER || yflag == USER || zflag == USER)
- comm->layout = LAYOUT_NONUNIFORM;
- } else if (comm->layout == LAYOUT_NONUNIFORM) {
+ comm->layout = Comm::LAYOUT_NONUNIFORM;
+ } else if (comm->layout == Comm::LAYOUT_NONUNIFORM) {
if (xflag == UNIFORM && yflag == UNIFORM && zflag == UNIFORM)
- comm->layout = LAYOUT_UNIFORM;
- } else if (comm->layout == LAYOUT_TILED) {
+ comm->layout = Comm::LAYOUT_UNIFORM;
+ } else if (comm->layout == Comm::LAYOUT_TILED) {
if (xflag == UNIFORM && yflag == UNIFORM && zflag == UNIFORM)
- comm->layout = LAYOUT_UNIFORM;
- else comm->layout = LAYOUT_NONUNIFORM;
+ comm->layout = Comm::LAYOUT_UNIFORM;
+ else comm->layout = Comm::LAYOUT_NONUNIFORM;
}
if (xflag == UNIFORM) {
@@ -333,7 +331,7 @@ void Balance::command(int narg, char **arg)
// style SHIFT = adjust cutting planes of logical 3d grid
if (style == SHIFT) {
- comm->layout = LAYOUT_NONUNIFORM;
+ comm->layout = Comm::LAYOUT_NONUNIFORM;
shift_setup_static(bstr);
niter = shift();
}
@@ -341,7 +339,7 @@ void Balance::command(int narg, char **arg)
// style BISECTION = recursive coordinate bisectioning
if (style == BISECTION) {
- comm->layout = LAYOUT_TILED;
+ comm->layout = Comm::LAYOUT_TILED;
bisection(1);
}
@@ -745,7 +743,7 @@ void Balance::shift_setup_static(char *str)
// if current layout is TILED, set initial uniform splits in Comm
// this gives starting point to subsequent shift balancing
- if (comm->layout == LAYOUT_TILED) {
+ if (comm->layout == Comm::LAYOUT_TILED) {
int *procgrid = comm->procgrid;
double *xsplit = comm->xsplit;
double *ysplit = comm->ysplit;
diff --git a/src/comm.cpp b/src/comm.cpp
index b53ce251f27bce20a5c964015f52286e9a266a74..f355a562fc03e3089c35bce57c0ceae4c9de95a3 100644
--- a/src/comm.cpp
+++ b/src/comm.cpp
@@ -40,11 +40,8 @@ using namespace LAMMPS_NS;
#define BUFMIN 1000 // also in comm styles
-enum{SINGLE,MULTI}; // same as in Comm sub-styles
-enum{MULTIPLE}; // same as in ProcMap
enum{ONELEVEL,TWOLEVEL,NUMA,CUSTOM};
enum{CART,CARTREORDER,XYZ};
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
/* ---------------------------------------------------------------------- */
@@ -244,14 +241,14 @@ void Comm::modify_params(int narg, char **arg)
if (iarg+2 > narg) error->all(FLERR,"Illegal comm_modify command");
if (strcmp(arg[iarg+1],"single") == 0) {
// need to reset cutghostuser when switching comm mode
- if (mode == MULTI) cutghostuser = 0.0;
+ if (mode == Comm::MULTI) cutghostuser = 0.0;
memory->destroy(cutusermulti);
cutusermulti = NULL;
- mode = SINGLE;
+ mode = Comm::SINGLE;
} else if (strcmp(arg[iarg+1],"multi") == 0) {
// need to reset cutghostuser when switching comm mode
- if (mode == SINGLE) cutghostuser = 0.0;
- mode = MULTI;
+ if (mode == Comm::SINGLE) cutghostuser = 0.0;
+ mode = Comm::MULTI;
} else error->all(FLERR,"Illegal comm_modify command");
iarg += 2;
} else if (strcmp(arg[iarg],"group") == 0) {
@@ -265,7 +262,7 @@ void Comm::modify_params(int narg, char **arg)
iarg += 2;
} else if (strcmp(arg[iarg],"cutoff") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal comm_modify command");
- if (mode == MULTI)
+ if (mode == Comm::MULTI)
error->all(FLERR,
"Use cutoff/multi keyword to set cutoff in multi mode");
cutghostuser = force->numeric(FLERR,arg[iarg+1]);
@@ -275,7 +272,7 @@ void Comm::modify_params(int narg, char **arg)
} else if (strcmp(arg[iarg],"cutoff/multi") == 0) {
int i,nlo,nhi;
double cut;
- if (mode == SINGLE)
+ if (mode == Comm::SINGLE)
error->all(FLERR,"Use cutoff keyword to set cutoff in single mode");
if (domain->box_exist == 0)
error->all(FLERR,
@@ -415,7 +412,7 @@ void Comm::set_processors(int narg, char **arg)
if (strcmp(arg[iarg+3],"multiple") == 0) {
if (universe->iworld == irecv-1) {
otherflag = 1;
- other_style = MULTIPLE;
+ other_style = Comm::MULTIPLE;
}
} else error->all(FLERR,"Illegal processors command");
iarg += 4;
@@ -607,7 +604,7 @@ int Comm::coord2proc(double *x, int &igx, int &igy, int &igz)
triclinic = domain->triclinic;
- if (layout == LAYOUT_UNIFORM) {
+ if (layout == Comm::LAYOUT_UNIFORM) {
if (triclinic == 0) {
igx = static_cast<int> (procgrid[0] * (x[0]-boxlo[0]) / prd[0]);
igy = static_cast<int> (procgrid[1] * (x[1]-boxlo[1]) / prd[1]);
@@ -618,7 +615,7 @@ int Comm::coord2proc(double *x, int &igx, int &igy, int &igz)
igz = static_cast<int> (procgrid[2] * x[2]);
}
- } else if (layout == LAYOUT_NONUNIFORM) {
+ } else if (layout == Comm::LAYOUT_NONUNIFORM) {
if (triclinic == 0) {
igx = binary((x[0]-boxlo[0])/prd[0],procgrid[0],xsplit);
igy = binary((x[1]-boxlo[1])/prd[1],procgrid[1],ysplit);
diff --git a/src/comm.h b/src/comm.h
index aefcb14391c71df20f2157e11bd71e76c4b467ef..298f435c3d94f2e04af2077cfadb21ef3372896c 100644
--- a/src/comm.h
+++ b/src/comm.h
@@ -24,7 +24,9 @@ class Comm : protected Pointers {
int layout; // LAYOUT_UNIFORM = equal-sized bricks
// LAYOUT_NONUNIFORM = logical bricks, but diff sizes via LB
// LAYOUT_TILED = general tiling, due to RCB LB
+ enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED};
int mode; // 0 = single cutoff, 1 = multi-type cutoff
+ enum{SINGLE,MULTI};
int me,nprocs; // proc info
int ghost_velocity; // 1 if ghost atoms have velocity, 0 if not
@@ -137,6 +139,8 @@ class Comm : protected Pointers {
int ncores; // # of cores per node
int coregrid[3]; // 3d grid of cores within a node
int user_coregrid[3]; // user request for cores in each dim
+ public:
+ enum{MULTIPLE};
};
}
diff --git a/src/comm_brick.cpp b/src/comm_brick.cpp
index 070c93bd3cb765dffe29ecab6dd6cfa42742579e..3aa6d7748ce96b2c095857ed5e9f6c2acbc9314d 100644
--- a/src/comm_brick.cpp
+++ b/src/comm_brick.cpp
@@ -46,9 +46,6 @@ using namespace LAMMPS_NS;
#define BUFEXTRA 1000
#define BIG 1.0e20
-enum{SINGLE,MULTI}; // same as in Comm
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
/* ---------------------------------------------------------------------- */
CommBrick::CommBrick(LAMMPS *lmp) :
@@ -61,7 +58,7 @@ CommBrick::CommBrick(LAMMPS *lmp) :
sendlist(NULL), maxsendlist(NULL), buf_send(NULL), buf_recv(NULL)
{
style = 0;
- layout = LAYOUT_UNIFORM;
+ layout = Comm::LAYOUT_UNIFORM;
pbc_flag = NULL;
init_buffers();
}
@@ -71,7 +68,7 @@ CommBrick::CommBrick(LAMMPS *lmp) :
CommBrick::~CommBrick()
{
free_swap();
- if (mode == MULTI) {
+ if (mode == Comm::MULTI) {
free_multi();
memory->destroy(cutghostmulti);
}
@@ -93,7 +90,7 @@ CommBrick::~CommBrick()
CommBrick::CommBrick(LAMMPS *lmp, Comm *oldcomm) : Comm(*oldcomm)
{
- if (oldcomm->layout == LAYOUT_TILED)
+ if (oldcomm->layout == Comm::LAYOUT_TILED)
error->all(FLERR,"Cannot change to comm_style brick from tiled layout");
style = 0;
@@ -144,11 +141,11 @@ void CommBrick::init()
// memory for multi-style communication
- if (mode == MULTI && multilo == NULL) {
+ if (mode == Comm::MULTI && multilo == NULL) {
allocate_multi(maxswap);
memory->create(cutghostmulti,atom->ntypes+1,3,"comm:cutghostmulti");
}
- if (mode == SINGLE && multilo) {
+ if (mode == Comm::SINGLE && multilo) {
free_multi();
memory->destroy(cutghostmulti);
}
@@ -184,7 +181,7 @@ void CommBrick::setup()
subhi = domain->subhi;
cutghost[0] = cutghost[1] = cutghost[2] = cut;
- if (mode == MULTI) {
+ if (mode == Comm::MULTI) {
double *cuttype = neighbor->cuttype;
for (i = 1; i <= ntypes; i++) {
cut = 0.0;
@@ -208,7 +205,7 @@ void CommBrick::setup()
length2 = h_inv[2];
cutghost[2] = cut * length2;
- if (mode == MULTI) {
+ if (mode == Comm::MULTI) {
double *cuttype = neighbor->cuttype;
for (i = 1; i <= ntypes; i++) {
cut = 0.0;
@@ -240,7 +237,7 @@ void CommBrick::setup()
int *periodicity = domain->periodicity;
int left,right;
- if (layout == LAYOUT_UNIFORM) {
+ if (layout == Comm::LAYOUT_UNIFORM) {
maxneed[0] = static_cast<int> (cutghost[0] * procgrid[0] / prd[0]) + 1;
maxneed[1] = static_cast<int> (cutghost[1] * procgrid[1] / prd[1]) + 1;
maxneed[2] = static_cast<int> (cutghost[2] * procgrid[2] / prd[2]) + 1;
@@ -357,7 +354,7 @@ void CommBrick::setup()
if (ineed % 2 == 0) {
sendproc[iswap] = procneigh[dim][0];
recvproc[iswap] = procneigh[dim][1];
- if (mode == SINGLE) {
+ if (mode == Comm::SINGLE) {
if (ineed < 2) slablo[iswap] = -BIG;
else slablo[iswap] = 0.5 * (sublo[dim] + subhi[dim]);
slabhi[iswap] = sublo[dim] + cutghost[dim];
@@ -380,7 +377,7 @@ void CommBrick::setup()
} else {
sendproc[iswap] = procneigh[dim][1];
recvproc[iswap] = procneigh[dim][0];
- if (mode == SINGLE) {
+ if (mode == Comm::SINGLE) {
slablo[iswap] = subhi[dim] - cutghost[dim];
if (ineed < 2) slabhi[iswap] = BIG;
else slabhi[iswap] = 0.5 * (sublo[dim] + subhi[dim]);
@@ -737,7 +734,7 @@ void CommBrick::borders()
// store sent atom indices in sendlist for use in future timesteps
x = atom->x;
- if (mode == SINGLE) {
+ if (mode == Comm::SINGLE) {
lo = slablo[iswap];
hi = slabhi[iswap];
} else {
@@ -766,7 +763,7 @@ void CommBrick::borders()
if (sendflag) {
if (!bordergroup || ineed >= 2) {
- if (mode == SINGLE) {
+ if (mode == Comm::SINGLE) {
for (i = nfirst; i < nlast; i++)
if (x[i][dim] >= lo && x[i][dim] <= hi) {
if (nsend == maxsendlist[iswap]) grow_list(iswap,nsend);
@@ -783,7 +780,7 @@ void CommBrick::borders()
}
} else {
- if (mode == SINGLE) {
+ if (mode == Comm::SINGLE) {
ngroup = atom->nfirst;
for (i = 0; i < ngroup; i++)
if (x[i][dim] >= lo && x[i][dim] <= hi) {
@@ -1396,7 +1393,7 @@ void CommBrick::grow_swap(int n)
{
free_swap();
allocate_swap(n);
- if (mode == MULTI) {
+ if (mode == Comm::MULTI) {
free_multi();
allocate_multi(n);
}
diff --git a/src/comm_tiled.cpp b/src/comm_tiled.cpp
index 12b10be2f7652f113c5e736a80d032a8d7a9942d..584be94bf3f32f8bc91451f545fef888e3f54bc6 100644
--- a/src/comm_tiled.cpp
+++ b/src/comm_tiled.cpp
@@ -38,15 +38,12 @@ using namespace LAMMPS_NS;
#define DELTA_PROCS 16
-enum{SINGLE,MULTI}; // same as in Comm
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
/* ---------------------------------------------------------------------- */
CommTiled::CommTiled(LAMMPS *lmp) : Comm(lmp)
{
style = 1;
- layout = LAYOUT_UNIFORM;
+ layout = Comm::LAYOUT_UNIFORM;
pbc_flag = NULL;
init_buffers();
}
@@ -115,7 +112,7 @@ void CommTiled::init()
if (triclinic)
error->all(FLERR,"Cannot yet use comm_style tiled with triclinic box");
- if (mode == MULTI)
+ if (mode == Comm::MULTI)
error->all(FLERR,"Cannot yet use comm_style tiled with multi-mode comm");
}
@@ -141,7 +138,7 @@ void CommTiled::setup()
// set function pointers
- if (layout != LAYOUT_TILED) {
+ if (layout != Comm::LAYOUT_TILED) {
box_drop = &CommTiled::box_drop_brick;
box_other = &CommTiled::box_other_brick;
box_touch = &CommTiled::box_touch_brick;
@@ -155,7 +152,7 @@ void CommTiled::setup()
// if RCB decomp exists and just changed, gather needed global RCB info
- if (layout == LAYOUT_TILED) coord2proc_setup();
+ if (layout == Comm::LAYOUT_TILED) coord2proc_setup();
// set cutoff for comm forward and comm reverse
// check that cutoff < any periodic box length
@@ -1807,7 +1804,7 @@ void CommTiled::coord2proc_setup()
int CommTiled::coord2proc(double *x, int &igx, int &igy, int &igz)
{
- if (layout != LAYOUT_TILED) return Comm::coord2proc(x,igx,igy,igz);
+ if (layout != Comm::LAYOUT_TILED) return Comm::coord2proc(x,igx,igy,igz);
return point_drop_tiled_recurse(x,0,nprocs-1);
}
diff --git a/src/compute_temp_deform.cpp b/src/compute_temp_deform.cpp
index 907f7090987a5cb4c0b3261835aaec74384479c7..8729204e700d5bd105ed4ac35f4468068dfb8547 100644
--- a/src/compute_temp_deform.cpp
+++ b/src/compute_temp_deform.cpp
@@ -32,8 +32,6 @@
using namespace LAMMPS_NS;
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-
/* ---------------------------------------------------------------------- */
ComputeTempDeform::ComputeTempDeform(LAMMPS *lmp, int narg, char **arg) :
@@ -71,7 +69,7 @@ void ComputeTempDeform::init()
for (i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
- if (((FixDeform *) modify->fix[i])->remapflag == X_REMAP &&
+ if (((FixDeform *) modify->fix[i])->remapflag == Domain::X_REMAP &&
comm->me == 0)
error->warning(FLERR,"Using compute temp/deform with inconsistent "
"fix deform remap option");
diff --git a/src/create_atoms.cpp b/src/create_atoms.cpp
index cdd2c5e37f4528da2f7844ba0462c608fad18b62..dcc36aa5d433e889d1baa5b7c09302aa66416ce9 100644
--- a/src/create_atoms.cpp
+++ b/src/create_atoms.cpp
@@ -44,7 +44,6 @@ using namespace MathConst;
enum{BOX,REGION,SINGLE,RANDOM};
enum{ATOM,MOLECULE};
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
/* ---------------------------------------------------------------------- */
@@ -314,7 +313,7 @@ void CreateAtoms::command(int narg, char **arg)
}
if (style == BOX || style == REGION) {
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (domain->xperiodic) {
if (comm->myloc[0] == 0) sublo[0] -= epsilon[0];
if (comm->myloc[0] == comm->procgrid[0]-1) subhi[0] -= 2.0*epsilon[0];
diff --git a/src/domain.cpp b/src/domain.cpp
index eae4f3926df72c1c587785e78db2f79efdb12bcf..4c1d12018b1b6076819e3938239e70129c61d989 100644
--- a/src/domain.cpp
+++ b/src/domain.cpp
@@ -44,10 +44,6 @@
using namespace LAMMPS_NS;
using namespace MathConst;
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-enum{IGNORE,WARN,ERROR}; // same as thermo.cpp
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
#define BIG 1.0e20
#define SMALL 1.0e-4
#define DELTAREGION 4
@@ -155,7 +151,7 @@ void Domain::init()
for (int i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
deform_flag = 1;
- if (((FixDeform *) modify->fix[i])->remapflag == V_REMAP) {
+ if (((FixDeform *) modify->fix[i])->remapflag == Domain::V_REMAP) {
deform_vremap = 1;
deform_groupbit = modify->fix[i]->groupbit;
}
@@ -276,7 +272,7 @@ void Domain::set_global_box()
void Domain::set_lamda_box()
{
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
int *myloc = comm->myloc;
double *xsplit = comm->xsplit;
double *ysplit = comm->ysplit;
@@ -313,7 +309,7 @@ void Domain::set_local_box()
{
if (triclinic) return;
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
int *myloc = comm->myloc;
int *procgrid = comm->procgrid;
double *xsplit = comm->xsplit;
@@ -762,7 +758,7 @@ void Domain::image_check()
if (k == -1) {
nmissing++;
- if (lostbond == ERROR)
+ if (lostbond == Thermo::ERROR)
error->one(FLERR,"Bond atom missing in image check");
continue;
}
@@ -785,7 +781,7 @@ void Domain::image_check()
if (flagall && comm->me == 0)
error->warning(FLERR,"Inconsistent image flags");
- if (lostbond == WARN) {
+ if (lostbond == Thermo::WARN) {
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
if (all && comm->me == 0)
@@ -861,7 +857,7 @@ void Domain::box_too_small_check()
if (k == -1) {
nmissing++;
- if (lostbond == ERROR)
+ if (lostbond == Thermo::ERROR)
error->one(FLERR,"Bond atom missing in box size check");
continue;
}
@@ -875,7 +871,7 @@ void Domain::box_too_small_check()
}
}
- if (lostbond == WARN) {
+ if (lostbond == Thermo::WARN) {
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
if (all && comm->me == 0)
diff --git a/src/domain.h b/src/domain.h
index fe2502c691476c32ebed7448ed035b4533f9d608..5581f9a45fc85eb2a5d79f353e5bc53e41ba2e89 100644
--- a/src/domain.h
+++ b/src/domain.h
@@ -93,6 +93,7 @@ class Domain : protected Pointers {
class Region **regions; // list of defined Regions
int copymode;
+ enum{NO_REMAP,X_REMAP,V_REMAP};
typedef Region *(*RegionCreator)(LAMMPS *,int,char**);
typedef std::map<std::string,RegionCreator> RegionCreatorMap;
diff --git a/src/dump.cpp b/src/dump.cpp
index 46e556600a8473ffa80feca8bc90a28e08a788f1..7c171015bb4e201f8a8602294ad53c2c9937fbf6 100644
--- a/src/dump.cpp
+++ b/src/dump.cpp
@@ -91,6 +91,11 @@ Dump::Dump(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp)
pbcflag = 0;
delay_flag = 0;
+ maxfiles = -1;
+ numfiles = 0;
+ fileidx = 0;
+ nameslist = NULL;
+
maxbuf = maxids = maxsort = maxproc = 0;
buf = bufsort = NULL;
ids = idsort = NULL;
@@ -187,6 +192,14 @@ Dump::~Dump()
if (multiproc) MPI_Comm_free(&clustercomm);
+ // delete storage for caching file names
+
+ if (maxfiles > 0) {
+ for (int idx=0; idx < numfiles; ++idx)
+ delete[] nameslist[idx];
+ delete[] nameslist;
+ }
+
// XTC style sets fp to NULL since it closes file in its destructor
if (multifile == 0 && fp != NULL) {
@@ -549,6 +562,19 @@ void Dump::openfile()
sprintf(filecurrent,pad,filestar,update->ntimestep,ptr+1);
}
*ptr = '*';
+ if (maxfiles > 0) {
+ if (numfiles < maxfiles) {
+ nameslist[numfiles] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[numfiles],filecurrent);
+ ++numfiles;
+ } else {
+ remove(nameslist[fileidx]);
+ delete[] nameslist[fileidx];
+ nameslist[fileidx] = new char[strlen(filecurrent)+1];
+ strcpy(nameslist[fileidx],filecurrent);
+ fileidx = (fileidx + 1) % maxfiles;
+ }
+ }
}
// each proc with filewriter = 1 opens a file
@@ -1003,6 +1029,27 @@ void Dump::modify_params(int narg, char **arg)
iarg += n;
}
+ } else if (strcmp(arg[iarg],"maxfiles") == 0) {
+ if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
+ if (!multifile)
+ error->all(FLERR,"Cannot use dump_modify maxfiles "
+ "without * in dump file name");
+ // wipe out existing storage
+ if (maxfiles > 0) {
+ for (int idx=0; idx < numfiles; ++idx)
+ delete[] nameslist[idx];
+ delete[] nameslist;
+ }
+ maxfiles = force->inumeric(FLERR,arg[iarg+1]);
+ if (maxfiles == 0) error->all(FLERR,"Illegal dump_modify command");
+ if (maxfiles > 0) {
+ nameslist = new char*[maxfiles];
+ numfiles = 0;
+ for (int idx=0; idx < maxfiles; ++idx)
+ nameslist[idx] = NULL;
+ fileidx = 0;
+ }
+ iarg += 2;
} else if (strcmp(arg[iarg],"nfile") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
if (!multiproc)
diff --git a/src/dump.h b/src/dump.h
index b726b9c3e73d764cbe1bfbaebbc2791f3fc52347..1c6a131f76088e060cf5cb697e7ef744d660d95e 100644
--- a/src/dump.h
+++ b/src/dump.h
@@ -94,6 +94,7 @@ class Dump : protected Pointers {
char *format_int_user;
char *format_bigint_user;
char **format_column_user;
+ enum{INT,DOUBLE,STRING,BIGINT};
FILE *fp; // file to write dump to
int size_one; // # of quantities for one atom
@@ -105,6 +106,11 @@ class Dump : protected Pointers {
double boxzlo,boxzhi;
double boxxy,boxxz,boxyz;
+ int maxfiles; // max number of files created, -1 == infinite
+ int numfiles; // number of files in names list
+ int fileidx; // index of file in names list
+ char **nameslist; // list of history file names
+
bigint ntotal; // total # of per-atom lines in snapshot
int reorderflag; // 1 if OK to reorder instead of sort
int ntotal_reorder; // # of atoms that must be in snapshot
diff --git a/src/dump_cfg.cpp b/src/dump_cfg.cpp
index 15f64688464ad35e20952de8adc746f787a42302..ddd662c8a65c3ae0365672ed3e82f27c7587be52 100644
--- a/src/dump_cfg.cpp
+++ b/src/dump_cfg.cpp
@@ -33,8 +33,6 @@
using namespace LAMMPS_NS;
-enum{INT,DOUBLE,STRING,BIGINT}; // same as in DumpCustom
-
#define UNWRAPEXPAND 10.0
#define ONEFIELD 32
#define DELTA 1048576
@@ -181,15 +179,15 @@ int DumpCFG::convert_string(int n, double *mybuf)
} else if (j == 1) {
offset += sprintf(&sbuf[offset],"%s \n",typenames[(int) mybuf[m]]);
} else if (j >= 2) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
offset +=
sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
offset +=
sprintf(&sbuf[offset],vformat[j],typenames[(int) mybuf[m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
offset +=
sprintf(&sbuf[offset],vformat[j],static_cast<bigint> (mybuf[m]));
}
@@ -216,15 +214,15 @@ int DumpCFG::convert_string(int n, double *mybuf)
unwrap_coord = (mybuf[m] - 0.5)/UNWRAPEXPAND + 0.5;
offset += sprintf(&sbuf[offset],vformat[j],unwrap_coord);
} else if (j >= 5 ) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
offset +=
sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
offset +=
sprintf(&sbuf[offset],vformat[j],typenames[(int) mybuf[m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
offset +=
sprintf(&sbuf[offset],vformat[j],static_cast<bigint> (mybuf[m]));
}
@@ -266,13 +264,13 @@ void DumpCFG::write_lines(int n, double *mybuf)
} else if (j == 1) {
fprintf(fp,"%s \n",typenames[(int) mybuf[m]]);
} else if (j >= 2) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
fprintf(fp,vformat[j],static_cast<int> (mybuf[m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
fprintf(fp,vformat[j],mybuf[m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
fprintf(fp,vformat[j],typenames[(int) mybuf[m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
fprintf(fp,vformat[j],static_cast<bigint> (mybuf[m]));
}
m++;
@@ -292,13 +290,13 @@ void DumpCFG::write_lines(int n, double *mybuf)
unwrap_coord = (mybuf[m] - 0.5)/UNWRAPEXPAND + 0.5;
fprintf(fp,vformat[j],unwrap_coord);
} else if (j >= 5 ) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
fprintf(fp,vformat[j],static_cast<int> (mybuf[m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
fprintf(fp,vformat[j],mybuf[m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
fprintf(fp,vformat[j],typenames[(int) mybuf[m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
fprintf(fp,vformat[j],static_cast<bigint> (mybuf[m]));
}
m++;
diff --git a/src/dump_custom.cpp b/src/dump_custom.cpp
index 174fcd73c875ac04f308d612d81ce25835712f4e..40d74a1d746837209696e13f6aecfd9ea6deaf54 100644
--- a/src/dump_custom.cpp
+++ b/src/dump_custom.cpp
@@ -45,7 +45,6 @@ enum{ID,MOL,PROC,PROCP1,TYPE,ELEMENT,MASS,
TQX,TQY,TQZ,
COMPUTE,FIX,VARIABLE,INAME,DNAME};
enum{LT,LE,GT,GE,EQ,NEQ,XOR};
-enum{INT,DOUBLE,STRING,BIGINT}; // same as in DumpCFG
#define INVOKED_PERATOM 8
#define ONEFIELD 32
@@ -60,7 +59,7 @@ DumpCustom::DumpCustom(LAMMPS *lmp, int narg, char **arg) :
earg(NULL), vtype(NULL), vformat(NULL), columns(NULL), choose(NULL),
dchoose(NULL), clist(NULL), field2index(NULL), argindex(NULL), id_compute(NULL),
compute(NULL), id_fix(NULL), fix(NULL), id_variable(NULL), variable(NULL),
- vbuf(NULL), id_custom(NULL), flag_custom(NULL), typenames(NULL),
+ vbuf(NULL), id_custom(NULL), flag_custom(NULL), typenames(NULL),
pack_choice(NULL)
{
if (narg == 5) error->all(FLERR,"No dump custom arguments specified");
@@ -165,10 +164,10 @@ DumpCustom::DumpCustom(LAMMPS *lmp, int narg, char **arg) :
format_default[0] = '\0';
for (int i = 0; i < size_one; i++) {
- if (vtype[i] == INT) strcat(format_default,"%d ");
- else if (vtype[i] == DOUBLE) strcat(format_default,"%g ");
- else if (vtype[i] == STRING) strcat(format_default,"%s ");
- else if (vtype[i] == BIGINT) strcat(format_default,BIGINT_FORMAT " ");
+ if (vtype[i] == Dump::INT) strcat(format_default,"%d ");
+ else if (vtype[i] == Dump::DOUBLE) strcat(format_default,"%g ");
+ else if (vtype[i] == Dump::STRING) strcat(format_default,"%s ");
+ else if (vtype[i] == Dump::BIGINT) strcat(format_default,BIGINT_FORMAT " ");
vformat[i] = NULL;
}
@@ -288,13 +287,13 @@ void DumpCustom::init_style()
if (format_column_user[i]) {
vformat[i] = new char[strlen(format_column_user[i]) + 2];
strcpy(vformat[i],format_column_user[i]);
- } else if (vtype[i] == INT && format_int_user) {
+ } else if (vtype[i] == Dump::INT && format_int_user) {
vformat[i] = new char[strlen(format_int_user) + 2];
strcpy(vformat[i],format_int_user);
- } else if (vtype[i] == DOUBLE && format_float_user) {
+ } else if (vtype[i] == Dump::DOUBLE && format_float_user) {
vformat[i] = new char[strlen(format_float_user) + 2];
strcpy(vformat[i],format_float_user);
- } else if (vtype[i] == BIGINT && format_bigint_user) {
+ } else if (vtype[i] == Dump::BIGINT && format_bigint_user) {
vformat[i] = new char[strlen(format_bigint_user) + 2];
strcpy(vformat[i],format_bigint_user);
} else {
@@ -1088,13 +1087,13 @@ int DumpCustom::convert_string(int n, double *mybuf)
}
for (j = 0; j < size_one; j++) {
- if (vtype[j] == INT)
+ if (vtype[j] == Dump::INT)
offset += sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
- else if (vtype[j] == DOUBLE)
+ else if (vtype[j] == Dump::DOUBLE)
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
- else if (vtype[j] == STRING)
+ else if (vtype[j] == Dump::STRING)
offset += sprintf(&sbuf[offset],vformat[j],typenames[(int) mybuf[m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
offset += sprintf(&sbuf[offset],vformat[j],
static_cast<bigint> (mybuf[m]));
m++;
@@ -1137,11 +1136,11 @@ void DumpCustom::write_lines(int n, double *mybuf)
int m = 0;
for (i = 0; i < n; i++) {
for (j = 0; j < size_one; j++) {
- if (vtype[j] == INT) fprintf(fp,vformat[j],static_cast<int> (mybuf[m]));
- else if (vtype[j] == DOUBLE) fprintf(fp,vformat[j],mybuf[m]);
- else if (vtype[j] == STRING)
+ if (vtype[j] == Dump::INT) fprintf(fp,vformat[j],static_cast<int> (mybuf[m]));
+ else if (vtype[j] == Dump::DOUBLE) fprintf(fp,vformat[j],mybuf[m]);
+ else if (vtype[j] == Dump::STRING)
fprintf(fp,vformat[j],typenames[(int) mybuf[m]]);
- else if (vtype[j] == BIGINT)
+ else if (vtype[j] == Dump::BIGINT)
fprintf(fp,vformat[j],static_cast<bigint> (mybuf[m]));
m++;
}
@@ -1161,192 +1160,192 @@ int DumpCustom::parse_fields(int narg, char **arg)
if (strcmp(arg[iarg],"id") == 0) {
pack_choice[i] = &DumpCustom::pack_id;
- if (sizeof(tagint) == sizeof(smallint)) vtype[i] = INT;
- else vtype[i] = BIGINT;
+ if (sizeof(tagint) == sizeof(smallint)) vtype[i] = Dump::INT;
+ else vtype[i] = Dump::BIGINT;
} else if (strcmp(arg[iarg],"mol") == 0) {
if (!atom->molecule_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_molecule;
- if (sizeof(tagint) == sizeof(smallint)) vtype[i] = INT;
- else vtype[i] = BIGINT;
+ if (sizeof(tagint) == sizeof(smallint)) vtype[i] = Dump::INT;
+ else vtype[i] = Dump::BIGINT;
} else if (strcmp(arg[iarg],"proc") == 0) {
pack_choice[i] = &DumpCustom::pack_proc;
- vtype[i] = INT;
+ vtype[i] = Dump::INT;
} else if (strcmp(arg[iarg],"procp1") == 0) {
pack_choice[i] = &DumpCustom::pack_procp1;
- vtype[i] = INT;
+ vtype[i] = Dump::INT;
} else if (strcmp(arg[iarg],"type") == 0) {
pack_choice[i] = &DumpCustom::pack_type;
- vtype[i] = INT;
+ vtype[i] = Dump::INT;
} else if (strcmp(arg[iarg],"element") == 0) {
pack_choice[i] = &DumpCustom::pack_type;
- vtype[i] = STRING;
+ vtype[i] = Dump::STRING;
} else if (strcmp(arg[iarg],"mass") == 0) {
pack_choice[i] = &DumpCustom::pack_mass;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"x") == 0) {
pack_choice[i] = &DumpCustom::pack_x;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"y") == 0) {
pack_choice[i] = &DumpCustom::pack_y;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"z") == 0) {
pack_choice[i] = &DumpCustom::pack_z;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"xs") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_xs_triclinic;
else pack_choice[i] = &DumpCustom::pack_xs;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"ys") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_ys_triclinic;
else pack_choice[i] = &DumpCustom::pack_ys;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"zs") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_zs_triclinic;
else pack_choice[i] = &DumpCustom::pack_zs;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"xu") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_xu_triclinic;
else pack_choice[i] = &DumpCustom::pack_xu;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"yu") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_yu_triclinic;
else pack_choice[i] = &DumpCustom::pack_yu;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"zu") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_zu_triclinic;
else pack_choice[i] = &DumpCustom::pack_zu;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"xsu") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_xsu_triclinic;
else pack_choice[i] = &DumpCustom::pack_xsu;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"ysu") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_ysu_triclinic;
else pack_choice[i] = &DumpCustom::pack_ysu;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"zsu") == 0) {
if (domain->triclinic) pack_choice[i] = &DumpCustom::pack_zsu_triclinic;
else pack_choice[i] = &DumpCustom::pack_zsu;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"ix") == 0) {
pack_choice[i] = &DumpCustom::pack_ix;
- vtype[i] = INT;
+ vtype[i] = Dump::INT;
} else if (strcmp(arg[iarg],"iy") == 0) {
pack_choice[i] = &DumpCustom::pack_iy;
- vtype[i] = INT;
+ vtype[i] = Dump::INT;
} else if (strcmp(arg[iarg],"iz") == 0) {
pack_choice[i] = &DumpCustom::pack_iz;
- vtype[i] = INT;
+ vtype[i] = Dump::INT;
} else if (strcmp(arg[iarg],"vx") == 0) {
pack_choice[i] = &DumpCustom::pack_vx;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"vy") == 0) {
pack_choice[i] = &DumpCustom::pack_vy;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"vz") == 0) {
pack_choice[i] = &DumpCustom::pack_vz;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"fx") == 0) {
pack_choice[i] = &DumpCustom::pack_fx;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"fy") == 0) {
pack_choice[i] = &DumpCustom::pack_fy;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"fz") == 0) {
pack_choice[i] = &DumpCustom::pack_fz;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"q") == 0) {
if (!atom->q_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_q;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"mux") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_mux;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"muy") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_muy;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"muz") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_muz;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"mu") == 0) {
if (!atom->mu_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_mu;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"radius") == 0) {
if (!atom->radius_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_radius;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"diameter") == 0) {
if (!atom->radius_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_diameter;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"omegax") == 0) {
if (!atom->omega_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_omegax;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"omegay") == 0) {
if (!atom->omega_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_omegay;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"omegaz") == 0) {
if (!atom->omega_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_omegaz;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"angmomx") == 0) {
if (!atom->angmom_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_angmomx;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"angmomy") == 0) {
if (!atom->angmom_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_angmomy;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"angmomz") == 0) {
if (!atom->angmom_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_angmomz;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"tqx") == 0) {
if (!atom->torque_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_tqx;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"tqy") == 0) {
if (!atom->torque_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_tqy;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
} else if (strcmp(arg[iarg],"tqz") == 0) {
if (!atom->torque_flag)
error->all(FLERR,"Dumping an atom property that isn't allocated");
pack_choice[i] = &DumpCustom::pack_tqz;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
// compute value = c_ID
// if no trailing [], then arg is set to 0, else arg is int between []
} else if (strncmp(arg[iarg],"c_",2) == 0) {
pack_choice[i] = &DumpCustom::pack_compute;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1382,7 +1381,7 @@ int DumpCustom::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"f_",2) == 0) {
pack_choice[i] = &DumpCustom::pack_fix;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1415,7 +1414,7 @@ int DumpCustom::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"v_",2) == 0) {
pack_choice[i] = &DumpCustom::pack_variable;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1435,7 +1434,7 @@ int DumpCustom::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"d_",2) == 0) {
pack_choice[i] = &DumpCustom::pack_custom;
- vtype[i] = DOUBLE;
+ vtype[i] = Dump::DOUBLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1457,7 +1456,7 @@ int DumpCustom::parse_fields(int narg, char **arg)
} else if (strncmp(arg[iarg],"i_",2) == 0) {
pack_choice[i] = &DumpCustom::pack_custom;
- vtype[i] = INT;
+ vtype[i] = Dump::INT;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
@@ -1676,7 +1675,7 @@ int DumpCustom::modify_param(int narg, char **arg)
if (strcmp(arg[0],"refresh") == 0) {
if (narg < 2) error->all(FLERR,"Illegal dump_modify command");
- if (strncmp(arg[1],"c_",2) != 0)
+ if (strncmp(arg[1],"c_",2) != 0)
error->all(FLERR,"Illegal dump_modify command");
if (refreshflag) error->all(FLERR,"Dump modify can only have one refresh");
diff --git a/src/error.cpp b/src/error.cpp
index 3c8d8bb34cda52a888342a6165f7a19424f91991..2ab816c992ac910fd0a2dd789a1c08d3f9da2576 100644
--- a/src/error.cpp
+++ b/src/error.cpp
@@ -21,6 +21,20 @@
using namespace LAMMPS_NS;
+// helper function to truncate a string to a segment starting with "src/";
+
+static const char *truncpath(const char *path)
+{
+ if (path) {
+ int len = strlen(path);
+ for (int i = len-4; i > 0; --i) {
+ if (strncmp("src/",path+i,4) == 0)
+ return path+i;
+ }
+ }
+ return path;
+}
+
/* ---------------------------------------------------------------------- */
Error::Error(LAMMPS *lmp) : Pointers(lmp) {
@@ -42,9 +56,9 @@ void Error::universe_all(const char *file, int line, const char *str)
if (universe->me == 0) {
if (universe->uscreen) fprintf(universe->uscreen,
- "ERROR: %s (%s:%d)\n",str,file,line);
+ "ERROR: %s (%s:%d)\n",str,truncpath(file),line);
if (universe->ulogfile) fprintf(universe->ulogfile,
- "ERROR: %s (%s:%d)\n",str,file,line);
+ "ERROR: %s (%s:%d)\n",str,truncpath(file),line);
}
if (output) delete output;
@@ -73,7 +87,7 @@ void Error::universe_one(const char *file, int line, const char *str)
{
if (universe->uscreen)
fprintf(universe->uscreen,"ERROR on proc %d: %s (%s:%d)\n",
- universe->me,str,file,line);
+ universe->me,str,truncpath(file),line);
#ifdef LAMMPS_EXCEPTIONS
char msg[100];
@@ -93,7 +107,7 @@ void Error::universe_warn(const char *file, int line, const char *str)
{
if (universe->uscreen)
fprintf(universe->uscreen,"WARNING on proc %d: %s (%s:%d)\n",
- universe->me,str,file,line);
+ universe->me,str,truncpath(file),line);
}
/* ----------------------------------------------------------------------
@@ -116,10 +130,10 @@ void Error::all(const char *file, int line, const char *str)
if (input && input->line) lastcmd = input->line;
if (screen) fprintf(screen,"ERROR: %s (%s:%d)\n"
"Last command: %s\n",
- str,file,line,lastcmd);
+ str,truncpath(file),line,lastcmd);
if (logfile) fprintf(logfile,"ERROR: %s (%s:%d)\n"
"Last command: %s\n",
- str,file,line,lastcmd);
+ str,truncpath(file),line,lastcmd);
}
#ifdef LAMMPS_EXCEPTIONS
@@ -158,15 +172,15 @@ void Error::one(const char *file, int line, const char *str)
if (input && input->line) lastcmd = input->line;
if (screen) fprintf(screen,"ERROR on proc %d: %s (%s:%d)\n"
"Last command: %s\n",
- me,str,file,line,lastcmd);
+ me,str,truncpath(file),line,lastcmd);
if (logfile) fprintf(logfile,"ERROR on proc %d: %s (%s:%d)\n"
"Last command: %s\n",
- me,str,file,line,lastcmd);
+ me,str,truncpath(file),line,lastcmd);
if (universe->nworlds > 1)
if (universe->uscreen)
fprintf(universe->uscreen,"ERROR on proc %d: %s (%s:%d)\n",
- universe->me,str,file,line);
+ universe->me,str,truncpath(file),line);
#ifdef LAMMPS_EXCEPTIONS
char msg[100];
@@ -184,9 +198,9 @@ void Error::one(const char *file, int line, const char *str)
void Error::warning(const char *file, int line, const char *str, int logflag)
{
- if (screen) fprintf(screen,"WARNING: %s (%s:%d)\n",str,file,line);
+ if (screen) fprintf(screen,"WARNING: %s (%s:%d)\n",str,truncpath(file),line);
if (logflag && logfile) fprintf(logfile,"WARNING: %s (%s:%d)\n",
- str,file,line);
+ str,truncpath(file),line);
}
/* ----------------------------------------------------------------------
@@ -196,8 +210,8 @@ void Error::warning(const char *file, int line, const char *str, int logflag)
void Error::message(const char *file, int line, const char *str, int logflag)
{
- if (screen) fprintf(screen,"%s (%s:%d)\n",str,file,line);
- if (logflag && logfile) fprintf(logfile,"%s (%s:%d)\n",str,file,line);
+ if (screen) fprintf(screen,"%s (%s:%d)\n",str,truncpath(file),line);
+ if (logflag && logfile) fprintf(logfile,"%s (%s:%d)\n",str,truncpath(file),line);
}
/* ----------------------------------------------------------------------
diff --git a/src/fix_balance.cpp b/src/fix_balance.cpp
index 29ea5360661a8779db624233f79594526f1d5a45..b2f545c73fbbefa4af2f2a46d244ea78be352a60 100644
--- a/src/fix_balance.cpp
+++ b/src/fix_balance.cpp
@@ -33,7 +33,6 @@ using namespace LAMMPS_NS;
using namespace FixConst;
enum{SHIFT,BISECTION};
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
/* ---------------------------------------------------------------------- */
@@ -264,10 +263,10 @@ void FixBalance::rebalance()
int *sendproc;
if (lbstyle == SHIFT) {
itercount = balance->shift();
- comm->layout = LAYOUT_NONUNIFORM;
+ comm->layout = Comm::LAYOUT_NONUNIFORM;
} else if (lbstyle == BISECTION) {
sendproc = balance->bisection();
- comm->layout = LAYOUT_TILED;
+ comm->layout = Comm::LAYOUT_TILED;
}
// output of new decomposition
diff --git a/src/fix_deform.cpp b/src/fix_deform.cpp
index 65e1bee63b47310bff3931ff0b7b7679e91a9c1f..c960c310f618fda3816cd658768505161611d2e5 100644
--- a/src/fix_deform.cpp
+++ b/src/fix_deform.cpp
@@ -40,10 +40,6 @@ using namespace MathConst;
enum{NONE=0,FINAL,DELTA,SCALE,VEL,ERATE,TRATE,VOLUME,WIGGLE,VARIABLE};
enum{ONE_FROM_ONE,ONE_FROM_TWO,TWO_FROM_ONE};
-// same as domain.cpp, fix_nvt_sllod.cpp, compute_temp_deform.cpp
-
-enum{NO_REMAP,X_REMAP,V_REMAP};
-
/* ---------------------------------------------------------------------- */
FixDeform::FixDeform(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg),
@@ -208,7 +204,7 @@ rfix(NULL), irregular(NULL), set(NULL)
// no x remap effectively moves atoms within box, so set restart_pbc
options(narg-iarg,&arg[iarg]);
- if (remapflag != X_REMAP) restart_pbc = 1;
+ if (remapflag != Domain::X_REMAP) restart_pbc = 1;
// setup dimflags used by other classes to check for volume-change conflicts
@@ -890,7 +886,7 @@ void FixDeform::end_of_step()
// convert atoms and rigid bodies to lamda coords
- if (remapflag == X_REMAP) {
+ if (remapflag == Domain::X_REMAP) {
double **x = atom->x;
int *mask = atom->mask;
int nlocal = atom->nlocal;
@@ -930,7 +926,7 @@ void FixDeform::end_of_step()
// convert atoms and rigid bodies back to box coords
- if (remapflag == X_REMAP) {
+ if (remapflag == Domain::X_REMAP) {
double **x = atom->x;
int *mask = atom->mask;
int nlocal = atom->nlocal;
@@ -992,7 +988,7 @@ void FixDeform::options(int narg, char **arg)
{
if (narg < 0) error->all(FLERR,"Illegal fix deform command");
- remapflag = X_REMAP;
+ remapflag = Domain::X_REMAP;
scaleflag = 1;
flipflag = 1;
@@ -1000,9 +996,9 @@ void FixDeform::options(int narg, char **arg)
while (iarg < narg) {
if (strcmp(arg[iarg],"remap") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal fix deform command");
- if (strcmp(arg[iarg+1],"x") == 0) remapflag = X_REMAP;
- else if (strcmp(arg[iarg+1],"v") == 0) remapflag = V_REMAP;
- else if (strcmp(arg[iarg+1],"none") == 0) remapflag = NO_REMAP;
+ if (strcmp(arg[iarg+1],"x") == 0) remapflag = Domain::X_REMAP;
+ else if (strcmp(arg[iarg+1],"v") == 0) remapflag = Domain::V_REMAP;
+ else if (strcmp(arg[iarg+1],"none") == 0) remapflag = Domain::NO_REMAP;
else error->all(FLERR,"Illegal fix deform command");
iarg += 2;
} else if (strcmp(arg[iarg],"units") == 0) {
diff --git a/src/fix_nvt_sllod.cpp b/src/fix_nvt_sllod.cpp
index 270d721a008552e1de3ac0096df02cbfd7b17529..012d67c3187c82b795ea6eaf60e4a8abed5e678a 100644
--- a/src/fix_nvt_sllod.cpp
+++ b/src/fix_nvt_sllod.cpp
@@ -31,8 +31,6 @@
using namespace LAMMPS_NS;
using namespace FixConst;
-enum{NO_REMAP,X_REMAP,V_REMAP}; // same as fix_deform.cpp
-
/* ---------------------------------------------------------------------- */
FixNVTSllod::FixNVTSllod(LAMMPS *lmp, int narg, char **arg) :
@@ -83,7 +81,7 @@ void FixNVTSllod::init()
int i;
for (i = 0; i < modify->nfix; i++)
if (strcmp(modify->fix[i]->style,"deform") == 0) {
- if (((FixDeform *) modify->fix[i])->remapflag != V_REMAP)
+ if (((FixDeform *) modify->fix[i])->remapflag != Domain::V_REMAP)
error->all(FLERR,"Using fix nvt/sllod with inconsistent fix deform "
"remap option");
break;
diff --git a/src/irregular.cpp b/src/irregular.cpp
index 451ce7af668b9550bcfcdf3dce670b90d9cf997b..9c15f135d0e4bfc13e40e799f69499a754af68bf 100644
--- a/src/irregular.cpp
+++ b/src/irregular.cpp
@@ -34,8 +34,6 @@ static int compare_standalone(const void *, const void *);
static int compare_standalone(const int, const int, void *);
#endif
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
#define BUFFACTOR 1.5
#define BUFMIN 1000
#define BUFEXTRA 1000
@@ -222,7 +220,7 @@ int Irregular::migrate_check()
// migrate required if comm layout is tiled
// cannot use myloc[] logic below
- if (comm->layout == LAYOUT_TILED) return 1;
+ if (comm->layout == Comm::LAYOUT_TILED) return 1;
// subbox bounds for orthogonal or triclinic box
diff --git a/src/nbin_standard.cpp b/src/nbin_standard.cpp
index e6941014f9f60e441bc99831656465f2bb2d7a0e..9a28121384bcfd15bf5c3fe13b8b406e6bc1fa18 100644
--- a/src/nbin_standard.cpp
+++ b/src/nbin_standard.cpp
@@ -22,8 +22,6 @@
using namespace LAMMPS_NS;
-enum{NSQ,BIN,MULTI}; // also in Neighbor
-
#define SMALL 1.0e-6
#define CUT2BIN_RATIO 100
@@ -94,7 +92,7 @@ void NBinStandard::setup_bins(int style)
double binsize_optimal;
if (binsizeflag) binsize_optimal = binsize_user;
- else if (style == BIN) binsize_optimal = 0.5*cutneighmax;
+ else if (style == Neighbor::BIN) binsize_optimal = 0.5*cutneighmax;
else binsize_optimal = 0.5*cutneighmin;
if (binsize_optimal == 0.0) binsize_optimal = bbox[0];
double binsizeinv = 1.0/binsize_optimal;
diff --git a/src/neigh_list.cpp b/src/neigh_list.cpp
index 2851bae14fa4fec3d1dfe233465c6d03bb002acb..6bdb9beae62e6924b201af00b3c60bb5086dcc9a 100644
--- a/src/neigh_list.cpp
+++ b/src/neigh_list.cpp
@@ -26,8 +26,6 @@ using namespace LAMMPS_NS;
#define PGDELTA 1
-enum{NSQ,BIN,MULTI}; // also in Neighbor
-
/* ---------------------------------------------------------------------- */
NeighList::NeighList(LAMMPS *lmp) : Pointers(lmp)
diff --git a/src/neighbor.cpp b/src/neighbor.cpp
index d5cb83b37d5fc8176b4aceb7527aa4f11e23ca34..912a636227b84b0fb81a549a7297698bbd98c1b0 100644
--- a/src/neighbor.cpp
+++ b/src/neighbor.cpp
@@ -57,7 +57,6 @@ using namespace NeighConst;
#define BIG 1.0e20
-enum{NSQ,BIN,MULTI}; // also in NBin, NeighList, NStencil
enum{NONE,ALL,PARTIAL,TEMPLATE};
static const char cite_neigh_multi[] =
@@ -85,7 +84,7 @@ pairclass(NULL), pairnames(NULL), pairmasks(NULL)
firsttime = 1;
- style = BIN;
+ style = Neighbor::BIN;
every = 1;
delay = 10;
dist_check = 1;
@@ -651,7 +650,7 @@ int Neighbor::init_pair()
// at time of binning when neighbor lists are rebuilt,
// similar to what vanilla Nbin::coord2atom() does now in atom2bin
- if (style == BIN) {
+ if (style == Neighbor::BIN) {
for (i = 0; i < nrequest; i++)
if (requests[i]->occasional && requests[i]->ghost)
error->all(FLERR,"Cannot request an occasional binned neighbor list "
@@ -1420,7 +1419,7 @@ void Neighbor::print_pairwise_info()
bbox[1] = bboxhi[1]-bboxlo[1];
bbox[2] = bboxhi[2]-bboxlo[2];
if (binsizeflag) binsize = binsize_user;
- else if (style == BIN) binsize = 0.5*cutneighmax;
+ else if (style == Neighbor::BIN) binsize = 0.5*cutneighmax;
else binsize = 0.5*cutneighmin;
if (binsize == 0.0) binsize = bbox[0];
@@ -1445,7 +1444,7 @@ void Neighbor::print_pairwise_info()
oneatom, pgsize);
fprintf(out," master list distance cutoff = %g\n",cutneighmax);
fprintf(out," ghost atom cutoff = %g\n",cutghost);
- if (style != NSQ)
+ if (style != Neighbor::NSQ)
fprintf(out," binsize = %g, bins = %g %g %g\n",binsize,
ceil(bbox[0]/binsize), ceil(bbox[1]/binsize),
ceil(bbox[2]/binsize));
@@ -1595,7 +1594,7 @@ int Neighbor::choose_bin(NeighRequest *rq)
{
// no binning needed
- if (style == NSQ) return 0;
+ if (style == Neighbor::NSQ) return 0;
if (rq->skip || rq->copy || rq->halffull) return 0;
// use request settings to match exactly one NBin class mask
@@ -1635,7 +1634,7 @@ int Neighbor::choose_stencil(NeighRequest *rq)
{
// no stencil creation needed
- if (style == NSQ) return 0;
+ if (style == Neighbor::NSQ) return 0;
if (rq->skip || rq->copy || rq->halffull) return 0;
// convert newton request to newtflag = on or off
@@ -1686,9 +1685,9 @@ int Neighbor::choose_stencil(NeighRequest *rq)
// neighbor style is BIN or MULTI and must match
- if (style == BIN) {
+ if (style == Neighbor::BIN) {
if (!(mask & NS_BIN)) continue;
- } else if (style == MULTI) {
+ } else if (style == Neighbor::MULTI) {
if (!(mask & NS_MULTI)) continue;
}
@@ -1817,11 +1816,11 @@ int Neighbor::choose_pair(NeighRequest *rq)
// neighbor style is one of NSQ,BIN,MULTI and must match
- if (style == NSQ) {
+ if (style == Neighbor::NSQ) {
if (!(mask & NP_NSQ)) continue;
- } else if (style == BIN) {
+ } else if (style == Neighbor::BIN) {
if (!(mask & NP_BIN)) continue;
- } else if (style == MULTI) {
+ } else if (style == Neighbor::MULTI) {
if (!(mask & NP_MULTI)) continue;
}
@@ -2063,7 +2062,7 @@ void Neighbor::build(int topoflag)
// if bin then, atoms may have moved outside of proc domain & bin extent,
// leading to errors or even a crash
- if (style != NSQ) {
+ if (style != Neighbor::NSQ) {
for (int i = 0; i < nbin; i++) {
neigh_bin[i]->bin_atoms_setup(nall);
neigh_bin[i]->bin_atoms();
@@ -2180,12 +2179,12 @@ void Neighbor::set(int narg, char **arg)
skin = force->numeric(FLERR,arg[0]);
if (skin < 0.0) error->all(FLERR,"Illegal neighbor command");
- if (strcmp(arg[1],"nsq") == 0) style = NSQ;
- else if (strcmp(arg[1],"bin") == 0) style = BIN;
- else if (strcmp(arg[1],"multi") == 0) style = MULTI;
+ if (strcmp(arg[1],"nsq") == 0) style = Neighbor::NSQ;
+ else if (strcmp(arg[1],"bin") == 0) style = Neighbor::BIN;
+ else if (strcmp(arg[1],"multi") == 0) style = Neighbor::MULTI;
else error->all(FLERR,"Illegal neighbor command");
- if (style == MULTI && lmp->citeme) lmp->citeme->add(cite_neigh_multi);
+ if (style == Neighbor::MULTI && lmp->citeme) lmp->citeme->add(cite_neigh_multi);
}
/* ----------------------------------------------------------------------
diff --git a/src/neighbor.h b/src/neighbor.h
index 12f90556c2512013e99663a9d7f92cebd7f5f5da..751beeae4baa47d40cdcf0e7c852e0d056207cd7 100644
--- a/src/neighbor.h
+++ b/src/neighbor.h
@@ -21,6 +21,7 @@ namespace LAMMPS_NS {
class Neighbor : protected Pointers {
public:
+ enum{NSQ,BIN,MULTI};
int style; // 0,1,2 = nsq, bin, multi
int every; // build every this many steps
int delay; // delay build for this many steps
diff --git a/src/nstencil.cpp b/src/nstencil.cpp
index 18e815d0c998c37f814441ee4333dfbddc5a8b8d..431906e89856a9fe7f03767b0635dc11d5f5566c 100644
--- a/src/nstencil.cpp
+++ b/src/nstencil.cpp
@@ -22,8 +22,6 @@
using namespace LAMMPS_NS;
-enum{NSQ,BIN,MULTI}; // also in Neighbor
-
/* ----------------------------------------------------------------------
NStencil classes
each has method to create a stencil = list of bin offsets
@@ -161,7 +159,7 @@ void NStencil::create_setup()
// reallocate stencil structs if necessary
// for BIN and MULTI styles
- if (neighstyle == BIN) {
+ if (neighstyle == Neighbor::BIN) {
if (smax > maxstencil) {
maxstencil = smax;
memory->destroy(stencil);
@@ -227,10 +225,10 @@ double NStencil::bin_distance(int i, int j, int k)
bigint NStencil::memory_usage()
{
bigint bytes = 0;
- if (neighstyle == BIN) {
+ if (neighstyle == Neighbor::BIN) {
bytes += memory->usage(stencil,maxstencil);
bytes += memory->usage(stencilxyz,maxstencil,3);
- } else if (neighstyle == MULTI) {
+ } else if (neighstyle == Neighbor::MULTI) {
bytes += atom->ntypes*maxstencil_multi * sizeof(int);
bytes += atom->ntypes*maxstencil_multi * sizeof(double);
}
diff --git a/src/ntopo.h b/src/ntopo.h
index 9512606ba49a3ebd45ca87fbab91f2513f625f22..b115b329654c7468d7d41057d03c3a83fb89912b 100644
--- a/src/ntopo.h
+++ b/src/ntopo.h
@@ -31,8 +31,6 @@ class NTopo : protected Pointers {
bigint memory_usage();
protected:
- enum{IGNORE,WARN,ERROR}; // same as thermo.cpp
-
int me,nprocs;
int maxbond,maxangle,maxdihedral,maximproper;
int cluster_check; // copy from Neighbor
diff --git a/src/ntopo_angle_all.cpp b/src/ntopo_angle_all.cpp
index 3a079ab4677533c46902cf337702c7dec60d74b5..2a358c8ce574e4291ec8a492e7904402d7622a90 100644
--- a/src/ntopo_angle_all.cpp
+++ b/src/ntopo_angle_all.cpp
@@ -58,7 +58,7 @@ void NTopoAngleAll::build()
atom3 = atom->map(angle_atom3[i][m]);
if (atom1 == -1 || atom2 == -1 || atom3 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Angle atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT
@@ -86,7 +86,7 @@ void NTopoAngleAll::build()
}
if (cluster_check) angle_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_angle_partial.cpp b/src/ntopo_angle_partial.cpp
index f1d668a3ba2503e0c3992078b518a3239d3d24b6..c82110cda5309f01f3a223c1c4eabc6d3e615051 100644
--- a/src/ntopo_angle_partial.cpp
+++ b/src/ntopo_angle_partial.cpp
@@ -59,7 +59,7 @@ void NTopoAnglePartial::build()
atom3 = atom->map(angle_atom3[i][m]);
if (atom1 == -1 || atom2 == -1 || atom3 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Angle atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT
@@ -87,7 +87,7 @@ void NTopoAnglePartial::build()
}
if (cluster_check) angle_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_angle_template.cpp b/src/ntopo_angle_template.cpp
index 05d5de28a473389bb7c791ce56cf0095b70a2a63..15a8b658f328357e0a27368a22d888bc6b9f9387 100644
--- a/src/ntopo_angle_template.cpp
+++ b/src/ntopo_angle_template.cpp
@@ -76,7 +76,7 @@ void NTopoAngleTemplate::build()
atom3 = atom->map(angle_atom3[iatom][m]+tagprev);
if (atom1 == -1 || atom2 == -1 || atom3 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Angle atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT
@@ -106,7 +106,7 @@ void NTopoAngleTemplate::build()
}
if (cluster_check) angle_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_bond_all.cpp b/src/ntopo_bond_all.cpp
index 03cb2ad86bc83027190c34386adbc5a1d062b503..42e9e2303d7af91b2bf3086f14defc3465bf48b0 100644
--- a/src/ntopo_bond_all.cpp
+++ b/src/ntopo_bond_all.cpp
@@ -55,7 +55,7 @@ void NTopoBondAll::build()
atom1 = atom->map(bond_atom[i][m]);
if (atom1 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Bond atoms " TAGINT_FORMAT " " TAGINT_FORMAT
" missing on proc %d at step " BIGINT_FORMAT,
@@ -78,7 +78,7 @@ void NTopoBondAll::build()
}
if (cluster_check) bond_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_bond_partial.cpp b/src/ntopo_bond_partial.cpp
index cda4bdcf09cbd0150dacd076392e8e5630839c37..5a1acd9191b7a1485d435382fcc1c909622d336d 100644
--- a/src/ntopo_bond_partial.cpp
+++ b/src/ntopo_bond_partial.cpp
@@ -56,7 +56,7 @@ void NTopoBondPartial::build()
atom1 = atom->map(bond_atom[i][m]);
if (atom1 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Bond atoms " TAGINT_FORMAT " " TAGINT_FORMAT
" missing on proc %d at step " BIGINT_FORMAT,
@@ -79,7 +79,7 @@ void NTopoBondPartial::build()
}
if (cluster_check) bond_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_bond_template.cpp b/src/ntopo_bond_template.cpp
index de16d78585cf65175babe696dafb14aa4bd5a2d6..fa98658ad0d640023e80149bdd2b878906120930 100644
--- a/src/ntopo_bond_template.cpp
+++ b/src/ntopo_bond_template.cpp
@@ -72,7 +72,7 @@ void NTopoBondTemplate::build()
atom1 = atom->map(bond_atom[iatom][m]+tagprev);
if (atom1 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Bond atoms " TAGINT_FORMAT " " TAGINT_FORMAT
" missing on proc %d at step " BIGINT_FORMAT,
@@ -96,7 +96,7 @@ void NTopoBondTemplate::build()
}
if (cluster_check) bond_check();
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_dihedral_all.cpp b/src/ntopo_dihedral_all.cpp
index 7a5b350fa018096a7bd46cea5226003d5248ca83..9c94fb10f959e93a9f56db10474039f62feb39b1 100644
--- a/src/ntopo_dihedral_all.cpp
+++ b/src/ntopo_dihedral_all.cpp
@@ -60,7 +60,7 @@ void NTopoDihedralAll::build()
atom4 = atom->map(dihedral_atom4[i][m]);
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Dihedral atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " "
@@ -93,7 +93,7 @@ void NTopoDihedralAll::build()
}
if (cluster_check) dihedral_check(ndihedrallist,dihedrallist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_dihedral_partial.cpp b/src/ntopo_dihedral_partial.cpp
index 603c81e68c674487f48fdd241c93516e9cc47439..14749e6511d74d26fc09ee55cedfab9b272c0443 100644
--- a/src/ntopo_dihedral_partial.cpp
+++ b/src/ntopo_dihedral_partial.cpp
@@ -62,7 +62,7 @@ void NTopoDihedralPartial::build()
atom4 = atom->map(dihedral_atom4[i][m]);
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Dihedral atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " "
@@ -95,7 +95,7 @@ void NTopoDihedralPartial::build()
}
if (cluster_check) dihedral_check(ndihedrallist,dihedrallist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_dihedral_template.cpp b/src/ntopo_dihedral_template.cpp
index 38d319fb9643b6ed5c18cb185f2eea22fbffed5a..8ea860c2e21fdb9b5ee4bc548866b13c2d91a6c0 100644
--- a/src/ntopo_dihedral_template.cpp
+++ b/src/ntopo_dihedral_template.cpp
@@ -78,7 +78,7 @@ void NTopoDihedralTemplate::build()
atom4 = atom->map(dihedral_atom4[iatom][m]+tagprev);
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Dihedral atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " "
@@ -114,7 +114,7 @@ void NTopoDihedralTemplate::build()
}
if (cluster_check) dihedral_check(ndihedrallist,dihedrallist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_improper_all.cpp b/src/ntopo_improper_all.cpp
index ada3927e79c7f70cd2188b46bd40886412e74035..6c478dec527c8e3d75e6c1bdfe0a3dfe5c7396a6 100644
--- a/src/ntopo_improper_all.cpp
+++ b/src/ntopo_improper_all.cpp
@@ -60,7 +60,7 @@ void NTopoImproperAll::build()
atom4 = atom->map(improper_atom4[i][m]);
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Improper atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " "
@@ -93,7 +93,7 @@ void NTopoImproperAll::build()
}
if (cluster_check) dihedral_check(nimproperlist,improperlist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_improper_partial.cpp b/src/ntopo_improper_partial.cpp
index 072a2aa76790d0f762e070c8034dc78e0b73607c..2c37668ca850eb80bcb661ad7857c89bdc185162 100644
--- a/src/ntopo_improper_partial.cpp
+++ b/src/ntopo_improper_partial.cpp
@@ -62,7 +62,7 @@ void NTopoImproperPartial::build()
atom4 = atom->map(improper_atom4[i][m]);
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Improper atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " "
@@ -95,7 +95,7 @@ void NTopoImproperPartial::build()
}
if (cluster_check) dihedral_check(nimproperlist,improperlist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/ntopo_improper_template.cpp b/src/ntopo_improper_template.cpp
index cb15a077d5a30e928cbd35fc20fdac0ce1322a84..953e010d89a0bd5f1e2f862349216edfdaf97a91 100644
--- a/src/ntopo_improper_template.cpp
+++ b/src/ntopo_improper_template.cpp
@@ -78,7 +78,7 @@ void NTopoImproperTemplate::build()
atom4 = atom->map(improper_atom4[iatom][m]+tagprev);
if (atom1 == -1 || atom2 == -1 || atom3 == -1 || atom4 == -1) {
nmissing++;
- if (lostbond == ERROR) {
+ if (lostbond == Thermo::ERROR) {
char str[128];
sprintf(str,"Improper atoms "
TAGINT_FORMAT " " TAGINT_FORMAT " "
@@ -114,7 +114,7 @@ void NTopoImproperTemplate::build()
}
if (cluster_check) dihedral_check(nimproperlist,improperlist);
- if (lostbond == IGNORE) return;
+ if (lostbond == Thermo::IGNORE) return;
int all;
MPI_Allreduce(&nmissing,&all,1,MPI_INT,MPI_SUM,world);
diff --git a/src/pair.h b/src/pair.h
index 3ce567876cf38cd97093abcd26cb0baea3abd758..26488103a05e2044d537f31dcf8a4b8b2fa0073b 100644
--- a/src/pair.h
+++ b/src/pair.h
@@ -96,6 +96,8 @@ class Pair : protected Pointers {
// public so external driver can check
int compute_flag; // 0 if skip compute()
+ enum{GEOMETRIC,ARITHMETIC,SIXTHPOWER}; // mixing options
+
// KOKKOS host/device flag and data masks
ExecutionSpace execution_space;
@@ -191,8 +193,6 @@ class Pair : protected Pointers {
protected:
int instance_me; // which Pair class instantiation I am
- enum{GEOMETRIC,ARITHMETIC,SIXTHPOWER}; // mixing options
-
int special_lj[4]; // copied from force->special_lj for Kokkos
int suffix_flag; // suffix compatibility flag
diff --git a/src/procmap.cpp b/src/procmap.cpp
index f54ac1915c4b4c934d18d3ae0e848c57df856808..9d1ed83e7382a75dd418b311f6ad9bc10c0ca8a0 100644
--- a/src/procmap.cpp
+++ b/src/procmap.cpp
@@ -17,6 +17,7 @@
#include "procmap.h"
#include "universe.h"
+#include "comm.h"
#include "domain.h"
#include "math_extra.h"
#include "memory.h"
@@ -29,8 +30,6 @@ using namespace LAMMPS_NS;
#define MAXLINE 128
-enum{MULTIPLE}; // same as in Comm
-
/* ---------------------------------------------------------------------- */
ProcMap::ProcMap(LAMMPS *lmp) : Pointers(lmp) {}
@@ -811,7 +810,7 @@ int ProcMap::cull_other(int n, int **factors, int m,
{
int i = 0;
while (i < n) {
- if (other_style == MULTIPLE) {
+ if (other_style == Comm::MULTIPLE) {
int flag = 0;
if ((other_procgrid[0]/other_coregrid[0]) % factors[i][0]) flag = 1;
if ((other_procgrid[1]/other_coregrid[1]) % factors[i][1]) flag = 1;
diff --git a/src/replicate.cpp b/src/replicate.cpp
index 7954e4291416b6ea666a18dc639d442093efe610..cdadf1fd1f6e1636f02ef227d2a4a912153f81ba 100644
--- a/src/replicate.cpp
+++ b/src/replicate.cpp
@@ -30,8 +30,6 @@ using namespace LAMMPS_NS;
#define LB_FACTOR 1.1
#define EPSILON 1.0e-6
-enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED}; // several files
-
/* ---------------------------------------------------------------------- */
Replicate::Replicate(LAMMPS *lmp) : Pointers(lmp) {}
@@ -267,7 +265,7 @@ void Replicate::command(int narg, char **arg)
sublo[2] = domain->sublo_lamda[2]; subhi[2] = domain->subhi_lamda[2];
}
- if (comm->layout != LAYOUT_TILED) {
+ if (comm->layout != Comm::LAYOUT_TILED) {
if (domain->xperiodic) {
if (comm->myloc[0] == 0) sublo[0] -= epsilon[0];
if (comm->myloc[0] == comm->procgrid[0]-1) subhi[0] += epsilon[0];
diff --git a/src/thermo.cpp b/src/thermo.cpp
index fb3d1abc485df6cd975a68835dd99d771247b287..ade7a3c333ab893fdd14a17adaa94a966328b64a 100644
--- a/src/thermo.cpp
+++ b/src/thermo.cpp
@@ -72,7 +72,6 @@ using namespace MathConst;
#define ONE "step temp epair emol etotal press"
#define MULTI "etotal ke temp pe ebond eangle edihed eimp evdwl ecoul elong press"
-enum{IGNORE,WARN,ERROR}; // same as several files
enum{ONELINE,MULTILINE};
enum{INT,FLOAT,BIGINT};
enum{SCALAR,VECTOR,ARRAY};
@@ -98,7 +97,7 @@ Thermo::Thermo(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp)
modified = 0;
normuserflag = 0;
lineflag = ONELINE;
- lostflag = lostbond = ERROR;
+ lostflag = lostbond = Thermo::ERROR;
lostbefore = 0;
flushflag = 0;
@@ -427,14 +426,14 @@ bigint Thermo::lost_check()
if (ntotal == atom->natoms) return ntotal;
// if not checking or already warned, just return
- if (lostflag == IGNORE) return ntotal;
- if (lostflag == WARN && lostbefore == 1) {
+ if (lostflag == Thermo::IGNORE) return ntotal;
+ if (lostflag == Thermo::WARN && lostbefore == 1) {
return ntotal;
}
// error message
- if (lostflag == ERROR) {
+ if (lostflag == Thermo::ERROR) {
char str[64];
sprintf(str,
"Lost atoms: original " BIGINT_FORMAT " current " BIGINT_FORMAT,
@@ -536,17 +535,17 @@ void Thermo::modify_params(int narg, char **arg)
} else if (strcmp(arg[iarg],"lost") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal thermo_modify command");
- if (strcmp(arg[iarg+1],"ignore") == 0) lostflag = IGNORE;
- else if (strcmp(arg[iarg+1],"warn") == 0) lostflag = WARN;
- else if (strcmp(arg[iarg+1],"error") == 0) lostflag = ERROR;
+ if (strcmp(arg[iarg+1],"ignore") == 0) lostflag = Thermo::IGNORE;
+ else if (strcmp(arg[iarg+1],"warn") == 0) lostflag = Thermo::WARN;
+ else if (strcmp(arg[iarg+1],"error") == 0) lostflag = Thermo::ERROR;
else error->all(FLERR,"Illegal thermo_modify command");
iarg += 2;
} else if (strcmp(arg[iarg],"lost/bond") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal thermo_modify command");
- if (strcmp(arg[iarg+1],"ignore") == 0) lostbond = IGNORE;
- else if (strcmp(arg[iarg+1],"warn") == 0) lostbond = WARN;
- else if (strcmp(arg[iarg+1],"error") == 0) lostbond = ERROR;
+ if (strcmp(arg[iarg+1],"ignore") == 0) lostbond = Thermo::IGNORE;
+ else if (strcmp(arg[iarg+1],"warn") == 0) lostbond = Thermo::WARN;
+ else if (strcmp(arg[iarg+1],"error") == 0) lostbond = Thermo::ERROR;
else error->all(FLERR,"Illegal thermo_modify command");
iarg += 2;
diff --git a/src/thermo.h b/src/thermo.h
index 8023a8867cf337b17b73059775d7e3314b6559d6..8c32f24d3c2142f4953126e97dbd6bdd38921928 100644
--- a/src/thermo.h
+++ b/src/thermo.h
@@ -33,6 +33,8 @@ class Thermo : protected Pointers {
int lostflag; // IGNORE,WARN,ERROR
int lostbond; // ditto for atoms in bonds
+ enum {IGNORE,WARN,ERROR};
+
Thermo(class LAMMPS *, int, char **);
~Thermo();
void init();
diff --git a/src/version.h b/src/version.h
index eda788fad6773bb16161275945120f9ce15f5383..981ec8eeb8c8e47627ab4d756be509e7b7fe2872 100644
--- a/src/version.h
+++ b/src/version.h
@@ -1 +1 @@
-#define LAMMPS_VERSION "20 Apr 2018"
+#define LAMMPS_VERSION "11 May 2018"
diff --git a/src/write_data.cpp b/src/write_data.cpp
index 85fbe6e5c7cd3f4fa309f8ee470feb1ae87214cc..96bf081157eaf13096dd215477950a5ee0792201 100644
--- a/src/write_data.cpp
+++ b/src/write_data.cpp
@@ -36,7 +36,6 @@
using namespace LAMMPS_NS;
-enum{IGNORE,WARN,ERROR}; // same as thermo.cpp
enum{II,IJ};
/* ---------------------------------------------------------------------- */
@@ -153,7 +152,7 @@ void WriteData::write(char *file)
bigint nblocal = atom->nlocal;
bigint natoms;
MPI_Allreduce(&nblocal,&natoms,1,MPI_LMP_BIGINT,MPI_SUM,world);
- if (natoms != atom->natoms && output->thermo->lostflag == ERROR)
+ if (natoms != atom->natoms && output->thermo->lostflag == Thermo::ERROR)
error->all(FLERR,"Atom count is inconsistent, cannot write data file");
// sum up bond,angle,dihedral,improper counts
diff --git a/src/write_restart.cpp b/src/write_restart.cpp
index c82791e2c1c7c1d0f8a50e6faa695c39873e8c8a..69b731870d0c28fe6f389da0183f1905664a0400 100644
--- a/src/write_restart.cpp
+++ b/src/write_restart.cpp
@@ -63,8 +63,6 @@ enum{VERSION,SMALLINT,TAGINT,BIGINT,
ATOM_ID,ATOM_MAP_STYLE,ATOM_MAP_USER,ATOM_SORTFREQ,ATOM_SORTBIN,
COMM_MODE,COMM_CUTOFF,COMM_VEL,NO_PAIR};
-enum{IGNORE,WARN,ERROR}; // same as thermo.cpp
-
/* ---------------------------------------------------------------------- */
WriteRestart::WriteRestart(LAMMPS *lmp) : Pointers(lmp)
@@ -252,7 +250,7 @@ void WriteRestart::write(char *file)
bigint nblocal = atom->nlocal;
MPI_Allreduce(&nblocal,&natoms,1,MPI_LMP_BIGINT,MPI_SUM,world);
- if (natoms != atom->natoms && output->thermo->lostflag == ERROR)
+ if (natoms != atom->natoms && output->thermo->lostflag == Thermo::ERROR)
error->all(FLERR,"Atom count is inconsistent, cannot write restart file");
// open single restart file or base file for multiproc case