jianmu
/
lammps
forked from lijiext/lammps
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge pull request #926 from ndtrung81/gpu-updates 

some updates to the GPU package
This commit is contained in:
Steve Plimpton 2018-06-13 12:57:01 -06:00 committed by GitHub
parent 29e555213c 3c781afa6c
commit 661848139c
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
23 changed files with 4757 additions and 75 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
cmake/CMakeLists.txt
View File

@ -760,7 +760,7 @@ if(PKG_GPU)
add_library(gpu STATIC ${GPU_LIB_SOURCES})
target_link_libraries(gpu ${OpenCL_LIBRARIES})
target_include_directories(gpu PRIVATE ${CMAKE_CURRENT_BINARY_DIR}/gpu ${OpenCL_INCLUDE_DIRS})
target_compile_definitions(gpu PRIVATE -D_${GPU_PREC} -DMPI_GERYON -DUCL_NO_EXIT)
target_compile_definitions(gpu PRIVATE -D_${GPU_PREC} -D${OCL_TUNE}_OCL -DMPI_GERYON -DUCL_NO_EXIT)
target_compile_definitions(gpu PRIVATE -DUSE_OPENCL)
list(APPEND LAMMPS_LINK_LIBS gpu)

53
lib/gpu/Makefile.linux_multi Normal file
View File

@ -0,0 +1,53 @@
# /* ----------------------------------------------------------------------
# Generic Linux Makefile for CUDA
# - Change CUDA_ARCH for your GPU
# ------------------------------------------------------------------------- */
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.standard
ifeq ($(CUDA_HOME),)
CUDA_HOME = /usr/local/cuda
endif
NVCC = nvcc
# Kepler CUDA
#CUDA_ARCH = -arch=sm_35
# newer CUDA
#CUDA_ARCH = -arch=sm_13
# older CUDA
#CUDA_ARCH = -arch=sm_10 -DCUDA_PRE_THREE
CUDA_ARCH = -gencode arch=compute_60,code=sm_60 -gencode arch=compute_61,code=sm_61
# this setting should match LAMMPS Makefile
# one of LAMMPS_SMALLBIG (default), LAMMPS_BIGBIG and LAMMPS_SMALLSMALL
LMP_INC = -DLAMMPS_SMALLBIG
# precision for GPU calculations
# -D_SINGLE_SINGLE # Single precision for all calculations
# -D_DOUBLE_DOUBLE # Double precision for all calculations
# -D_SINGLE_DOUBLE # Accumulation of forces, etc. in double
CUDA_PRECISION = -D_SINGLE_DOUBLE
CUDA_INCLUDE = -I$(CUDA_HOME)/include
CUDA_LIB = -L$(CUDA_HOME)/lib64
CUDA_OPTS = -DUNIX -O3 -Xptxas -v --use_fast_math $(LMP_INC) -Xcompiler "-fPIC -std=c++98"
CUDR_CPP = mpicxx -DMPI_GERYON -DUCL_NO_EXIT -DMPICH_IGNORE_CXX_SEEK -DOMPI_SKIP_MPICXX=1 -fPIC
CUDR_OPTS = -O2 $(LMP_INC) # -xHost -no-prec-div -ansi-alias
BIN_DIR = ./
OBJ_DIR = ./
LIB_DIR = ./
AR = ar
BSH = /bin/sh
CUDPP_OPT = -DUSE_CUDPP -Icudpp_mini
include Nvidia.makefile_multi

32
lib/gpu/Nvidia.makefile
View File

@ -77,7 +77,9 @@ OBJS = $(OBJ_DIR)/lal_atom.o $(OBJ_DIR)/lal_ans.o \
$(OBJ_DIR)/lal_coul_debye.o $(OBJ_DIR)/lal_coul_debye_ext.o \
$(OBJ_DIR)/lal_zbl.o $(OBJ_DIR)/lal_zbl_ext.o \
$(OBJ_DIR)/lal_lj_cubic.o $(OBJ_DIR)/lal_lj_cubic_ext.o \
$(OBJ_DIR)/lal_ufm.o $(OBJ_DIR)/lal_ufm_ext.o
$(OBJ_DIR)/lal_ufm.o $(OBJ_DIR)/lal_ufm_ext.o \
$(OBJ_DIR)/lal_dipole_long_lj.o $(OBJ_DIR)/lal_dipole_long_lj_ext.o \
$(OBJ_DIR)/lal_lj_expand_coul_long.o $(OBJ_DIR)/lal_lj_expand_coul_long_ext.o
CBNS = $(OBJ_DIR)/device.cubin $(OBJ_DIR)/device_cubin.h \
$(OBJ_DIR)/atom.cubin $(OBJ_DIR)/atom_cubin.h \
@ -133,7 +135,9 @@ CBNS = $(OBJ_DIR)/device.cubin $(OBJ_DIR)/device_cubin.h \
$(OBJ_DIR)/coul_debye.cubin $(OBJ_DIR)/coul_debye_cubin.h \
$(OBJ_DIR)/zbl.cubin $(OBJ_DIR)/zbl_cubin.h \
$(OBJ_DIR)/lj_cubic.cubin $(OBJ_DIR)/lj_cubic_cubin.h \
$(OBJ_DIR)/ufm.cubin $(OBJ_DIR)/ufm_cubin.h
$(OBJ_DIR)/ufm.cubin $(OBJ_DIR)/ufm_cubin.h \
$(OBJ_DIR)/dipole_long_lj.cubin $(OBJ_DIR)/dipole_long_lj_cubin.h \
$(OBJ_DIR)/lj_expand_coul_long.cubin $(OBJ_DIR)/lj_expand_coul_long_cubin.h
all: $(OBJ_DIR) $(GPU_LIB) $(EXECS)
@ -809,6 +813,30 @@ $(OBJ_DIR)/lal_lj_cubic.o: $(ALL_H) lal_lj_cubic.h lal_lj_cubic.cpp $(OBJ_DIR)/l
$(OBJ_DIR)/lal_lj_cubic_ext.o: $(ALL_H) lal_lj_cubic.h lal_lj_cubic_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_lj_cubic_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/dipole_long_lj.cubin: lal_dipole_long_lj.cu lal_precision.h lal_preprocessor.h
$(CUDA) --cubin -DNV_KERNEL -o $@ lal_dipole_long_lj.cu
$(OBJ_DIR)/dipole_long_lj_cubin.h: $(OBJ_DIR)/dipole_long_lj.cubin $(OBJ_DIR)/dipole_long_lj.cubin
$(BIN2C) -c -n dipole_long_lj $(OBJ_DIR)/dipole_long_lj.cubin > $(OBJ_DIR)/dipole_long_lj_cubin.h
$(OBJ_DIR)/lal_dipole_long_lj.o: $(ALL_H) lal_dipole_long_lj.h lal_dipole_long_lj.cpp $(OBJ_DIR)/dipole_long_lj_cubin.h $(OBJ_DIR)/lal_base_dipole.o
$(CUDR) -o $@ -c lal_dipole_long_lj.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_dipole_long_lj_ext.o: $(ALL_H) lal_dipole_long_lj.h lal_dipole_long_lj_ext.cpp lal_base_dipole.h
$(CUDR) -o $@ -c lal_dipole_long_lj_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_expand_coul_long.cubin: lal_lj_expand_coul_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --cubin -DNV_KERNEL -o $@ lal_lj_expand_coul_long.cu
$(OBJ_DIR)/lj_expand_coul_long_cubin.h: $(OBJ_DIR)/lj_expand_coul_long.cubin $(OBJ_DIR)/lj_expand_coul_long.cubin
$(BIN2C) -c -n lj_expand_coul_long $(OBJ_DIR)/lj_expand_coul_long.cubin > $(OBJ_DIR)/lj_expand_coul_long_cubin.h
$(OBJ_DIR)/lal_lj_expand_coul_long.o: $(ALL_H) lal_lj_expand_coul_long.h lal_lj_expand_coul_long.cpp $(OBJ_DIR)/lj_expand_coul_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_expand_coul_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_expand_coul_long_ext.o: $(ALL_H) lal_lj_expand_coul_long.h lal_lj_expand_coul_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_expand_coul_long_ext.cpp -I$(OBJ_DIR)
$(BIN_DIR)/nvc_get_devices: ./geryon/ucl_get_devices.cpp $(NVD_H)
$(CUDR) -o $@ ./geryon/ucl_get_devices.cpp -DUCL_CUDADR $(CUDA_LIB) -lcuda

854
lib/gpu/Nvidia.makefile_multi Normal file
View File

@ -0,0 +1,854 @@
CUDA = $(NVCC) $(CUDA_INCLUDE) $(CUDA_OPTS) -Icudpp_mini $(CUDA_ARCH) \
$(CUDA_PRECISION)
CUDR = $(CUDR_CPP) $(CUDR_OPTS) $(CUDA_PRECISION) $(CUDA_INCLUDE) \
$(CUDPP_OPT)
CUDA_LINK = $(CUDA_LIB) -lcudart
BIN2C = $(CUDA_HOME)/bin/bin2c
GPU_LIB = $(LIB_DIR)/libgpu.a
# Headers for Geryon
UCL_H = $(wildcard ./geryon/ucl*.h)
NVC_H = $(wildcard ./geryon/nvc*.h) $(UCL_H)
NVD_H = $(wildcard ./geryon/nvd*.h) $(UCL_H) lal_preprocessor.h
# Headers for Pair Stuff
PAIR_H = lal_atom.h lal_answer.h lal_neighbor_shared.h \
lal_neighbor.h lal_precision.h lal_device.h \
lal_balance.h lal_pppm.h
ALL_H = $(NVD_H) $(PAIR_H)
EXECS = $(BIN_DIR)/nvc_get_devices
ifdef CUDPP_OPT
CUDPP = $(OBJ_DIR)/cudpp.o $(OBJ_DIR)/cudpp_plan.o \
$(OBJ_DIR)/cudpp_maximal_launch.o $(OBJ_DIR)/cudpp_plan_manager.o \
$(OBJ_DIR)/radixsort_app.cu_o $(OBJ_DIR)/scan_app.cu_o
endif
OBJS = $(OBJ_DIR)/lal_atom.o $(OBJ_DIR)/lal_ans.o \
$(OBJ_DIR)/lal_neighbor.o $(OBJ_DIR)/lal_neighbor_shared.o \
$(OBJ_DIR)/lal_device.o $(OBJ_DIR)/lal_base_atomic.o \
$(OBJ_DIR)/lal_base_charge.o $(OBJ_DIR)/lal_base_ellipsoid.o \
$(OBJ_DIR)/lal_base_dipole.o $(OBJ_DIR)/lal_base_three.o \
$(OBJ_DIR)/lal_base_dpd.o \
$(OBJ_DIR)/lal_pppm.o $(OBJ_DIR)/lal_pppm_ext.o \
$(OBJ_DIR)/lal_gayberne.o $(OBJ_DIR)/lal_gayberne_ext.o \
$(OBJ_DIR)/lal_re_squared.o $(OBJ_DIR)/lal_re_squared_ext.o \
$(OBJ_DIR)/lal_lj.o $(OBJ_DIR)/lal_lj_ext.o \
$(OBJ_DIR)/lal_lj96.o $(OBJ_DIR)/lal_lj96_ext.o \
$(OBJ_DIR)/lal_lj_expand.o $(OBJ_DIR)/lal_lj_expand_ext.o \
$(OBJ_DIR)/lal_lj_coul.o $(OBJ_DIR)/lal_lj_coul_ext.o \
$(OBJ_DIR)/lal_lj_coul_long.o $(OBJ_DIR)/lal_lj_coul_long_ext.o \
$(OBJ_DIR)/lal_lj_dsf.o $(OBJ_DIR)/lal_lj_dsf_ext.o \
$(OBJ_DIR)/lal_lj_class2_long.o $(OBJ_DIR)/lal_lj_class2_long_ext.o \
$(OBJ_DIR)/lal_coul_long.o $(OBJ_DIR)/lal_coul_long_ext.o \
$(OBJ_DIR)/lal_morse.o $(OBJ_DIR)/lal_morse_ext.o \
$(OBJ_DIR)/lal_charmm_long.o $(OBJ_DIR)/lal_charmm_long_ext.o \
$(OBJ_DIR)/lal_lj_sdk.o $(OBJ_DIR)/lal_lj_sdk_ext.o \
$(OBJ_DIR)/lal_lj_sdk_long.o $(OBJ_DIR)/lal_lj_sdk_long_ext.o \
$(OBJ_DIR)/lal_eam.o $(OBJ_DIR)/lal_eam_ext.o \
$(OBJ_DIR)/lal_eam_fs_ext.o $(OBJ_DIR)/lal_eam_alloy_ext.o \
$(OBJ_DIR)/lal_buck.o $(OBJ_DIR)/lal_buck_ext.o \
$(OBJ_DIR)/lal_buck_coul.o $(OBJ_DIR)/lal_buck_coul_ext.o \
$(OBJ_DIR)/lal_buck_coul_long.o $(OBJ_DIR)/lal_buck_coul_long_ext.o \
$(OBJ_DIR)/lal_table.o $(OBJ_DIR)/lal_table_ext.o \
$(OBJ_DIR)/lal_yukawa.o $(OBJ_DIR)/lal_yukawa_ext.o \
$(OBJ_DIR)/lal_born.o $(OBJ_DIR)/lal_born_ext.o \
$(OBJ_DIR)/lal_born_coul_wolf.o $(OBJ_DIR)/lal_born_coul_wolf_ext.o \
$(OBJ_DIR)/lal_born_coul_long.o $(OBJ_DIR)/lal_born_coul_long_ext.o \
$(OBJ_DIR)/lal_dipole_lj.o $(OBJ_DIR)/lal_dipole_lj_ext.o \
$(OBJ_DIR)/lal_dipole_lj_sf.o $(OBJ_DIR)/lal_dipole_lj_sf_ext.o \
$(OBJ_DIR)/lal_colloid.o $(OBJ_DIR)/lal_colloid_ext.o \
$(OBJ_DIR)/lal_gauss.o $(OBJ_DIR)/lal_gauss_ext.o \
$(OBJ_DIR)/lal_yukawa_colloid.o $(OBJ_DIR)/lal_yukawa_colloid_ext.o \
$(OBJ_DIR)/lal_lj_coul_debye.o $(OBJ_DIR)/lal_lj_coul_debye_ext.o \
$(OBJ_DIR)/lal_coul_dsf.o $(OBJ_DIR)/lal_coul_dsf_ext.o \
$(OBJ_DIR)/lal_sw.o $(OBJ_DIR)/lal_sw_ext.o \
$(OBJ_DIR)/lal_vashishta.o $(OBJ_DIR)/lal_vashishta_ext.o \
$(OBJ_DIR)/lal_beck.o $(OBJ_DIR)/lal_beck_ext.o \
$(OBJ_DIR)/lal_mie.o $(OBJ_DIR)/lal_mie_ext.o \
$(OBJ_DIR)/lal_soft.o $(OBJ_DIR)/lal_soft_ext.o \
$(OBJ_DIR)/lal_lj_coul_msm.o $(OBJ_DIR)/lal_lj_coul_msm_ext.o \
$(OBJ_DIR)/lal_lj_gromacs.o $(OBJ_DIR)/lal_lj_gromacs_ext.o \
$(OBJ_DIR)/lal_dpd.o $(OBJ_DIR)/lal_dpd_ext.o \
$(OBJ_DIR)/lal_tersoff.o $(OBJ_DIR)/lal_tersoff_ext.o \
$(OBJ_DIR)/lal_tersoff_zbl.o $(OBJ_DIR)/lal_tersoff_zbl_ext.o \
$(OBJ_DIR)/lal_tersoff_mod.o $(OBJ_DIR)/lal_tersoff_mod_ext.o \
$(OBJ_DIR)/lal_coul.o $(OBJ_DIR)/lal_coul_ext.o \
$(OBJ_DIR)/lal_coul_debye.o $(OBJ_DIR)/lal_coul_debye_ext.o \
$(OBJ_DIR)/lal_zbl.o $(OBJ_DIR)/lal_zbl_ext.o \
$(OBJ_DIR)/lal_lj_cubic.o $(OBJ_DIR)/lal_lj_cubic_ext.o \
$(OBJ_DIR)/lal_ufm.o $(OBJ_DIR)/lal_ufm_ext.o \
$(OBJ_DIR)/lal_dipole_long_lj.o $(OBJ_DIR)/lal_dipole_long_lj_ext.o \
$(OBJ_DIR)/lal_lj_expand_coul_long.o $(OBJ_DIR)/lal_lj_expand_coul_long_ext.o
CBNS = $(OBJ_DIR)/device.cubin $(OBJ_DIR)/device_cubin.h \
$(OBJ_DIR)/atom.cubin $(OBJ_DIR)/atom_cubin.h \
$(OBJ_DIR)/neighbor_cpu.cubin $(OBJ_DIR)/neighbor_cpu_cubin.h \
$(OBJ_DIR)/neighbor_gpu.cubin $(OBJ_DIR)/neighbor_gpu_cubin.h \
$(OBJ_DIR)/pppm_f.cubin $(OBJ_DIR)/pppm_f_cubin.h \
$(OBJ_DIR)/pppm_d.cubin $(OBJ_DIR)/pppm_d_cubin.h \
$(OBJ_DIR)/ellipsoid_nbor.cubin $(OBJ_DIR)/ellipsoid_nbor_cubin.h \
$(OBJ_DIR)/gayberne.cubin $(OBJ_DIR)/gayberne_lj.cubin \
$(OBJ_DIR)/gayberne_cubin.h $(OBJ_DIR)/gayberne_lj_cubin.h \
$(OBJ_DIR)/re_squared.cubin $(OBJ_DIR)/re_squared_lj.cubin \
$(OBJ_DIR)/re_squared_cubin.h $(OBJ_DIR)/re_squared_lj_cubin.h \
$(OBJ_DIR)/lj.cubin $(OBJ_DIR)/lj_cubin.h \
$(OBJ_DIR)/lj96.cubin $(OBJ_DIR)/lj96_cubin.h \
$(OBJ_DIR)/lj_expand.cubin $(OBJ_DIR)/lj_expand_cubin.h \
$(OBJ_DIR)/lj_coul.cubin $(OBJ_DIR)/lj_coul_cubin.h \
$(OBJ_DIR)/lj_coul_long.cubin $(OBJ_DIR)/lj_coul_long_cubin.h \
$(OBJ_DIR)/lj_dsf.cubin $(OBJ_DIR)/lj_dsf_cubin.h \
$(OBJ_DIR)/lj_class2_long.cubin $(OBJ_DIR)/lj_class2_long_cubin.h \
$(OBJ_DIR)/coul_long.cubin $(OBJ_DIR)/coul_long_cubin.h \
$(OBJ_DIR)/morse.cubin $(OBJ_DIR)/morse_cubin.h \
$(OBJ_DIR)/charmm_long.cubin $(OBJ_DIR)/charmm_long_cubin.h \
$(OBJ_DIR)/lj_sdk.cubin $(OBJ_DIR)/lj_sdk_cubin.h \
$(OBJ_DIR)/lj_sdk_long.cubin $(OBJ_DIR)/lj_sdk_long_cubin.h \
$(OBJ_DIR)/eam.cubin $(OBJ_DIR)/eam_cubin.h \
$(OBJ_DIR)/buck.cubin $(OBJ_DIR)/buck_cubin.h \
$(OBJ_DIR)/buck_coul_long.cubin $(OBJ_DIR)/buck_coul_long_cubin.h \
$(OBJ_DIR)/buck_coul.cubin $(OBJ_DIR)/buck_coul_cubin.h \
$(OBJ_DIR)/table.cubin $(OBJ_DIR)/table_cubin.h \
$(OBJ_DIR)/yukawa.cubin $(OBJ_DIR)/yukawa_cubin.h \
$(OBJ_DIR)/born.cubin $(OBJ_DIR)/born_cubin.h \
$(OBJ_DIR)/born_coul_wolf.cubin $(OBJ_DIR)/born_coul_wolf_cubin.h \
$(OBJ_DIR)/born_coul_long.cubin $(OBJ_DIR)/born_coul_long_cubin.h \
$(OBJ_DIR)/dipole_lj.cubin $(OBJ_DIR)/dipole_lj_cubin.h \
$(OBJ_DIR)/dipole_lj_sf.cubin $(OBJ_DIR)/dipole_lj_sf_cubin.h \
$(OBJ_DIR)/colloid.cubin $(OBJ_DIR)/colloid_cubin.h \
$(OBJ_DIR)/gauss.cubin $(OBJ_DIR)/gauss_cubin.h \
$(OBJ_DIR)/yukawa_colloid.cubin $(OBJ_DIR)/yukawa_colloid_cubin.h \
$(OBJ_DIR)/lj_coul_debye.cubin $(OBJ_DIR)/lj_coul_debye_cubin.h \
$(OBJ_DIR)/coul_dsf.cubin $(OBJ_DIR)/coul_dsf_cubin.h \
$(OBJ_DIR)/sw.cubin $(OBJ_DIR)/sw_cubin.h \
$(OBJ_DIR)/vashishta.cubin $(OBJ_DIR)/vashishta_cubin.h \
$(OBJ_DIR)/beck.cubin $(OBJ_DIR)/beck_cubin.h \
$(OBJ_DIR)/mie.cubin $(OBJ_DIR)/mie_cubin.h \
$(OBJ_DIR)/soft.cubin $(OBJ_DIR)/soft_cubin.h \
$(OBJ_DIR)/lj_coul_msm.cubin $(OBJ_DIR)/lj_coul_msm_cubin.h \
$(OBJ_DIR)/lj_gromacs.cubin $(OBJ_DIR)/lj_gromacs_cubin.h \
$(OBJ_DIR)/dpd.cubin $(OBJ_DIR)/dpd_cubin.h \
$(OBJ_DIR)/tersoff.cubin $(OBJ_DIR)/tersoff_cubin.h \
$(OBJ_DIR)/tersoff_zbl.cubin $(OBJ_DIR)/tersoff_zbl_cubin.h \
$(OBJ_DIR)/tersoff_mod.cubin $(OBJ_DIR)/tersoff_mod_cubin.h \
$(OBJ_DIR)/coul.cubin $(OBJ_DIR)/coul_cubin.h \
$(OBJ_DIR)/coul_debye.cubin $(OBJ_DIR)/coul_debye_cubin.h \
$(OBJ_DIR)/zbl.cubin $(OBJ_DIR)/zbl_cubin.h \
$(OBJ_DIR)/lj_cubic.cubin $(OBJ_DIR)/lj_cubic_cubin.h \
$(OBJ_DIR)/ufm.cubin $(OBJ_DIR)/ufm_cubin.h \
$(OBJ_DIR)/dipole_long_lj.cubin $(OBJ_DIR)/dipole_long_lj_cubin.h \
$(OBJ_DIR)/lj_expand_coul_long.cubin $(OBJ_DIR)/lj_expand_coul_long_cubin.h
all: $(OBJ_DIR) $(GPU_LIB) $(EXECS)
$(OBJ_DIR):
mkdir -p $@
$(OBJ_DIR)/cudpp.o: cudpp_mini/cudpp.cpp
$(CUDR) -o $@ -c cudpp_mini/cudpp.cpp -Icudpp_mini
$(OBJ_DIR)/cudpp_plan.o: cudpp_mini/cudpp_plan.cpp
$(CUDR) -o $@ -c cudpp_mini/cudpp_plan.cpp -Icudpp_mini
$(OBJ_DIR)/cudpp_maximal_launch.o: cudpp_mini/cudpp_maximal_launch.cpp
$(CUDR) -o $@ -c cudpp_mini/cudpp_maximal_launch.cpp -Icudpp_mini
$(OBJ_DIR)/cudpp_plan_manager.o: cudpp_mini/cudpp_plan_manager.cpp
$(CUDR) -o $@ -c cudpp_mini/cudpp_plan_manager.cpp -Icudpp_mini
$(OBJ_DIR)/radixsort_app.cu_o: cudpp_mini/radixsort_app.cu
$(CUDA) -o $@ -c cudpp_mini/radixsort_app.cu
$(OBJ_DIR)/scan_app.cu_o: cudpp_mini/scan_app.cu
$(CUDA) -o $@ -c cudpp_mini/scan_app.cu
$(OBJ_DIR)/atom.cubin: lal_atom.cu lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_atom.cu
$(OBJ_DIR)/atom_cubin.h: $(OBJ_DIR)/atom.cubin
$(BIN2C) -c -n atom $(OBJ_DIR)/atom.cubin > $(OBJ_DIR)/atom_cubin.h
$(OBJ_DIR)/lal_atom.o: lal_atom.cpp lal_atom.h $(NVD_H) $(OBJ_DIR)/atom_cubin.h
$(CUDR) -o $@ -c lal_atom.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_ans.o: lal_answer.cpp lal_answer.h $(NVD_H)
$(CUDR) -o $@ -c lal_answer.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/neighbor_cpu.cubin: lal_neighbor_cpu.cu lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_neighbor_cpu.cu
$(OBJ_DIR)/neighbor_cpu_cubin.h: $(OBJ_DIR)/neighbor_cpu.cubin
$(BIN2C) -c -n neighbor_cpu $(OBJ_DIR)/neighbor_cpu.cubin > $(OBJ_DIR)/neighbor_cpu_cubin.h
$(OBJ_DIR)/neighbor_gpu.cubin: lal_neighbor_gpu.cu lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_neighbor_gpu.cu
$(OBJ_DIR)/neighbor_gpu_cubin.h: $(OBJ_DIR)/neighbor_gpu.cubin
$(BIN2C) -c -n neighbor_gpu $(OBJ_DIR)/neighbor_gpu.cubin > $(OBJ_DIR)/neighbor_gpu_cubin.h
$(OBJ_DIR)/lal_neighbor_shared.o: lal_neighbor_shared.cpp lal_neighbor_shared.h $(OBJ_DIR)/neighbor_cpu_cubin.h $(OBJ_DIR)/neighbor_gpu_cubin.h $(NVD_H)
$(CUDR) -o $@ -c lal_neighbor_shared.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_neighbor.o: lal_neighbor.cpp lal_neighbor.h lal_neighbor_shared.h $(NVD_H)
$(CUDR) -o $@ -c lal_neighbor.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/device.cubin: lal_device.cu lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_device.cu
$(OBJ_DIR)/device_cubin.h: $(OBJ_DIR)/device.cubin
$(BIN2C) -c -n device $(OBJ_DIR)/device.cubin > $(OBJ_DIR)/device_cubin.h
$(OBJ_DIR)/lal_device.o: lal_device.cpp lal_device.h $(ALL_H) $(OBJ_DIR)/device_cubin.h
$(CUDR) -o $@ -c lal_device.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_base_atomic.o: $(ALL_H) lal_base_atomic.h lal_base_atomic.cpp
$(CUDR) -o $@ -c lal_base_atomic.cpp
$(OBJ_DIR)/lal_base_charge.o: $(ALL_H) lal_base_charge.h lal_base_charge.cpp
$(CUDR) -o $@ -c lal_base_charge.cpp
$(OBJ_DIR)/lal_base_ellipsoid.o: $(ALL_H) lal_base_ellipsoid.h lal_base_ellipsoid.cpp $(OBJ_DIR)/ellipsoid_nbor_cubin.h
$(CUDR) -o $@ -c lal_base_ellipsoid.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_base_dipole.o: $(ALL_H) lal_base_dipole.h lal_base_dipole.cpp
$(CUDR) -o $@ -c lal_base_dipole.cpp
$(OBJ_DIR)/lal_base_three.o: $(ALL_H) lal_base_three.h lal_base_three.cpp
$(CUDR) -o $@ -c lal_base_three.cpp
$(OBJ_DIR)/lal_base_dpd.o: $(ALL_H) lal_base_dpd.h lal_base_dpd.cpp
$(CUDR) -o $@ -c lal_base_dpd.cpp
$(OBJ_DIR)/pppm_f.cubin: lal_pppm.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -Dgrdtyp=float -Dgrdtyp4=float4 -o $@ lal_pppm.cu
$(OBJ_DIR)/pppm_f_cubin.h: $(OBJ_DIR)/pppm_f.cubin
$(BIN2C) -c -n pppm_f $(OBJ_DIR)/pppm_f.cubin > $(OBJ_DIR)/pppm_f_cubin.h
$(OBJ_DIR)/pppm_d.cubin: lal_pppm.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -Dgrdtyp=double -Dgrdtyp4=double4 -o $@ lal_pppm.cu
$(OBJ_DIR)/pppm_d_cubin.h: $(OBJ_DIR)/pppm_d.cubin
$(BIN2C) -c -n pppm_d $(OBJ_DIR)/pppm_d.cubin > $(OBJ_DIR)/pppm_d_cubin.h
$(OBJ_DIR)/lal_pppm.o: $(ALL_H) lal_pppm.h lal_pppm.cpp $(OBJ_DIR)/pppm_f_cubin.h $(OBJ_DIR)/pppm_d_cubin.h
$(CUDR) -o $@ -c lal_pppm.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_pppm_ext.o: $(ALL_H) lal_pppm.h lal_pppm_ext.cpp
$(CUDR) -o $@ -c lal_pppm_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/ellipsoid_nbor.cubin: lal_ellipsoid_nbor.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_ellipsoid_nbor.cu
$(OBJ_DIR)/ellipsoid_nbor_cubin.h: $(OBJ_DIR)/ellipsoid_nbor.cubin
$(BIN2C) -c -n ellipsoid_nbor $(OBJ_DIR)/ellipsoid_nbor.cubin > $(OBJ_DIR)/ellipsoid_nbor_cubin.h
$(OBJ_DIR)/gayberne.cubin: lal_gayberne.cu lal_precision.h lal_ellipsoid_extra.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_gayberne.cu
$(OBJ_DIR)/gayberne_lj.cubin: lal_gayberne_lj.cu lal_precision.h lal_ellipsoid_extra.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_gayberne_lj.cu
$(OBJ_DIR)/gayberne_cubin.h: $(OBJ_DIR)/gayberne.cubin
$(BIN2C) -c -n gayberne $(OBJ_DIR)/gayberne.cubin > $(OBJ_DIR)/gayberne_cubin.h
$(OBJ_DIR)/gayberne_lj_cubin.h: $(OBJ_DIR)/gayberne_lj.cubin
$(BIN2C) -c -n gayberne_lj $(OBJ_DIR)/gayberne_lj.cubin > $(OBJ_DIR)/gayberne_lj_cubin.h
$(OBJ_DIR)/lal_gayberne.o: $(ALL_H) lal_gayberne.h lal_gayberne.cpp $(OBJ_DIR)/gayberne_cubin.h $(OBJ_DIR)/gayberne_lj_cubin.h $(OBJ_DIR)/lal_base_ellipsoid.o
$(CUDR) -o $@ -c lal_gayberne.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_gayberne_ext.o: $(ALL_H) $(OBJ_DIR)/lal_gayberne.o lal_gayberne_ext.cpp
$(CUDR) -o $@ -c lal_gayberne_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/re_squared.cubin: lal_re_squared.cu lal_precision.h lal_ellipsoid_extra.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_re_squared.cu
$(OBJ_DIR)/re_squared_lj.cubin: lal_re_squared_lj.cu lal_precision.h lal_ellipsoid_extra.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_re_squared_lj.cu
$(OBJ_DIR)/re_squared_cubin.h: $(OBJ_DIR)/re_squared.cubin
$(BIN2C) -c -n re_squared $(OBJ_DIR)/re_squared.cubin > $(OBJ_DIR)/re_squared_cubin.h
$(OBJ_DIR)/re_squared_lj_cubin.h: $(OBJ_DIR)/re_squared_lj.cubin
$(BIN2C) -c -n re_squared_lj $(OBJ_DIR)/re_squared_lj.cubin > $(OBJ_DIR)/re_squared_lj_cubin.h
$(OBJ_DIR)/lal_re_squared.o: $(ALL_H) lal_re_squared.h lal_re_squared.cpp $(OBJ_DIR)/re_squared_cubin.h $(OBJ_DIR)/re_squared_lj_cubin.h $(OBJ_DIR)/lal_base_ellipsoid.o
$(CUDR) -o $@ -c lal_re_squared.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_re_squared_ext.o: $(ALL_H) $(OBJ_DIR)/lal_re_squared.o lal_re_squared_ext.cpp
$(CUDR) -o $@ -c lal_re_squared_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj.cubin: lal_lj.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj.cu
$(OBJ_DIR)/lj_cubin.h: $(OBJ_DIR)/lj.cubin $(OBJ_DIR)/lj.cubin
$(BIN2C) -c -n lj $(OBJ_DIR)/lj.cubin > $(OBJ_DIR)/lj_cubin.h
$(OBJ_DIR)/lal_lj.o: $(ALL_H) lal_lj.h lal_lj.cpp $(OBJ_DIR)/lj_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_lj.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_ext.o: $(ALL_H) lal_lj.h lal_lj_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_lj_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_coul.cubin: lal_lj_coul.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_coul.cu
$(OBJ_DIR)/lj_coul_cubin.h: $(OBJ_DIR)/lj_coul.cubin $(OBJ_DIR)/lj_coul.cubin
$(BIN2C) -c -n lj_coul $(OBJ_DIR)/lj_coul.cubin > $(OBJ_DIR)/lj_coul_cubin.h
$(OBJ_DIR)/lal_lj_coul.o: $(ALL_H) lal_lj_coul.h lal_lj_coul.cpp $(OBJ_DIR)/lj_coul_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_coul.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_coul_ext.o: $(ALL_H) lal_lj_coul.h lal_lj_coul_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_coul_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_class2_long.cubin: lal_lj_class2_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_class2_long.cu
$(OBJ_DIR)/lj_class2_long_cubin.h: $(OBJ_DIR)/lj_class2_long.cubin $(OBJ_DIR)/lj_class2_long.cubin
$(BIN2C) -c -n lj_class2_long $(OBJ_DIR)/lj_class2_long.cubin > $(OBJ_DIR)/lj_class2_long_cubin.h
$(OBJ_DIR)/lal_lj_class2_long.o: $(ALL_H) lal_lj_class2_long.h lal_lj_class2_long.cpp $(OBJ_DIR)/lj_class2_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_class2_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_class2_long_ext.o: $(ALL_H) lal_lj_class2_long.h lal_lj_class2_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_class2_long_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/coul_long.cubin: lal_coul_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_coul_long.cu
$(OBJ_DIR)/coul_long_cubin.h: $(OBJ_DIR)/coul_long.cubin $(OBJ_DIR)/coul_long.cubin
$(BIN2C) -c -n coul_long $(OBJ_DIR)/coul_long.cubin > $(OBJ_DIR)/coul_long_cubin.h
$(OBJ_DIR)/lal_coul_long.o: $(ALL_H) lal_coul_long.h lal_coul_long.cpp $(OBJ_DIR)/coul_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_coul_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_coul_long_ext.o: $(ALL_H) lal_coul_long.h lal_coul_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_coul_long_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_coul_long.cubin: lal_lj_coul_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_coul_long.cu
$(OBJ_DIR)/lj_coul_long_cubin.h: $(OBJ_DIR)/lj_coul_long.cubin $(OBJ_DIR)/lj_coul_long.cubin
$(BIN2C) -c -n lj_coul_long $(OBJ_DIR)/lj_coul_long.cubin > $(OBJ_DIR)/lj_coul_long_cubin.h
$(OBJ_DIR)/lal_lj_coul_long.o: $(ALL_H) lal_lj_coul_long.h lal_lj_coul_long.cpp $(OBJ_DIR)/lj_coul_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_coul_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_coul_long_ext.o: $(ALL_H) lal_lj_coul_long.h lal_lj_coul_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_coul_long_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_dsf.cubin: lal_lj_dsf.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_dsf.cu
$(OBJ_DIR)/lj_dsf_cubin.h: $(OBJ_DIR)/lj_dsf.cubin $(OBJ_DIR)/lj_dsf.cubin
$(BIN2C) -c -n lj_dsf $(OBJ_DIR)/lj_dsf.cubin > $(OBJ_DIR)/lj_dsf_cubin.h
$(OBJ_DIR)/lal_lj_dsf.o: $(ALL_H) lal_lj_dsf.h lal_lj_dsf.cpp $(OBJ_DIR)/lj_dsf_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_dsf.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_dsf_ext.o: $(ALL_H) lal_lj_dsf.h lal_lj_dsf_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_dsf_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/morse.cubin: lal_morse.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_morse.cu
$(OBJ_DIR)/morse_cubin.h: $(OBJ_DIR)/morse.cubin $(OBJ_DIR)/morse.cubin
$(BIN2C) -c -n morse $(OBJ_DIR)/morse.cubin > $(OBJ_DIR)/morse_cubin.h
$(OBJ_DIR)/lal_morse.o: $(ALL_H) lal_morse.h lal_morse.cpp $(OBJ_DIR)/morse_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_morse.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_morse_ext.o: $(ALL_H) lal_morse.h lal_morse_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_morse_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/charmm_long.cubin: lal_charmm_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_charmm_long.cu
$(OBJ_DIR)/charmm_long_cubin.h: $(OBJ_DIR)/charmm_long.cubin $(OBJ_DIR)/charmm_long.cubin
$(BIN2C) -c -n charmm_long $(OBJ_DIR)/charmm_long.cubin > $(OBJ_DIR)/charmm_long_cubin.h
$(OBJ_DIR)/lal_charmm_long.o: $(ALL_H) lal_charmm_long.h lal_charmm_long.cpp $(OBJ_DIR)/charmm_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_charmm_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_charmm_long_ext.o: $(ALL_H) lal_charmm_long.h lal_charmm_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_charmm_long_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj96.cubin: lal_lj96.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj96.cu
$(OBJ_DIR)/lj96_cubin.h: $(OBJ_DIR)/lj96.cubin $(OBJ_DIR)/lj96.cubin
$(BIN2C) -c -n lj96 $(OBJ_DIR)/lj96.cubin > $(OBJ_DIR)/lj96_cubin.h
$(OBJ_DIR)/lal_lj96.o: $(ALL_H) lal_lj96.h lal_lj96.cpp $(OBJ_DIR)/lj96_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_lj96.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj96_ext.o: $(ALL_H) lal_lj96.h lal_lj96_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_lj96_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_expand.cubin: lal_lj_expand.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_expand.cu
$(OBJ_DIR)/lj_expand_cubin.h: $(OBJ_DIR)/lj_expand.cubin $(OBJ_DIR)/lj_expand.cubin
$(BIN2C) -c -n lj_expand $(OBJ_DIR)/lj_expand.cubin > $(OBJ_DIR)/lj_expand_cubin.h
$(OBJ_DIR)/lal_lj_expand.o: $(ALL_H) lal_lj_expand.h lal_lj_expand.cpp $(OBJ_DIR)/lj_expand_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_lj_expand.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_expand_ext.o: $(ALL_H) lal_lj_expand.h lal_lj_expand_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_lj_expand_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_sdk.cubin: lal_lj_sdk.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_sdk.cu
$(OBJ_DIR)/lj_sdk_cubin.h: $(OBJ_DIR)/lj_sdk.cubin $(OBJ_DIR)/lj_sdk.cubin
$(BIN2C) -c -n lj_sdk $(OBJ_DIR)/lj_sdk.cubin > $(OBJ_DIR)/lj_sdk_cubin.h
$(OBJ_DIR)/lal_lj_sdk.o: $(ALL_H) lal_lj_sdk.h lal_lj_sdk.cpp $(OBJ_DIR)/lj_sdk_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_lj_sdk.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_sdk_ext.o: $(ALL_H) lal_lj_sdk.h lal_lj_sdk_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_lj_sdk_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_sdk_long.cubin: lal_lj_sdk_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_sdk_long.cu
$(OBJ_DIR)/lj_sdk_long_cubin.h: $(OBJ_DIR)/lj_sdk_long.cubin $(OBJ_DIR)/lj_sdk_long.cubin
$(BIN2C) -c -n lj_sdk_long $(OBJ_DIR)/lj_sdk_long.cubin > $(OBJ_DIR)/lj_sdk_long_cubin.h
$(OBJ_DIR)/lal_lj_sdk_long.o: $(ALL_H) lal_lj_sdk_long.h lal_lj_sdk_long.cpp $(OBJ_DIR)/lj_sdk_long_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_lj_sdk_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_sdk_long_ext.o: $(ALL_H) lal_lj_sdk_long.h lal_lj_sdk_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_sdk_long_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/eam.cubin: lal_eam.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_eam.cu
$(OBJ_DIR)/eam_cubin.h: $(OBJ_DIR)/eam.cubin $(OBJ_DIR)/eam.cubin
$(BIN2C) -c -n eam $(OBJ_DIR)/eam.cubin > $(OBJ_DIR)/eam_cubin.h
$(OBJ_DIR)/lal_eam.o: $(ALL_H) lal_eam.h lal_eam.cpp $(OBJ_DIR)/eam_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_eam.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_eam_ext.o: $(ALL_H) lal_eam.h lal_eam_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_eam_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_eam_fs_ext.o: $(ALL_H) lal_eam.h lal_eam_fs_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_eam_fs_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_eam_alloy_ext.o: $(ALL_H) lal_eam.h lal_eam_alloy_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_eam_alloy_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/buck.cubin: lal_buck.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_buck.cu
$(OBJ_DIR)/buck_cubin.h: $(OBJ_DIR)/buck.cubin $(OBJ_DIR)/buck.cubin
$(BIN2C) -c -n buck $(OBJ_DIR)/buck.cubin > $(OBJ_DIR)/buck_cubin.h
$(OBJ_DIR)/lal_buck.o: $(ALL_H) lal_buck.h lal_buck.cpp $(OBJ_DIR)/buck_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_buck.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_buck_ext.o: $(ALL_H) lal_buck.h lal_buck_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_buck_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/buck_coul.cubin: lal_buck_coul.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_buck_coul.cu
$(OBJ_DIR)/buck_coul_cubin.h: $(OBJ_DIR)/buck_coul.cubin $(OBJ_DIR)/buck_coul.cubin
$(BIN2C) -c -n buck_coul $(OBJ_DIR)/buck_coul.cubin > $(OBJ_DIR)/buck_coul_cubin.h
$(OBJ_DIR)/lal_buck_coul.o: $(ALL_H) lal_buck_coul.h lal_buck_coul.cpp $(OBJ_DIR)/buck_coul_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_buck_coul.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_buck_coul_ext.o: $(ALL_H) lal_buck_coul.h lal_buck_coul_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_buck_coul_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/buck_coul_long.cubin: lal_buck_coul_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_buck_coul_long.cu
$(OBJ_DIR)/buck_coul_long_cubin.h: $(OBJ_DIR)/buck_coul_long.cubin $(OBJ_DIR)/buck_coul_long.cubin
$(BIN2C) -c -n buck_coul_long $(OBJ_DIR)/buck_coul_long.cubin > $(OBJ_DIR)/buck_coul_long_cubin.h
$(OBJ_DIR)/lal_buck_coul_long.o: $(ALL_H) lal_buck_coul_long.h lal_buck_coul_long.cpp $(OBJ_DIR)/buck_coul_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_buck_coul_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_buck_coul_long_ext.o: $(ALL_H) lal_buck_coul_long.h lal_buck_coul_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_buck_coul_long_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/table.cubin: lal_table.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_table.cu
$(OBJ_DIR)/table_cubin.h: $(OBJ_DIR)/table.cubin $(OBJ_DIR)/table.cubin
$(BIN2C) -c -n table $(OBJ_DIR)/table.cubin > $(OBJ_DIR)/table_cubin.h
$(OBJ_DIR)/lal_table.o: $(ALL_H) lal_table.h lal_table.cpp $(OBJ_DIR)/table_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_table.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_table_ext.o: $(ALL_H) lal_table.h lal_table_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_table_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/yukawa.cubin: lal_yukawa.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_yukawa.cu
$(OBJ_DIR)/yukawa_cubin.h: $(OBJ_DIR)/yukawa.cubin $(OBJ_DIR)/yukawa.cubin
$(BIN2C) -c -n yukawa $(OBJ_DIR)/yukawa.cubin > $(OBJ_DIR)/yukawa_cubin.h
$(OBJ_DIR)/lal_yukawa.o: $(ALL_H) lal_yukawa.h lal_yukawa.cpp $(OBJ_DIR)/yukawa_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_yukawa.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_yukawa_ext.o: $(ALL_H) lal_yukawa.h lal_yukawa_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_yukawa_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/born.cubin: lal_born.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_born.cu
$(OBJ_DIR)/born_cubin.h: $(OBJ_DIR)/born.cubin $(OBJ_DIR)/born.cubin
$(BIN2C) -c -n born $(OBJ_DIR)/born.cubin > $(OBJ_DIR)/born_cubin.h
$(OBJ_DIR)/lal_born.o: $(ALL_H) lal_born.h lal_born.cpp $(OBJ_DIR)/born_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_born.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_born_ext.o: $(ALL_H) lal_born.h lal_born_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_born_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/born_coul_wolf.cubin: lal_born_coul_wolf.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_born_coul_wolf.cu
$(OBJ_DIR)/born_coul_wolf_cubin.h: $(OBJ_DIR)/born_coul_wolf.cubin $(OBJ_DIR)/born_coul_wolf.cubin
$(BIN2C) -c -n born_coul_wolf $(OBJ_DIR)/born_coul_wolf.cubin > $(OBJ_DIR)/born_coul_wolf_cubin.h
$(OBJ_DIR)/lal_born_coul_wolf.o: $(ALL_H) lal_born_coul_wolf.h lal_born_coul_wolf.cpp $(OBJ_DIR)/born_coul_wolf_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_born_coul_wolf.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_born_coul_wolf_ext.o: $(ALL_H) lal_born_coul_wolf.h lal_born_coul_wolf_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_born_coul_wolf_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/born_coul_long.cubin: lal_born_coul_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_born_coul_long.cu
$(OBJ_DIR)/born_coul_long_cubin.h: $(OBJ_DIR)/born_coul_long.cubin $(OBJ_DIR)/born_coul_long.cubin
$(BIN2C) -c -n born_coul_long $(OBJ_DIR)/born_coul_long.cubin > $(OBJ_DIR)/born_coul_long_cubin.h
$(OBJ_DIR)/lal_born_coul_long.o: $(ALL_H) lal_born_coul_long.h lal_born_coul_long.cpp $(OBJ_DIR)/born_coul_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_born_coul_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_born_coul_long_ext.o: $(ALL_H) lal_born_coul_long.h lal_born_coul_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_born_coul_long_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/dipole_lj.cubin: lal_dipole_lj.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_dipole_lj.cu
$(OBJ_DIR)/dipole_lj_cubin.h: $(OBJ_DIR)/dipole_lj.cubin $(OBJ_DIR)/dipole_lj.cubin
$(BIN2C) -c -n dipole_lj $(OBJ_DIR)/dipole_lj.cubin > $(OBJ_DIR)/dipole_lj_cubin.h
$(OBJ_DIR)/lal_dipole_lj.o: $(ALL_H) lal_dipole_lj.h lal_dipole_lj.cpp $(OBJ_DIR)/dipole_lj_cubin.h $(OBJ_DIR)/lal_base_dipole.o
$(CUDR) -o $@ -c lal_dipole_lj.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_dipole_lj_ext.o: $(ALL_H) lal_dipole_lj.h lal_dipole_lj_ext.cpp lal_base_dipole.h
$(CUDR) -o $@ -c lal_dipole_lj_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/dipole_lj_sf.cubin: lal_dipole_lj_sf.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_dipole_lj_sf.cu
$(OBJ_DIR)/dipole_lj_sf_cubin.h: $(OBJ_DIR)/dipole_lj_sf.cubin $(OBJ_DIR)/dipole_lj_sf.cubin
$(BIN2C) -c -n dipole_lj_sf $(OBJ_DIR)/dipole_lj_sf.cubin > $(OBJ_DIR)/dipole_lj_sf_cubin.h
$(OBJ_DIR)/lal_dipole_lj_sf.o: $(ALL_H) lal_dipole_lj_sf.h lal_dipole_lj_sf.cpp $(OBJ_DIR)/dipole_lj_sf_cubin.h $(OBJ_DIR)/lal_base_dipole.o
$(CUDR) -o $@ -c lal_dipole_lj_sf.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_dipole_lj_sf_ext.o: $(ALL_H) lal_dipole_lj_sf.h lal_dipole_lj_sf_ext.cpp lal_base_dipole.h
$(CUDR) -o $@ -c lal_dipole_lj_sf_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/colloid.cubin: lal_colloid.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_colloid.cu
$(OBJ_DIR)/colloid_cubin.h: $(OBJ_DIR)/colloid.cubin $(OBJ_DIR)/colloid.cubin
$(BIN2C) -c -n colloid $(OBJ_DIR)/colloid.cubin > $(OBJ_DIR)/colloid_cubin.h
$(OBJ_DIR)/lal_colloid.o: $(ALL_H) lal_colloid.h lal_colloid.cpp $(OBJ_DIR)/colloid_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_colloid.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_colloid_ext.o: $(ALL_H) lal_colloid.h lal_colloid_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_colloid_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/gauss.cubin: lal_gauss.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_gauss.cu
$(OBJ_DIR)/gauss_cubin.h: $(OBJ_DIR)/gauss.cubin $(OBJ_DIR)/gauss.cubin
$(BIN2C) -c -n gauss $(OBJ_DIR)/gauss.cubin > $(OBJ_DIR)/gauss_cubin.h
$(OBJ_DIR)/lal_gauss.o: $(ALL_H) lal_gauss.h lal_gauss.cpp $(OBJ_DIR)/gauss_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_gauss.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_gauss_ext.o: $(ALL_H) lal_gauss.h lal_gauss_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_gauss_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/yukawa_colloid.cubin: lal_yukawa_colloid.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_yukawa_colloid.cu
$(OBJ_DIR)/yukawa_colloid_cubin.h: $(OBJ_DIR)/yukawa_colloid.cubin $(OBJ_DIR)/yukawa_colloid.cubin
$(BIN2C) -c -n yukawa_colloid $(OBJ_DIR)/yukawa_colloid.cubin > $(OBJ_DIR)/yukawa_colloid_cubin.h
$(OBJ_DIR)/lal_yukawa_colloid.o: $(ALL_H) lal_yukawa_colloid.h lal_yukawa_colloid.cpp $(OBJ_DIR)/yukawa_colloid_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_yukawa_colloid.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_yukawa_colloid_ext.o: $(ALL_H) lal_yukawa_colloid.h lal_yukawa_colloid_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_yukawa_colloid_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_coul_debye.cubin: lal_lj_coul_debye.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_coul_debye.cu
$(OBJ_DIR)/lj_coul_debye_cubin.h: $(OBJ_DIR)/lj_coul_debye.cubin $(OBJ_DIR)/lj_coul_debye.cubin
$(BIN2C) -c -n lj_coul_debye $(OBJ_DIR)/lj_coul_debye.cubin > $(OBJ_DIR)/lj_coul_debye_cubin.h
$(OBJ_DIR)/lal_lj_coul_debye.o: $(ALL_H) lal_lj_coul_debye.h lal_lj_coul_debye.cpp $(OBJ_DIR)/lj_coul_debye_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_coul_debye.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_coul_debye_ext.o: $(ALL_H) lal_lj_coul_debye.h lal_lj_coul_debye_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_coul_debye_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/coul_dsf.cubin: lal_coul_dsf.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_coul_dsf.cu
$(OBJ_DIR)/coul_dsf_cubin.h: $(OBJ_DIR)/coul_dsf.cubin $(OBJ_DIR)/coul_dsf.cubin
$(BIN2C) -c -n coul_dsf $(OBJ_DIR)/coul_dsf.cubin > $(OBJ_DIR)/coul_dsf_cubin.h
$(OBJ_DIR)/lal_coul_dsf.o: $(ALL_H) lal_coul_dsf.h lal_coul_dsf.cpp $(OBJ_DIR)/coul_dsf_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_coul_dsf.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_coul_dsf_ext.o: $(ALL_H) lal_coul_dsf.h lal_coul_dsf_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_coul_dsf_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/vashishta.cubin: lal_vashishta.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_vashishta.cu
$(OBJ_DIR)/vashishta_cubin.h: $(OBJ_DIR)/vashishta.cubin $(OBJ_DIR)/vashishta.cubin
$(BIN2C) -c -n vashishta $(OBJ_DIR)/vashishta.cubin > $(OBJ_DIR)/vashishta_cubin.h
$(OBJ_DIR)/lal_vashishta.o: $(ALL_H) lal_vashishta.h lal_vashishta.cpp $(OBJ_DIR)/vashishta_cubin.h $(OBJ_DIR)/lal_base_three.o
$(CUDR) -o $@ -c lal_vashishta.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_vashishta_ext.o: $(ALL_H) lal_vashishta.h lal_vashishta_ext.cpp lal_base_three.h
$(CUDR) -o $@ -c lal_vashishta_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/sw.cubin: lal_sw.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_sw.cu
$(OBJ_DIR)/sw_cubin.h: $(OBJ_DIR)/sw.cubin $(OBJ_DIR)/sw.cubin
$(BIN2C) -c -n sw $(OBJ_DIR)/sw.cubin > $(OBJ_DIR)/sw_cubin.h
$(OBJ_DIR)/lal_sw.o: $(ALL_H) lal_sw.h lal_sw.cpp $(OBJ_DIR)/sw_cubin.h $(OBJ_DIR)/lal_base_three.o
$(CUDR) -o $@ -c lal_sw.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_sw_ext.o: $(ALL_H) lal_sw.h lal_sw_ext.cpp lal_base_three.h
$(CUDR) -o $@ -c lal_sw_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/beck.cubin: lal_beck.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_beck.cu
$(OBJ_DIR)/beck_cubin.h: $(OBJ_DIR)/beck.cubin $(OBJ_DIR)/beck.cubin
$(BIN2C) -c -n beck $(OBJ_DIR)/beck.cubin > $(OBJ_DIR)/beck_cubin.h
$(OBJ_DIR)/lal_beck.o: $(ALL_H) lal_beck.h lal_beck.cpp $(OBJ_DIR)/beck_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_beck.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_beck_ext.o: $(ALL_H) lal_beck.h lal_beck_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_beck_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/mie.cubin: lal_mie.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_mie.cu
$(OBJ_DIR)/mie_cubin.h: $(OBJ_DIR)/mie.cubin $(OBJ_DIR)/mie.cubin
$(BIN2C) -c -n mie $(OBJ_DIR)/mie.cubin > $(OBJ_DIR)/mie_cubin.h
$(OBJ_DIR)/lal_mie.o: $(ALL_H) lal_mie.h lal_mie.cpp $(OBJ_DIR)/mie_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_mie.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_mie_ext.o: $(ALL_H) lal_mie.h lal_mie_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_mie_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/soft.cubin: lal_soft.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_soft.cu
$(OBJ_DIR)/soft_cubin.h: $(OBJ_DIR)/soft.cubin $(OBJ_DIR)/soft.cubin
$(BIN2C) -c -n soft $(OBJ_DIR)/soft.cubin > $(OBJ_DIR)/soft_cubin.h
$(OBJ_DIR)/lal_soft.o: $(ALL_H) lal_soft.h lal_soft.cpp $(OBJ_DIR)/soft_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_soft.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_soft_ext.o: $(ALL_H) lal_soft.h lal_soft_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_soft_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_coul_msm.cubin: lal_lj_coul_msm.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_coul_msm.cu
$(OBJ_DIR)/lj_coul_msm_cubin.h: $(OBJ_DIR)/lj_coul_msm.cubin $(OBJ_DIR)/lj_coul_msm.cubin
$(BIN2C) -c -n lj_coul_msm $(OBJ_DIR)/lj_coul_msm.cubin > $(OBJ_DIR)/lj_coul_msm_cubin.h
$(OBJ_DIR)/lal_lj_coul_msm.o: $(ALL_H) lal_lj_coul_msm.h lal_lj_coul_msm.cpp $(OBJ_DIR)/lj_coul_msm_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_coul_msm.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_coul_msm_ext.o: $(ALL_H) lal_lj_coul_msm.h lal_lj_coul_msm_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_coul_msm_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_gromacs.cubin: lal_lj_gromacs.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_gromacs.cu
$(OBJ_DIR)/lj_gromacs_cubin.h: $(OBJ_DIR)/lj_gromacs.cubin $(OBJ_DIR)/lj_gromacs.cubin
$(BIN2C) -c -n lj_gromacs $(OBJ_DIR)/lj_gromacs.cubin > $(OBJ_DIR)/lj_gromacs_cubin.h
$(OBJ_DIR)/lal_lj_gromacs.o: $(ALL_H) lal_lj_gromacs.h lal_lj_gromacs.cpp $(OBJ_DIR)/lj_gromacs_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_lj_gromacs.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_gromacs_ext.o: $(ALL_H) lal_lj_gromacs.h lal_lj_gromacs_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_lj_gromacs_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/dpd.cubin: lal_dpd.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_dpd.cu
$(OBJ_DIR)/dpd_cubin.h: $(OBJ_DIR)/dpd.cubin $(OBJ_DIR)/dpd.cubin
$(BIN2C) -c -n dpd $(OBJ_DIR)/dpd.cubin > $(OBJ_DIR)/dpd_cubin.h
$(OBJ_DIR)/ufm.cubin: lal_ufm.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_ufm.cu
$(OBJ_DIR)/ufm_cubin.h: $(OBJ_DIR)/ufm.cubin $(OBJ_DIR)/ufm.cubin
$(BIN2C) -c -n ufm $(OBJ_DIR)/ufm.cubin > $(OBJ_DIR)/ufm_cubin.h
$(OBJ_DIR)/lal_ufm.o: $(ALL_H) lal_ufm.h lal_ufm.cpp $(OBJ_DIR)/ufm_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_ufm.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_ufm_ext.o: $(ALL_H) lal_ufm.h lal_ufm_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_ufm_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_dpd.o: $(ALL_H) lal_dpd.h lal_dpd.cpp $(OBJ_DIR)/dpd_cubin.h $(OBJ_DIR)/lal_base_dpd.o
$(CUDR) -o $@ -c lal_dpd.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_dpd_ext.o: $(ALL_H) lal_dpd.h lal_dpd_ext.cpp lal_base_dpd.h
$(CUDR) -o $@ -c lal_dpd_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/tersoff.cubin: lal_tersoff.cu lal_precision.h lal_tersoff_extra.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_tersoff.cu
$(OBJ_DIR)/tersoff_cubin.h: $(OBJ_DIR)/tersoff.cubin $(OBJ_DIR)/tersoff.cubin
$(BIN2C) -c -n tersoff $(OBJ_DIR)/tersoff.cubin > $(OBJ_DIR)/tersoff_cubin.h
$(OBJ_DIR)/lal_tersoff.o: $(ALL_H) lal_tersoff.h lal_tersoff.cpp $(OBJ_DIR)/tersoff_cubin.h $(OBJ_DIR)/lal_base_three.o
$(CUDR) -o $@ -c lal_tersoff.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_tersoff_ext.o: $(ALL_H) lal_tersoff.h lal_tersoff_ext.cpp lal_base_three.h
$(CUDR) -o $@ -c lal_tersoff_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/tersoff_zbl.cubin: lal_tersoff_zbl.cu lal_precision.h lal_tersoff_zbl_extra.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_tersoff_zbl.cu
$(OBJ_DIR)/tersoff_zbl_cubin.h: $(OBJ_DIR)/tersoff_zbl.cubin $(OBJ_DIR)/tersoff_zbl.cubin
$(BIN2C) -c -n tersoff_zbl $(OBJ_DIR)/tersoff_zbl.cubin > $(OBJ_DIR)/tersoff_zbl_cubin.h
$(OBJ_DIR)/lal_tersoff_zbl.o: $(ALL_H) lal_tersoff_zbl.h lal_tersoff_zbl.cpp $(OBJ_DIR)/tersoff_zbl_cubin.h $(OBJ_DIR)/lal_base_three.o
$(CUDR) -o $@ -c lal_tersoff_zbl.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_tersoff_zbl_ext.o: $(ALL_H) lal_tersoff_zbl.h lal_tersoff_zbl_ext.cpp lal_base_three.h
$(CUDR) -o $@ -c lal_tersoff_zbl_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/tersoff_mod.cubin: lal_tersoff_mod.cu lal_precision.h lal_tersoff_mod_extra.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_tersoff_mod.cu
$(OBJ_DIR)/tersoff_mod_cubin.h: $(OBJ_DIR)/tersoff_mod.cubin $(OBJ_DIR)/tersoff_mod.cubin
$(BIN2C) -c -n tersoff_mod $(OBJ_DIR)/tersoff_mod.cubin > $(OBJ_DIR)/tersoff_mod_cubin.h
$(OBJ_DIR)/lal_tersoff_mod.o: $(ALL_H) lal_tersoff_mod.h lal_tersoff_mod.cpp $(OBJ_DIR)/tersoff_mod_cubin.h $(OBJ_DIR)/lal_base_three.o
$(CUDR) -o $@ -c lal_tersoff_mod.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_tersoff_mod_ext.o: $(ALL_H) lal_tersoff_mod.h lal_tersoff_mod_ext.cpp lal_base_three.h
$(CUDR) -o $@ -c lal_tersoff_mod_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/coul.cubin: lal_coul.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_coul.cu
$(OBJ_DIR)/coul_cubin.h: $(OBJ_DIR)/coul.cubin $(OBJ_DIR)/coul.cubin
$(BIN2C) -c -n coul $(OBJ_DIR)/coul.cubin > $(OBJ_DIR)/coul_cubin.h
$(OBJ_DIR)/lal_coul.o: $(ALL_H) lal_coul.h lal_coul.cpp $(OBJ_DIR)/coul_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_coul.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_coul_ext.o: $(ALL_H) lal_coul.h lal_coul_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_coul_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/coul_debye.cubin: lal_coul_debye.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_coul_debye.cu
$(OBJ_DIR)/coul_debye_cubin.h: $(OBJ_DIR)/coul_debye.cubin $(OBJ_DIR)/coul_debye.cubin
$(BIN2C) -c -n coul_debye $(OBJ_DIR)/coul_debye.cubin > $(OBJ_DIR)/coul_debye_cubin.h
$(OBJ_DIR)/lal_coul_debye.o: $(ALL_H) lal_coul_debye.h lal_coul_debye.cpp $(OBJ_DIR)/coul_debye_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_coul_debye.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_coul_debye_ext.o: $(ALL_H) lal_coul_debye.h lal_coul_debye_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_coul_debye_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/zbl.cubin: lal_zbl.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_zbl.cu
$(OBJ_DIR)/zbl_cubin.h: $(OBJ_DIR)/zbl.cubin $(OBJ_DIR)/zbl.cubin
$(BIN2C) -c -n zbl $(OBJ_DIR)/zbl.cubin > $(OBJ_DIR)/zbl_cubin.h
$(OBJ_DIR)/lal_zbl.o: $(ALL_H) lal_zbl.h lal_zbl.cpp $(OBJ_DIR)/zbl_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_zbl.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_zbl_ext.o: $(ALL_H) lal_zbl.h lal_zbl_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_zbl_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_cubic.cubin: lal_lj_cubic.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_cubic.cu
$(OBJ_DIR)/lj_cubic_cubin.h: $(OBJ_DIR)/lj_cubic.cubin $(OBJ_DIR)/lj_cubic.cubin
$(BIN2C) -c -n lj_cubic $(OBJ_DIR)/lj_cubic.cubin > $(OBJ_DIR)/lj_cubic_cubin.h
$(OBJ_DIR)/lal_lj_cubic.o: $(ALL_H) lal_lj_cubic.h lal_lj_cubic.cpp $(OBJ_DIR)/lj_cubic_cubin.h $(OBJ_DIR)/lal_base_atomic.o
$(CUDR) -o $@ -c lal_lj_cubic.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_cubic_ext.o: $(ALL_H) lal_lj_cubic.h lal_lj_cubic_ext.cpp lal_base_atomic.h
$(CUDR) -o $@ -c lal_lj_cubic_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/dipole_long_lj.cubin: lal_dipole_long_lj.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_dipole_long_lj.cu
$(OBJ_DIR)/dipole_long_lj_cubin.h: $(OBJ_DIR)/dipole_long_lj.cubin $(OBJ_DIR)/dipole_long_lj.cubin
$(BIN2C) -c -n dipole_long_lj $(OBJ_DIR)/dipole_long_lj.cubin > $(OBJ_DIR)/dipole_long_lj_cubin.h
$(OBJ_DIR)/lal_dipole_long_lj.o: $(ALL_H) lal_dipole_long_lj.h lal_dipole_long_lj.cpp $(OBJ_DIR)/dipole_long_lj_cubin.h $(OBJ_DIR)/lal_base_dipole.o
$(CUDR) -o $@ -c lal_dipole_long_lj.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_dipole_long_lj_ext.o: $(ALL_H) lal_dipole_long_lj.h lal_dipole_long_lj_ext.cpp lal_base_dipole.h
$(CUDR) -o $@ -c lal_dipole_long_lj_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_expand_coul_long.cubin: lal_lj_expand_coul_long.cu lal_precision.h lal_preprocessor.h
$(CUDA) --fatbin -DNV_KERNEL -o $@ lal_lj_expand_coul_long.cu
$(OBJ_DIR)/lj_expand_coul_long_cubin.h: $(OBJ_DIR)/lj_expand_coul_long.cubin $(OBJ_DIR)/lj_expand_coul_long.cubin
$(BIN2C) -c -n lj_expand_coul_long $(OBJ_DIR)/lj_expand_coul_long.cubin > $(OBJ_DIR)/lj_expand_coul_long_cubin.h
$(OBJ_DIR)/lal_lj_expand_coul_long.o: $(ALL_H) lal_lj_expand_coul_long.h lal_lj_expand_coul_long.cpp $(OBJ_DIR)/lj_expand_coul_long_cubin.h $(OBJ_DIR)/lal_base_charge.o
$(CUDR) -o $@ -c lal_lj_expand_coul_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_expand_coul_long_ext.o: $(ALL_H) lal_lj_expand_coul_long.h lal_lj_expand_coul_long_ext.cpp lal_base_charge.h
$(CUDR) -o $@ -c lal_lj_expand_coul_long_ext.cpp -I$(OBJ_DIR)
$(BIN_DIR)/nvc_get_devices: ./geryon/ucl_get_devices.cpp $(NVD_H)
$(CUDR) -o $@ ./geryon/ucl_get_devices.cpp -DUCL_CUDADR $(CUDA_LIB) -lcuda
$(GPU_LIB): $(OBJS) $(CUDPP)
$(AR) -crusv $(GPU_LIB) $(OBJS) $(CUDPP)
@cp $(EXTRAMAKE) Makefile.lammps
clean:
-rm -f $(EXECS) $(GPU_LIB) $(OBJS) $(CUDPP) $(CBNS) *.linkinfo
veryclean: clean
-rm -rf *~ *.linkinfo
cleanlib:
-rm -f $(EXECS) $(GPU_LIB) $(OBJS) $(CBNS) *.linkinfo

27
lib/gpu/Opencl.makefile
View File

@ -66,7 +66,9 @@ OBJS = $(OBJ_DIR)/lal_atom.o $(OBJ_DIR)/lal_answer.o \
$(OBJ_DIR)/lal_coul_debye.o $(OBJ_DIR)/lal_coul_debye_ext.o \
$(OBJ_DIR)/lal_zbl.o $(OBJ_DIR)/lal_zbl_ext.o \
$(OBJ_DIR)/lal_lj_cubic.o $(OBJ_DIR)/lal_lj_cubic_ext.o \
$(OBJ_DIR)/lal_ufm.o $(OBJ_DIR)/lal_ufm_ext.o
$(OBJ_DIR)/lal_ufm.o $(OBJ_DIR)/lal_ufm_ext.o \
$(OBJ_DIR)/lal_dipole_long_lj.o $(OBJ_DIR)/lal_dipole_long_lj_ext.o \
$(OBJ_DIR)/lal_lj_expand_coul_long.o $(OBJ_DIR)/lal_lj_expand_coul_long_ext.o
KERS = $(OBJ_DIR)/device_cl.h $(OBJ_DIR)/atom_cl.h \
$(OBJ_DIR)/neighbor_cpu_cl.h $(OBJ_DIR)/pppm_cl.h \
@ -95,7 +97,8 @@ KERS = $(OBJ_DIR)/device_cl.h $(OBJ_DIR)/atom_cl.h \
$(OBJ_DIR)/tersoff_mod_cl.h $(OBJ_DIR)/coul_cl.h \
$(OBJ_DIR)/coul_debye_cl.h $(OBJ_DIR)/zbl_cl.h \
$(OBJ_DIR)/lj_cubic_cl.h $(OBJ_DIR)/vashishta_cl.h \
$(OBJ_DIR)/ufm_cl.h
$(OBJ_DIR)/ufm_cl.h $(OBJ_DIR)/dipole_long_lj_cl.h \
$(OBJ_DIR)/lj_expand_coul_long_cl.h
OCL_EXECS = $(BIN_DIR)/ocl_get_devices
@ -588,7 +591,25 @@ $(OBJ_DIR)/lal_ufm.o: $(ALL_H) lal_ufm.h lal_ufm.cpp $(OBJ_DIR)/ufm_cl.h $(OBJ_
$(OBJ_DIR)/lal_ufm_ext.o: $(ALL_H) lal_ufm.h lal_ufm_ext.cpp lal_base_atomic.h
$(OCL) -o $@ -c lal_ufm_ext.cpp -I$(OBJ_DIR)
$(BIN_DIR)/ocl_get_devices: ./geryon/ucl_get_devices.cpp
$(OBJ_DIR)/dipole_long_lj_cl.h: lal_dipole_long_lj.cu $(PRE1_H)
$(BSH) ./geryon/file_to_cstr.sh dipole_long_lj $(PRE1_H) lal_dipole_long_lj.cu $(OBJ_DIR)/dipole_long_lj_cl.h;
$(OBJ_DIR)/lal_dipole_long_lj.o: $(ALL_H) lal_dipole_long_lj.h lal_dipole_long_lj.cpp $(OBJ_DIR)/dipole_long_lj_cl.h $(OBJ_DIR)/lj_expand_coul_long_cl.h $(OBJ_DIR)/lal_base_charge.o
$(OCL) -o $@ -c lal_dipole_long_lj.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_dipole_long_lj_ext.o: $(ALL_H) lal_dipole_long_lj.h lal_dipole_long_lj_ext.cpp lal_base_dipole.h
$(OCL) -o $@ -c lal_dipole_long_lj_ext.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lj_expand_coul_long_cl.h: lal_lj_expand_coul_long.cu $(PRE1_H)
$(BSH) ./geryon/file_to_cstr.sh lj_expand_coul_long $(PRE1_H) lal_lj_expand_coul_long.cu $(OBJ_DIR)/lj_expand_coul_long_cl.h;
$(OBJ_DIR)/lal_lj_expand_coul_long.o: $(ALL_H) lal_lj_expand_coul_long.h lal_lj_expand_coul_long.cpp $(OBJ_DIR)/lj_expand_coul_long_cl.h $(OBJ_DIR)/lj_expand_coul_long_cl.h $(OBJ_DIR)/lal_base_charge.o
$(OCL) -o $@ -c lal_lj_expand_coul_long.cpp -I$(OBJ_DIR)
$(OBJ_DIR)/lal_lj_expand_coul_long_ext.o: $(ALL_H) lal_lj_expand_coul_long.h lal_lj_expand_coul_long_ext.cpp lal_base_charge.h
$(OCL) -o $@ -c lal_lj_expand_coul_long_ext.cpp -I$(OBJ_DIR)
$(BIN_DIR)/ocl_get_devices: ./geryon/ucl_get_devices.cpp $(OCL_H)
$(OCL) -o $@ ./geryon/ucl_get_devices.cpp -DUCL_OPENCL $(OCL_LINK)
$(OCL_LIB): $(OBJS) $(PTXS)

42
lib/gpu/geryon/ocl_device.h
View File

@ -280,6 +280,9 @@ class UCL_Device {
/// Return the OpenCL type for the device
inline cl_device_id & cl_device() { return _cl_device; }
/// Select the platform that has accelerators
inline void set_platform_accelerator(int pid=-1);
private:
int _num_platforms; // Number of platforms
int _platform; // UCL_Device ID for current platform
@ -311,8 +314,8 @@ UCL_Device::UCL_Device() {
return;
} else
_num_platforms=static_cast<int>(nplatforms);
set_platform(0);
// note that platform 0 may not necessarily be associated with accelerators
set_platform_accelerator();
}
UCL_Device::~UCL_Device() {
@ -320,6 +323,7 @@ UCL_Device::~UCL_Device() {
}
void UCL_Device::clear() {
_properties.clear();
if (_device>-1) {
for (size_t i=0; i<_cq.size(); i++) {
CL_DESTRUCT_CALL(clReleaseCommandQueue(_cq.back()));
@ -529,8 +533,15 @@ int UCL_Device::set(int num) {
return create_context();
}
// List all devices along with all properties
// List all devices from all platforms along with all properties
void UCL_Device::print_all(std::ostream &out) {
// --- loop through the platforms
for (int n=0; n<_num_platforms; n++) {
set_platform(n);
out << "\nPlatform " << n << ":\n";
if (num_devices() == 0)
out << "There is no device supporting OpenCL\n";
for (int i=0; i<num_devices(); ++i) {
@ -597,7 +608,30 @@ void UCL_Device::print_all(std::ostream &out) {
<< max_sub_devices(i) << std::endl;
}
}
}
// Select the platform that is associated with accelerators
// if pid < 0, select the first platform
void UCL_Device::set_platform_accelerator(int pid) {
if (pid < 0) {
int found = 0;
for (int n=0; n<_num_platforms; n++) {
set_platform(n);
for (int i=0; i<num_devices(); i++) {
if (_properties[i].device_type==CL_DEVICE_TYPE_CPU ||
_properties[i].device_type==CL_DEVICE_TYPE_GPU ||
_properties[i].device_type==CL_DEVICE_TYPE_ACCELERATOR) {
found = 1;
break;
}
}
if (found) break;
}
} else {
set_platform(pid);
}
}
} // namespace ucl_opencl
#endif

2
lib/gpu/lal_device.cpp
View File

@ -80,7 +80,7 @@ int DeviceT::init_device(MPI_Comm world, MPI_Comm replica, const int first_gpu,
char node_name[MPI_MAX_PROCESSOR_NAME];
char *node_names = new char[MPI_MAX_PROCESSOR_NAME*_world_size];
MPI_Get_processor_name(node_name,&name_length);
MPI_Allgather(&node_name,MPI_MAX_PROCESSOR_NAME,MPI_CHAR,&node_names,
MPI_Allgather(&node_name,MPI_MAX_PROCESSOR_NAME,MPI_CHAR,&node_names[0],
MPI_MAX_PROCESSOR_NAME,MPI_CHAR,_comm_world);
std::string node_string=std::string(node_name);

173
lib/gpu/lal_dipole_long_lj.cpp Normal file
View File

@ -0,0 +1,173 @@
/***************************************************************************
dipole_lj.cpp
-------------------
Trung Dac Nguyen (ORNL)
Class for acceleration of the dipole/cut pair style.
__________________________________________________________________________
This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
__________________________________________________________________________
begin :
email : nguyentd@ornl.gov
***************************************************************************/
#ifdef USE_OPENCL
#include "dipole_long_lj_cl.h"
#elif defined(USE_CUDART)
const char *dipole_long_lj=0;
#else
#include "dipole_long_lj_cubin.h"
#endif
#include "lal_dipole_long_lj.h"
#include <cassert>
using namespace LAMMPS_AL;
#define DipoleLongLJT DipoleLongLJ<numtyp, acctyp>
extern Device<PRECISION,ACC_PRECISION> device;
template <class numtyp, class acctyp>
DipoleLongLJT::DipoleLongLJ() : BaseDipole<numtyp,acctyp>(),
_allocated(false) {
}
template <class numtyp, class acctyp>
DipoleLongLJT::~DipoleLongLJ() {
clear();
}
template <class numtyp, class acctyp>
int DipoleLongLJT::bytes_per_atom(const int max_nbors) const {
return this->bytes_per_atom_atomic(max_nbors);
}
template <class numtyp, class acctyp>
int DipoleLongLJT::init(const int ntypes,
double **host_cutsq, double **host_lj1,
double **host_lj2, double **host_lj3,
double **host_lj4, double **host_offset,
double *host_special_lj, const int nlocal,
const int nall, const int max_nbors,
const int maxspecial, const double cell_size,
const double gpu_split, FILE *_screen,
double **host_cut_ljsq, const double host_cut_coulsq,
double *host_special_coul, const double qqrd2e,
const double g_ewald) {
int success;
success=this->init_atomic(nlocal,nall,max_nbors,maxspecial,cell_size,gpu_split,
_screen,dipole_long_lj,"k_dipole_long_lj");
if (success!=0)
return success;
// If atom type constants fit in shared memory use fast kernel
int lj_types=ntypes;
shared_types=false;
int max_shared_types=this->device->max_shared_types();
if (lj_types<=max_shared_types && this->_block_size>=max_shared_types) {
lj_types=max_shared_types;
shared_types=true;
}
_lj_types=lj_types;
// Allocate a host write buffer for data initialization
UCL_H_Vec<numtyp> host_write(lj_types*lj_types*32,*(this->ucl_device),
UCL_WRITE_ONLY);
for (int i=0; i<lj_types*lj_types; i++)
host_write[i]=0.0;
lj1.alloc(lj_types*lj_types,*(this->ucl_device),UCL_READ_ONLY);
this->atom->type_pack4(ntypes,lj_types,lj1,host_write,host_lj1,host_lj2,
host_cut_ljsq);
lj3.alloc(lj_types*lj_types,*(this->ucl_device),UCL_READ_ONLY);
this->atom->type_pack4(ntypes,lj_types,lj3,host_write,host_lj3,host_lj4,
host_offset);
cutsq.alloc(lj_types*lj_types,*(this->ucl_device),UCL_READ_ONLY);
this->atom->type_pack1(ntypes,lj_types,cutsq,host_write,host_cutsq);
sp_lj.alloc(8,*(this->ucl_device),UCL_READ_ONLY);
for (int i=0; i<4; i++) {
host_write[i]=host_special_lj[i];
host_write[i+4]=host_special_coul[i];
}
ucl_copy(sp_lj,host_write,8,false);
_cut_coulsq=host_cut_coulsq;
_qqrd2e=qqrd2e;
_g_ewald=g_ewald;
_allocated=true;
this->_max_bytes=lj1.row_bytes()+lj3.row_bytes()+cutsq.row_bytes()+
sp_lj.row_bytes();
return 0;
}
template <class numtyp, class acctyp>
void DipoleLongLJT::clear() {
if (!_allocated)
return;
_allocated=false;
lj1.clear();
lj3.clear();
cutsq.clear();
sp_lj.clear();
this->clear_atomic();
}
template <class numtyp, class acctyp>
double DipoleLongLJT::host_memory_usage() const {
return this->host_memory_usage_atomic()+sizeof(DipoleLongLJ<numtyp,acctyp>);
}
// ---------------------------------------------------------------------------
// Calculate energies, forces, and torques
// ---------------------------------------------------------------------------
template <class numtyp, class acctyp>
void DipoleLongLJT::loop(const bool _eflag, const bool _vflag) {
// Compute the block size and grid size to keep all cores busy
const int BX=this->block_size();
int eflag, vflag;
if (_eflag)
eflag=1;
else
eflag=0;
if (_vflag)
vflag=1;
else
vflag=0;
int GX=static_cast<int>(ceil(static_cast<double>(this->ans->inum())/
(BX/this->_threads_per_atom)));
int ainum=this->ans->inum();
int nbor_pitch=this->nbor->nbor_pitch();
this->time_pair.start();
if (shared_types) {
this->k_pair_fast.set_size(GX,BX);
this->k_pair_fast.run(&this->atom->x, &lj1, &lj3, &sp_lj,
&this->nbor->dev_nbor,
&this->_nbor_data->begin(),
&this->ans->force, &this->ans->engv, &eflag, &vflag,
&ainum, &nbor_pitch, &this->atom->q,
&this->atom->quat, &cutsq, &_cut_coulsq,
&_qqrd2e, &_g_ewald, &this->_threads_per_atom);
} else {
this->k_pair.set_size(GX,BX);
this->k_pair.run(&this->atom->x, &lj1, &lj3,
&_lj_types, &sp_lj, &this->nbor->dev_nbor,
&this->_nbor_data->begin(), &this->ans->force,
&this->ans->engv, &eflag, &vflag, &ainum,
&nbor_pitch, &this->atom->q,
&this->atom->quat, &cutsq, &_cut_coulsq,
&_qqrd2e, &_g_ewald, &this->_threads_per_atom);
}
this->time_pair.stop();
}
template class DipoleLongLJ<PRECISION,ACC_PRECISION>;

640
lib/gpu/lal_dipole_long_lj.cu Normal file
View File

@ -0,0 +1,640 @@
// **************************************************************************
// dipole_lj.cu
// -------------------
// Trung Dac Nguyen (ORNL)
//
// Device code for acceleration of the dipole/cut pair style
//
// __________________________________________________________________________
// This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
// __________________________________________________________________________
//
// begin :
// email : nguyentd@ornl.gov
// ***************************************************************************/
#ifdef NV_KERNEL
#include "lal_aux_fun1.h"
#ifndef _DOUBLE_DOUBLE
texture<float4> pos_tex;
texture<float> q_tex;
texture<float4> mu_tex;
#else
texture<int4,1> pos_tex;
texture<int2> q_tex;
texture<int4,1> mu_tex;
#endif
#else
#define pos_tex x_
#define q_tex q_
#define mu_tex mu_
#endif
#if (ARCH < 300)
#define store_answers_tq(f, tor, energy, ecoul, virial, ii, inum, tid, \
t_per_atom, offset, eflag, vflag, ans, engv) \
if (t_per_atom>1) { \
__local acctyp red_acc[8][BLOCK_PAIR]; \
red_acc[0][tid]=f.x; \
red_acc[1][tid]=f.y; \
red_acc[2][tid]=f.z; \
red_acc[3][tid]=tor.x; \
red_acc[4][tid]=tor.y; \
red_acc[5][tid]=tor.z; \
for (unsigned int s=t_per_atom/2; s>0; s>>=1) { \
if (offset < s) { \
for (int r=0; r<6; r++) \
red_acc[r][tid] += red_acc[r][tid+s]; \
} \
} \
f.x=red_acc[0][tid]; \
f.y=red_acc[1][tid]; \
f.z=red_acc[2][tid]; \
tor.x=red_acc[3][tid]; \
tor.y=red_acc[4][tid]; \
tor.z=red_acc[5][tid]; \
if (eflag>0 || vflag>0) { \
for (int r=0; r<6; r++) \
red_acc[r][tid]=virial[r]; \
red_acc[6][tid]=energy; \
red_acc[7][tid]=ecoul; \
for (unsigned int s=t_per_atom/2; s>0; s>>=1) { \
if (offset < s) { \
for (int r=0; r<8; r++) \
red_acc[r][tid] += red_acc[r][tid+s]; \
} \
} \
for (int r=0; r<6; r++) \
virial[r]=red_acc[r][tid]; \
energy=red_acc[6][tid]; \
ecoul=red_acc[7][tid]; \
} \
} \
if (offset==0) { \
int ei=ii; \
if (eflag>0) { \
engv[ei]=energy*(acctyp)0.5; \
ei+=inum; \
engv[ei]=e_coul*(acctyp)0.5; \
ei+=inum; \
} \
if (vflag>0) { \
for (int i=0; i<6; i++) { \
engv[ei]=virial[i]*(acctyp)0.5; \
ei+=inum; \
} \
} \
ans[ii]=f; \
ans[ii+inum]=tor; \
}
#else
#define store_answers_tq(f, tor, energy, e_coul, virial, ii, inum, tid, \
t_per_atom, offset, eflag, vflag, ans, engv) \
if (t_per_atom>1) { \
for (unsigned int s=t_per_atom/2; s>0; s>>=1) { \
f.x += shfl_xor(f.x, s, t_per_atom); \
f.y += shfl_xor(f.y, s, t_per_atom); \
f.z += shfl_xor(f.z, s, t_per_atom); \
tor.x += shfl_xor(tor.x, s, t_per_atom); \
tor.y += shfl_xor(tor.y, s, t_per_atom); \
tor.z += shfl_xor(tor.z, s, t_per_atom); \
energy += shfl_xor(energy, s, t_per_atom); \
e_coul += shfl_xor(e_coul, s, t_per_atom); \
} \
if (vflag>0) { \
for (unsigned int s=t_per_atom/2; s>0; s>>=1) { \
for (int r=0; r<6; r++) \
virial[r] += shfl_xor(virial[r], s, t_per_atom); \
} \
} \
} \
if (offset==0) { \
int ei=ii; \
if (eflag>0) { \
engv[ei]=energy*(acctyp)0.5; \
ei+=inum; \
engv[ei]=e_coul*(acctyp)0.5; \
ei+=inum; \
} \
if (vflag>0) { \
for (int i=0; i<6; i++) { \
engv[ei]=virial[i]*(acctyp)0.5; \
ei+=inum; \
} \
} \
ans[ii]=f; \
ans[ii+inum]=tor; \
}
#endif
#define MY_PIS (acctyp)1.77245385090551602729
__kernel void k_dipole_long_lj(const __global numtyp4 *restrict x_,
const __global numtyp4 *restrict lj1,
const __global numtyp4 *restrict lj3,
const int lj_types,
const __global numtyp *restrict sp_lj_in,
const __global int *dev_nbor,
const __global int *dev_packed,
__global acctyp4 *restrict ans,
__global acctyp *restrict engv,
const int eflag, const int vflag, const int inum,
const int nbor_pitch,
const __global numtyp *restrict q_,
const __global numtyp4 *restrict mu_,
const __global numtyp *restrict cutsq,
const numtyp cut_coulsq, const numtyp qqrd2e,
const numtyp g_ewald, const int t_per_atom) {
int tid, ii, offset;
atom_info(t_per_atom,ii,tid,offset);
__local numtyp sp_lj[8];
sp_lj[0]=sp_lj_in[0];
sp_lj[1]=sp_lj_in[1];
sp_lj[2]=sp_lj_in[2];
sp_lj[3]=sp_lj_in[3];
sp_lj[4]=sp_lj_in[4];
sp_lj[5]=sp_lj_in[5];
sp_lj[6]=sp_lj_in[6];
sp_lj[7]=sp_lj_in[7];
acctyp energy=(acctyp)0;
acctyp e_coul=(acctyp)0;
acctyp4 f;
f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0;
acctyp4 tor;
tor.x=(acctyp)0;
tor.y=(acctyp)0;
tor.z=(acctyp)0;
acctyp virial[6];
for (int i=0; i<6; i++)
virial[i]=(acctyp)0;
numtyp pre1 = numtyp(2.0) * g_ewald / MY_PIS;
numtyp pre2 = numtyp(4.0) * (g_ewald*g_ewald*g_ewald) / MY_PIS;
numtyp pre3 = numtyp(8.0) * (g_ewald*g_ewald*g_ewald*g_ewald*g_ewald) / MY_PIS;
if (ii<inum) {
int nbor, nbor_end;
int i, numj;
__local int n_stride;
nbor_info(dev_nbor,dev_packed,nbor_pitch,t_per_atom,ii,offset,i,numj,
n_stride,nbor_end,nbor);
numtyp4 ix; fetch4(ix,i,pos_tex); //x_[i];
numtyp qtmp; fetch(qtmp,i,q_tex);
numtyp4 mui; fetch4(mui,i,mu_tex); //mu_[i];
int itype=ix.w;
for ( ; nbor<nbor_end; nbor+=n_stride) {
int j=dev_packed[nbor];
numtyp factor_lj, factor_coul;
factor_lj = sp_lj[sbmask(j)];
factor_coul = sp_lj[sbmask(j)+4];
j &= NEIGHMASK;
numtyp4 jx; fetch4(jx,j,pos_tex); //x_[j];
numtyp qj; fetch(qj,j,q_tex);
numtyp4 muj; fetch4(muj,j,mu_tex); //mu_[j];
int jtype=jx.w;
// Compute r12
numtyp delx = ix.x-jx.x;
numtyp dely = ix.y-jx.y;
numtyp delz = ix.z-jx.z;
numtyp rsq = delx*delx+dely*dely+delz*delz;
int mtype=itype*lj_types+jtype;
if (rsq<cutsq[mtype]) {
numtyp r2inv=ucl_recip(rsq);
numtyp force_lj,rinv,r6inv;
numtyp pdotp, pidotr, pjdotr, _erfc;
numtyp g0,g1,g2,b0,b1,b2,b3,d0,d1,d2,d3;
numtyp zdix,zdiy,zdiz,zdjx,zdjy,zdjz,zaix,zaiy,zaiz,zajx,zajy,zajz;
numtyp g0b1_g1b2_g2b3,g0d1_g1d2_g2d3,facm1;
numtyp fdx,fdy,fdz,fax,fay,faz;
acctyp4 forcecoul, ticoul;
acctyp4 force;
forcecoul.x = forcecoul.y = forcecoul.z = (acctyp)0;
ticoul.x = ticoul.y = ticoul.z = (acctyp)0;
if (rsq < lj1[mtype].z) {
r6inv = r2inv*r2inv*r2inv;
force_lj = factor_lj*r6inv*(lj1[mtype].x*r6inv-lj1[mtype].y)*r2inv;
} else force_lj = (numtyp)0.0;
if (rsq < cut_coulsq) {
rinv = ucl_rsqrt(rsq);
numtyp r = ucl_rsqrt(r2inv);
numtyp grij = g_ewald * r;
numtyp expm2 = ucl_exp(-grij*grij);
numtyp t = ucl_recip((numtyp)1.0 + EWALD_P*grij);
_erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
pdotp = mui.x*muj.x + mui.y*muj.y + mui.z*muj.z;
pidotr = mui.x*delx + mui.y*dely + mui.z*delz;
pjdotr = muj.x*delx + muj.y*dely + muj.z*delz;
g0 = qtmp*qj;
g1 = qtmp*pjdotr - qj*pidotr + pdotp;
g2 = -pidotr*pjdotr;
if (factor_coul > (numtyp)0.0) {
b0 = _erfc * rinv;
b1 = (b0 + pre1*expm2) * r2inv;
b2 = ((numtyp)3.0*b1 + pre2*expm2) * r2inv;
b3 = ((numtyp)5.0*b2 + pre3*expm2) * r2inv;
g0b1_g1b2_g2b3 = g0*b1 + g1*b2 + g2*b3;
fdx = delx * g0b1_g1b2_g2b3 -
b1 * (qtmp*muj.x - qj*mui.x) +
b2 * (pjdotr*mui.x + pidotr*muj.x);
fdy = dely * g0b1_g1b2_g2b3 -
b1 * (qtmp*muj.y - qj*mui.y) +
b2 * (pjdotr*mui.y + pidotr*muj.y);
fdz = delz * g0b1_g1b2_g2b3 -
b1 * (qtmp*muj.z - qj*mui.z) +
b2 * (pjdotr*mui.z + pidotr*muj.z);
zdix = delx * (qj*b1 + b2*pjdotr) - b1*muj.x;
zdiy = dely * (qj*b1 + b2*pjdotr) - b1*muj.y;
zdiz = delz * (qj*b1 + b2*pjdotr) - b1*muj.z;
zdjx = delx * (-qtmp*b1 + b2*pidotr) - b1*mui.x;
zdjy = dely * (-qtmp*b1 + b2*pidotr) - b1*mui.y;
zdjz = delz * (-qtmp*b1 + b2*pidotr) - b1*mui.z;
if (factor_coul < (numtyp)1.0) {
fdx *= factor_coul;
fdy *= factor_coul;
fdz *= factor_coul;
zdix *= factor_coul;
zdiy *= factor_coul;
zdiz *= factor_coul;
zdjx *= factor_coul;
zdjy *= factor_coul;
zdjz *= factor_coul;
}
} else {
fdx = fdy = fdz = (numtyp)0.0;
zdix = zdiy = zdiz = (numtyp)0.0;
zdjx = zdjy = zdjz = (numtyp)0.0;
}
if (factor_coul < (numtyp)1.0) {
d0 = (_erfc - (numtyp)1.0) * rinv;
d1 = (d0 + pre1*expm2) * r2inv;
d2 = ((numtyp)3.0*d1 + pre2*expm2) * r2inv;
d3 = ((numtyp)5.0*d2 + pre3*expm2) * r2inv;
g0d1_g1d2_g2d3 = g0*d1 + g1*d2 + g2*d3;
fax = delx * g0d1_g1d2_g2d3 -
d1 * (qtmp*muj.x - qj*mui.x) +
d2 * (pjdotr*mui.x + pidotr*muj.x);
fay = dely * g0d1_g1d2_g2d3 -
d1 * (qtmp*muj.y - qj*mui.y) +
d2 * (pjdotr*mui.y + pidotr*muj.y);
faz = delz * g0d1_g1d2_g2d3 -
d1 * (qtmp*muj.z - qj*mui.z) +
d2 * (pjdotr*mui.z + pidotr*muj.z);
zaix = delx * (qj*d1 + d2*pjdotr) - d1*muj.x;
zaiy = dely * (qj*d1 + d2*pjdotr) - d1*muj.y;
zaiz = delz * (qj*d1 + d2*pjdotr) - d1*muj.z;
zajx = delx * (-qtmp*d1 + d2*pidotr) - d1*mui.x;
zajy = dely * (-qtmp*d1 + d2*pidotr) - d1*mui.y;
zajz = delz * (-qtmp*d1 + d2*pidotr) - d1*mui.z;
if (factor_coul > (numtyp)0.0) {
facm1 = (numtyp)1.0 - factor_coul;
fax *= facm1;
fay *= facm1;
faz *= facm1;
zaix *= facm1;
zaiy *= facm1;
zaiz *= facm1;
zajx *= facm1;
zajy *= facm1;
zajz *= facm1;
}
} else {
fax = fay = faz = (numtyp)0.0;
zaix = zaiy = zaiz = (numtyp)0.0;
zajx = zajy = zajz = (numtyp)0.0;
}
forcecoul.x = fdx + fax;
forcecoul.y = fdy + fay;
forcecoul.z = fdz + faz;
ticoul.x = mui.y*(zdiz + zaiz) - mui.z*(zdiy + zaiy);
ticoul.y = mui.z*(zdix + zaix) - mui.x*(zdiz + zaiz);
ticoul.z = mui.x*(zdiy + zaiy) - mui.y*(zdix + zaix);
} else {
forcecoul.x = forcecoul.y = forcecoul.z = (numtyp)0.0;
ticoul.x = ticoul.y = ticoul.z = (numtyp)0.0;
}
force.x = qqrd2e*forcecoul.x + delx*force_lj;
force.y = qqrd2e*forcecoul.y + dely*force_lj;
force.z = qqrd2e*forcecoul.z + delz*force_lj;
f.x+=force.x;
f.y+=force.y;
f.z+=force.z;
tor.x+=qqrd2e*ticoul.x;
tor.y+=qqrd2e*ticoul.y;
tor.z+=qqrd2e*ticoul.z;
if (eflag>0) {
acctyp e = (acctyp)0.0;
if (rsq < cut_coulsq && factor_coul > (numtyp)0.0) {
e = qqrd2e*(b0*g0 + b1*g1 + b2*g2);
if (factor_coul < (numtyp)1.0) {
e_coul *= factor_coul;
e_coul += ((numtyp)1.0-factor_coul) * qqrd2e * (d0*g0 + d1*g1 + d2*g2);
}
} else e = (acctyp)0.0;
e_coul += e;
if (rsq < lj1[mtype].z) {
e=r6inv*(lj3[mtype].x*r6inv-lj3[mtype].y);
energy+=factor_lj*(e-lj3[mtype].z);
}
}
if (vflag>0) {
virial[0] += delx*force.x;
virial[1] += dely*force.y;
virial[2] += delz*force.z;
virial[3] += delx*force.y;
virial[4] += delx*force.z;
virial[5] += dely*force.z;
}
}
} // for nbor
store_answers_tq(f,tor,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,
vflag,ans,engv);
} // if ii
}
__kernel void k_dipole_long_lj_fast(const __global numtyp4 *restrict x_,
const __global numtyp4 *restrict lj1_in,
const __global numtyp4 *restrict lj3_in,
const __global numtyp *restrict sp_lj_in,
const __global int *dev_nbor,
const __global int *dev_packed,
__global acctyp4 *restrict ans,
__global acctyp *restrict engv,
const int eflag, const int vflag, const int inum,
const int nbor_pitch,
const __global numtyp *restrict q_,
const __global numtyp4 *restrict mu_,
const __global numtyp *restrict _cutsq,
const numtyp cut_coulsq, const numtyp qqrd2e,
const numtyp g_ewald, const int t_per_atom) {
int tid, ii, offset;
atom_info(t_per_atom,ii,tid,offset);
__local numtyp4 lj1[MAX_SHARED_TYPES*MAX_SHARED_TYPES];
__local numtyp4 lj3[MAX_SHARED_TYPES*MAX_SHARED_TYPES];
__local numtyp cutsq[MAX_SHARED_TYPES*MAX_SHARED_TYPES];
__local numtyp sp_lj[8];
if (tid<8)
sp_lj[tid]=sp_lj_in[tid];
if (tid<MAX_SHARED_TYPES*MAX_SHARED_TYPES) {
lj1[tid]=lj1_in[tid];
cutsq[tid]=_cutsq[tid];
if (eflag>0)
lj3[tid]=lj3_in[tid];
}
acctyp energy=(acctyp)0;
acctyp e_coul=(acctyp)0;
acctyp4 f;
f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0;
acctyp4 tor;
tor.x=(acctyp)0;
tor.y=(acctyp)0;
tor.z=(acctyp)0;
acctyp virial[6];
for (int i=0; i<6; i++)
virial[i]=(acctyp)0;
__syncthreads();
numtyp pre1 = numtyp(2.0) * g_ewald / MY_PIS;
numtyp pre2 = numtyp(4.0) * (g_ewald*g_ewald*g_ewald) / MY_PIS;
numtyp pre3 = numtyp(8.0) * (g_ewald*g_ewald*g_ewald*g_ewald*g_ewald) / MY_PIS;
if (ii<inum) {
int nbor, nbor_end;
int i, numj;
__local int n_stride;
nbor_info(dev_nbor,dev_packed,nbor_pitch,t_per_atom,ii,offset,i,numj,
n_stride,nbor_end,nbor);
numtyp4 ix; fetch4(ix,i,pos_tex); //x_[i];
numtyp qtmp; fetch(qtmp,i,q_tex);
numtyp4 mui; fetch4(mui,i,mu_tex); //mu_[i];
int iw=ix.w;
int itype=fast_mul((int)MAX_SHARED_TYPES,iw);
for ( ; nbor<nbor_end; nbor+=n_stride) {
int j=dev_packed[nbor];
numtyp factor_lj, factor_coul;
factor_lj = sp_lj[sbmask(j)];
factor_coul = sp_lj[sbmask(j)+4];
j &= NEIGHMASK;
numtyp4 jx; fetch4(jx,j,pos_tex); //x_[j];
numtyp qj; fetch(qj,j,q_tex);
numtyp4 muj; fetch4(muj,j,mu_tex); //mu_[j];
int mtype=itype+jx.w;
// Compute r12
numtyp delx = ix.x-jx.x;
numtyp dely = ix.y-jx.y;
numtyp delz = ix.z-jx.z;
numtyp rsq = delx*delx+dely*dely+delz*delz;
if (rsq<cutsq[mtype]) {
numtyp r2inv=ucl_recip(rsq);
numtyp force_lj,rinv,r6inv;
numtyp pdotp, pidotr, pjdotr, _erfc;
numtyp g0,g1,g2,b0,b1,b2,b3,d0,d1,d2,d3;
numtyp zdix,zdiy,zdiz,zdjx,zdjy,zdjz,zaix,zaiy,zaiz,zajx,zajy,zajz;
numtyp g0b1_g1b2_g2b3,g0d1_g1d2_g2d3,facm1;
numtyp fdx,fdy,fdz,fax,fay,faz;
acctyp4 forcecoul, ticoul;
acctyp4 force;
forcecoul.x = forcecoul.y = forcecoul.z = (acctyp)0;
ticoul.x = ticoul.y = ticoul.z = (acctyp)0;
if (rsq < lj1[mtype].z) {
r6inv = r2inv*r2inv*r2inv;
force_lj = factor_lj*r6inv*(lj1[mtype].x*r6inv-lj1[mtype].y)*r2inv;
} else force_lj = (numtyp)0.0;
if (rsq < cut_coulsq) {
rinv = ucl_rsqrt(rsq);
numtyp r = ucl_rsqrt(r2inv);
numtyp grij = g_ewald * r;
numtyp expm2 = ucl_exp(-grij*grij);
numtyp t = ucl_recip((numtyp)1.0 + EWALD_P*grij);
_erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
pdotp = mui.x*muj.x + mui.y*muj.y + mui.z*muj.z;
pidotr = mui.x*delx + mui.y*dely + mui.z*delz;
pjdotr = muj.x*delx + muj.y*dely + muj.z*delz;
g0 = qtmp*qj;
g1 = qtmp*pjdotr - qj*pidotr + pdotp;
g2 = -pidotr*pjdotr;
if (factor_coul > (numtyp)0.0) {
b0 = _erfc * rinv;
b1 = (b0 + pre1*expm2) * r2inv;
b2 = ((numtyp)3.0*b1 + pre2*expm2) * r2inv;
b3 = ((numtyp)5.0*b2 + pre3*expm2) * r2inv;
g0b1_g1b2_g2b3 = g0*b1 + g1*b2 + g2*b3;
fdx = delx * g0b1_g1b2_g2b3 -
b1 * (qtmp*muj.x - qj*mui.x) +
b2 * (pjdotr*mui.x + pidotr*muj.x);
fdy = dely * g0b1_g1b2_g2b3 -
b1 * (qtmp*muj.y - qj*mui.y) +
b2 * (pjdotr*mui.y + pidotr*muj.y);
fdz = delz * g0b1_g1b2_g2b3 -
b1 * (qtmp*muj.z - qj*mui.z) +
b2 * (pjdotr*mui.z + pidotr*muj.z);
zdix = delx * (qj*b1 + b2*pjdotr) - b1*muj.x;
zdiy = dely * (qj*b1 + b2*pjdotr) - b1*muj.y;
zdiz = delz * (qj*b1 + b2*pjdotr) - b1*muj.z;
zdjx = delx * (-qtmp*b1 + b2*pidotr) - b1*mui.x;
zdjy = dely * (-qtmp*b1 + b2*pidotr) - b1*mui.y;
zdjz = delz * (-qtmp*b1 + b2*pidotr) - b1*mui.z;
if (factor_coul < (numtyp)1.0) {
fdx *= factor_coul;
fdy *= factor_coul;
fdz *= factor_coul;
zdix *= factor_coul;
zdiy *= factor_coul;
zdiz *= factor_coul;
zdjx *= factor_coul;
zdjy *= factor_coul;
zdjz *= factor_coul;
}
} else {
fdx = fdy = fdz = (numtyp)0.0;
zdix = zdiy = zdiz = (numtyp)0.0;
zdjx = zdjy = zdjz = (numtyp)0.0;
}
if (factor_coul < (numtyp)1.0) {
d0 = (_erfc - (numtyp)1.0) * rinv;
d1 = (d0 + pre1*expm2) * r2inv;
d2 = ((numtyp)3.0*d1 + pre2*expm2) * r2inv;
d3 = ((numtyp)5.0*d2 + pre3*expm2) * r2inv;
g0d1_g1d2_g2d3 = g0*d1 + g1*d2 + g2*d3;
fax = delx * g0d1_g1d2_g2d3 -
d1 * (qtmp*muj.x - qj*mui.x) +
d2 * (pjdotr*mui.x + pidotr*muj.x);
fay = dely * g0d1_g1d2_g2d3 -
d1 * (qtmp*muj.y - qj*mui.y) +
d2 * (pjdotr*mui.y + pidotr*muj.y);
faz = delz * g0d1_g1d2_g2d3 -
d1 * (qtmp*muj.z - qj*mui.z) +
d2 * (pjdotr*mui.z + pidotr*muj.z);
zaix = delx * (qj*d1 + d2*pjdotr) - d1*muj.x;
zaiy = dely * (qj*d1 + d2*pjdotr) - d1*muj.y;
zaiz = delz * (qj*d1 + d2*pjdotr) - d1*muj.z;
zajx = delx * (-qtmp*d1 + d2*pidotr) - d1*mui.x;
zajy = dely * (-qtmp*d1 + d2*pidotr) - d1*mui.y;
zajz = delz * (-qtmp*d1 + d2*pidotr) - d1*mui.z;
if (factor_coul > (numtyp)0.0) {
facm1 = (numtyp)1.0 - factor_coul;
fax *= facm1;
fay *= facm1;
faz *= facm1;
zaix *= facm1;
zaiy *= facm1;
zaiz *= facm1;
zajx *= facm1;
zajy *= facm1;
zajz *= facm1;
}
} else {
fax = fay = faz = (numtyp)0.0;
zaix = zaiy = zaiz = (numtyp)0.0;
zajx = zajy = zajz = (numtyp)0.0;
}
forcecoul.x = fdx + fax;
forcecoul.y = fdy + fay;
forcecoul.z = fdz + faz;
ticoul.x = mui.y*(zdiz + zaiz) - mui.z*(zdiy + zaiy);
ticoul.y = mui.z*(zdix + zaix) - mui.x*(zdiz + zaiz);
ticoul.z = mui.x*(zdiy + zaiy) - mui.y*(zdix + zaix);
} else {
forcecoul.x = forcecoul.y = forcecoul.z = (numtyp)0.0;
ticoul.x = ticoul.y = ticoul.z = (numtyp)0.0;
}
force.x = qqrd2e*forcecoul.x + delx*force_lj;
force.y = qqrd2e*forcecoul.y + dely*force_lj;
force.z = qqrd2e*forcecoul.z + delz*force_lj;
f.x+=force.x;
f.y+=force.y;
f.z+=force.z;
tor.x+=qqrd2e*ticoul.x;
tor.y+=qqrd2e*ticoul.y;
tor.z+=qqrd2e*ticoul.z;
if (eflag>0) {
acctyp e = (acctyp)0.0;
if (rsq < cut_coulsq && factor_coul > (numtyp)0.0) {
e = qqrd2e*(b0*g0 + b1*g1 + b2*g2);
if (factor_coul < (numtyp)1.0) {
e_coul *= factor_coul;
e_coul += ((numtyp)1.0-factor_coul) * qqrd2e * (d0*g0 + d1*g1 + d2*g2);
}
} else e = (acctyp)0.0;
e_coul += e;
if (rsq < lj1[mtype].z) {
e=r6inv*(lj3[mtype].x*r6inv-lj3[mtype].y);
energy+=factor_lj*(e-lj3[mtype].z);
}
}
if (vflag>0) {
virial[0] += delx*force.x;
virial[1] += dely*force.y;
virial[2] += delz*force.z;
virial[3] += delx*force.y;
virial[4] += delx*force.z;
virial[5] += dely*force.z;
}
}
} // for nbor
store_answers_tq(f,tor,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,
vflag,ans,engv);
} // if ii
}

85
lib/gpu/lal_dipole_long_lj.h Normal file
View File

@ -0,0 +1,85 @@
/***************************************************************************
dipole_long_lj.h
-------------------
Trung Dac Nguyen (Northwestern)
Class for acceleration of the lj/cut/dipole/long pair style.
__________________________________________________________________________
This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
__________________________________________________________________________
begin :
email : ndactrung@gmail.com
***************************************************************************/
#ifndef LAL_DIPOLE_LONG_LJ_H
#define LAL_DIPOLE_LONG_LJ_H
#include "lal_base_dipole.h"
namespace LAMMPS_AL {
template <class numtyp, class acctyp>
class DipoleLongLJ : public BaseDipole<numtyp, acctyp> {
public:
DipoleLongLJ();
~DipoleLongLJ();
/// Clear any previous data and set up for a new LAMMPS run
/** \param max_nbors initial number of rows in the neighbor matrix
* \param cell_size cutoff + skin
* \param gpu_split fraction of particles handled by device
*
* Returns:
* - 0 if successfull
* - -1 if fix gpu not found
* - -3 if there is an out of memory error
* - -4 if the GPU library was not compiled for GPU
* - -5 Double precision is not supported on card **/
int init(const int ntypes, double **host_cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **host_offset, double *host_special_lj,
const int nlocal, const int nall, const int max_nbors,
const int maxspecial, const double cell_size,
const double gpu_split, FILE *screen, double **host_cut_ljsq,
const double host_cut_coulsq, double *host_special_coul,
const double qqrd2e, const double g_ewald);
/// Clear all host and device data
/** \note This is called at the beginning of the init() routine **/
void clear();
/// Returns memory usage on device per atom
int bytes_per_atom(const int max_nbors) const;
/// Total host memory used by library for pair style
double host_memory_usage() const;
// --------------------------- TYPE DATA --------------------------
/// lj1.x = lj1, lj1.y = lj2, lj1.z = cutsq_vdw, lj1.w = cutsq_coul
UCL_D_Vec<numtyp4> lj1;
/// lj3.x = lj3, lj3.y = lj4, lj3.z = offset
UCL_D_Vec<numtyp4> lj3;
/// cutsq
UCL_D_Vec<numtyp> cutsq;
/// Special LJ values [0-3] and Special Coul values [4-7]
UCL_D_Vec<numtyp> sp_lj;
/// If atom type constants fit in shared memory, use fast kernels
bool shared_types;
/// Number of atom types
int _lj_types;
numtyp _cut_coulsq, _qqrd2e, _g_ewald;
private:
bool _allocated;
void loop(const bool _eflag, const bool _vflag);
};
}
#endif

129
lib/gpu/lal_dipole_long_lj_ext.cpp Normal file
View File

@ -0,0 +1,129 @@
/***************************************************************************
dipole_long_lj_ext.cpp
-------------------
Trung Dac Nguyen (ORNL)
Functions for LAMMPS access to dipole/cut acceleration routines.
__________________________________________________________________________
This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
__________________________________________________________________________
begin :
email : nguyentd@ornl.gov
***************************************************************************/
#include <iostream>
#include <cassert>
#include <cmath>
#include "lal_dipole_long_lj.h"
using namespace std;
using namespace LAMMPS_AL;
static DipoleLongLJ<PRECISION,ACC_PRECISION> DPLJMF;
// ---------------------------------------------------------------------------
// Allocate memory on host and device and copy constants to device
// ---------------------------------------------------------------------------
int dplj_gpu_init(const int ntypes, double **cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **offset, double *special_lj, const int inum,
const int nall, const int max_nbors, const int maxspecial,
const double cell_size, int &gpu_mode, FILE *screen,
double **host_cut_ljsq, const double host_cut_coulsq,
double *host_special_coul, const double qqrd2e,
const double g_ewald) {
DPLJMF.clear();
gpu_mode=DPLJMF.device->gpu_mode();
double gpu_split=DPLJMF.device->particle_split();
int first_gpu=DPLJMF.device->first_device();
int last_gpu=DPLJMF.device->last_device();
int world_me=DPLJMF.device->world_me();
int gpu_rank=DPLJMF.device->gpu_rank();
int procs_per_gpu=DPLJMF.device->procs_per_gpu();
DPLJMF.device->init_message(screen,"lj/cut/dipole/long",first_gpu,last_gpu);
bool message=false;
if (DPLJMF.device->replica_me()==0 && screen)
message=true;
if (message) {
fprintf(screen,"Initializing Device and compiling on process 0...");
fflush(screen);
}
int init_ok=0;
if (world_me==0)
init_ok=DPLJMF.init(ntypes, cutsq, host_lj1, host_lj2, host_lj3,
host_lj4, offset, special_lj, inum, nall, 300,
maxspecial, cell_size, gpu_split, screen, host_cut_ljsq,
host_cut_coulsq, host_special_coul, qqrd2e, g_ewald);
DPLJMF.device->world_barrier();
if (message)
fprintf(screen,"Done.\n");
for (int i=0; i<procs_per_gpu; i++) {
if (message) {
if (last_gpu-first_gpu==0)
fprintf(screen,"Initializing Device %d on core %d...",first_gpu,i);
else
fprintf(screen,"Initializing Devices %d-%d on core %d...",first_gpu,
last_gpu,i);
fflush(screen);
}
if (gpu_rank==i && world_me!=0)
init_ok=DPLJMF.init(ntypes, cutsq, host_lj1, host_lj2, host_lj3, host_lj4,
offset, special_lj, inum, nall, 300, maxspecial,
cell_size, gpu_split, screen, host_cut_ljsq,
host_cut_coulsq, host_special_coul, qqrd2e, g_ewald);
DPLJMF.device->gpu_barrier();
if (message)
fprintf(screen,"Done.\n");
}
if (message)
fprintf(screen,"\n");
if (init_ok==0)
DPLJMF.estimate_gpu_overhead();
return init_ok;
}
void dplj_gpu_clear() {
DPLJMF.clear();
}
int** dplj_gpu_compute_n(const int ago, const int inum_full,
const int nall, double **host_x, int *host_type,
double *sublo, double *subhi, tagint *tag, int **nspecial,
tagint **special, const bool eflag, const bool vflag,
const bool eatom, const bool vatom, int &host_start,
int **ilist, int **jnum, const double cpu_time,
bool &success, double *host_q, double **host_mu,
double *boxlo, double *prd) {
return DPLJMF.compute(ago, inum_full, nall, host_x, host_type, sublo,
subhi, tag, nspecial, special, eflag, vflag, eatom,
vatom, host_start, ilist, jnum, cpu_time, success,
host_q, host_mu, boxlo, prd);
}
void dplj_gpu_compute(const int ago, const int inum_full, const int nall,
double **host_x, int *host_type, int *ilist, int *numj,
int **firstneigh, const bool eflag, const bool vflag,
const bool eatom, const bool vatom, int &host_start,
const double cpu_time, bool &success, double *host_q,
double **host_mu, const int nlocal, double *boxlo, double *prd) {
DPLJMF.compute(ago,inum_full,nall,host_x,host_type,ilist,numj,firstneigh,eflag,
vflag,eatom,vatom,host_start,cpu_time,success,host_q,host_mu,
nlocal,boxlo,prd);
}
double dplj_gpu_bytes() {
return DPLJMF.host_memory_usage();
}

183
lib/gpu/lal_lj_expand_coul_long.cpp Normal file
View File

@ -0,0 +1,183 @@
/***************************************************************************
lj_expand_coul_long.cpp
--------------------------
Trung Nguyen (Northwestern)
Class for acceleration of the lj/expand/coul/long pair style.
__________________________________________________________________________
This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
__________________________________________________________________________
begin :
email : trung.nguyen@northwestern.edu
***************************************************************************/
#if defined(USE_OPENCL)
#include "lj_expand_coul_long_cl.h"
#elif defined(USE_CUDART)
const char *lj_expand_coul_long=0;
#else
#include "lj_expand_coul_long_cubin.h"
#endif
#include "lal_lj_expand_coul_long.h"
#include <cassert>
using namespace LAMMPS_AL;
#define LJExpandCoulLongT LJExpandCoulLong<numtyp, acctyp>
extern Device<PRECISION,ACC_PRECISION> device;
template <class numtyp, class acctyp>
LJExpandCoulLongT::LJExpandCoulLong() : BaseCharge<numtyp,acctyp>(),
_allocated(false) {
}
template <class numtyp, class acctyp>
LJExpandCoulLongT::~LJExpandCoulLong() {
clear();
}
template <class numtyp, class acctyp>
int LJExpandCoulLongT::bytes_per_atom(const int max_nbors) const {
return this->bytes_per_atom_atomic(max_nbors);
}
template <class numtyp, class acctyp>
int LJExpandCoulLongT::init(const int ntypes,
double **host_cutsq, double **host_lj1,
double **host_lj2, double **host_lj3,
double **host_lj4, double **host_offset, double **host_shift,
double *host_special_lj, const int nlocal,
const int nall, const int max_nbors,
const int maxspecial, const double cell_size,
const double gpu_split, FILE *_screen,
double **host_cut_ljsq, const double host_cut_coulsq,
double *host_special_coul, const double qqrd2e,
const double g_ewald) {
int success;
success=this->init_atomic(nlocal,nall,max_nbors,maxspecial,cell_size,gpu_split,
_screen,lj_expand_coul_long,"k_lj_expand_coul_long");
if (success!=0)
return success;
// If atom type constants fit in shared memory use fast kernel
int lj_types=ntypes;
shared_types=false;
int max_shared_types=this->device->max_shared_types();
if (lj_types<=max_shared_types && this->_block_size>=max_shared_types) {
lj_types=max_shared_types;
shared_types=true;
}
_lj_types=lj_types;
// Allocate a host write buffer for data initialization
UCL_H_Vec<numtyp> host_write(lj_types*lj_types*32,*(this->ucl_device),
UCL_WRITE_ONLY);
for (int i=0; i<lj_types*lj_types; i++)
host_write[i]=0.0;
lj1.alloc(lj_types*lj_types,*(this->ucl_device),UCL_READ_ONLY);
this->atom->type_pack4(ntypes,lj_types,lj1,host_write,host_lj1,host_lj2,
host_cutsq, host_cut_ljsq);
lj3.alloc(lj_types*lj_types,*(this->ucl_device),UCL_READ_ONLY);
this->atom->type_pack4(ntypes,lj_types,lj3,host_write,host_lj3,host_lj4,
host_offset,host_shift);
sp_lj.alloc(8,*(this->ucl_device),UCL_READ_ONLY);
for (int i=0; i<4; i++) {
host_write[i]=host_special_lj[i];
host_write[i+4]=host_special_coul[i];
}
ucl_copy(sp_lj,host_write,8,false);
_cut_coulsq=host_cut_coulsq;
_qqrd2e=qqrd2e;
_g_ewald=g_ewald;
_allocated=true;
this->_max_bytes=lj1.row_bytes()+lj3.row_bytes()+sp_lj.row_bytes();
return 0;
}
template <class numtyp, class acctyp>
void LJExpandCoulLongT::reinit(const int ntypes, double **host_cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **host_offset, double **host_shift, double **host_cut_ljsq) {
// Allocate a host write buffer for data initialization
UCL_H_Vec<numtyp> host_write(_lj_types*_lj_types*32,*(this->ucl_device),
UCL_WRITE_ONLY);
for (int i=0; i<_lj_types*_lj_types; i++)
host_write[i]=0.0;
this->atom->type_pack4(ntypes,_lj_types,lj1,host_write,host_lj1,host_lj2,
host_cutsq, host_cut_ljsq);
this->atom->type_pack4(ntypes,_lj_types,lj3,host_write,host_lj3,host_lj4,
host_offset,host_shift);
}
template <class numtyp, class acctyp>
void LJExpandCoulLongT::clear() {
if (!_allocated)
return;
_allocated=false;
lj1.clear();
lj3.clear();
sp_lj.clear();
this->clear_atomic();
}
template <class numtyp, class acctyp>
double LJExpandCoulLongT::host_memory_usage() const {
return this->host_memory_usage_atomic()+sizeof(LJExpandCoulLong<numtyp,acctyp>);
}
// ---------------------------------------------------------------------------
// Calculate energies, forces, and torques
// ---------------------------------------------------------------------------
template <class numtyp, class acctyp>
void LJExpandCoulLongT::loop(const bool _eflag, const bool _vflag) {
// Compute the block size and grid size to keep all cores busy
const int BX=this->block_size();
int eflag, vflag;
if (_eflag)
eflag=1;
else
eflag=0;
if (_vflag)
vflag=1;
else
vflag=0;
int GX=static_cast<int>(ceil(static_cast<double>(this->ans->inum())/
(BX/this->_threads_per_atom)));
int ainum=this->ans->inum();
int nbor_pitch=this->nbor->nbor_pitch();
this->time_pair.start();
if (shared_types) {
this->k_pair_fast.set_size(GX,BX);
this->k_pair_fast.run(&this->atom->x, &lj1, &lj3, &sp_lj,
&this->nbor->dev_nbor, &this->_nbor_data->begin(),
&this->ans->force, &this->ans->engv, &eflag,
&vflag, &ainum, &nbor_pitch, &this->atom->q,
&_cut_coulsq, &_qqrd2e, &_g_ewald,
&this->_threads_per_atom);
} else {
this->k_pair.set_size(GX,BX);
this->k_pair.run(&this->atom->x, &lj1, &lj3,
&_lj_types, &sp_lj, &this->nbor->dev_nbor,
&this->_nbor_data->begin(), &this->ans->force,
&this->ans->engv, &eflag, &vflag, &ainum,
&nbor_pitch, &this->atom->q, &_cut_coulsq,
&_qqrd2e, &_g_ewald, &this->_threads_per_atom);
}
this->time_pair.stop();
}
template class LJExpandCoulLong<PRECISION,ACC_PRECISION>;

265
lib/gpu/lal_lj_expand_coul_long.cu Normal file
View File

@ -0,0 +1,265 @@
// **************************************************************************
// lj_coul_long.cu
// -------------------
// Trung Nguyen (Northwestern)
//
// Device code for acceleration of the lj/expand/coul/long pair style
//
// __________________________________________________________________________
// This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
// __________________________________________________________________________
//
// begin :
// email : trung.nguyen@northwestern.edu
// ***************************************************************************/
#ifdef NV_KERNEL
#include "lal_aux_fun1.h"
#ifndef _DOUBLE_DOUBLE
texture<float4> pos_tex;
texture<float> q_tex;
#else
texture<int4,1> pos_tex;
texture<int2> q_tex;
#endif
#else
#define pos_tex x_
#define q_tex q_
#endif
__kernel void k_lj_expand_coul_long(const __global numtyp4 *restrict x_,
const __global numtyp4 *restrict lj1,
const __global numtyp4 *restrict lj3,
const int lj_types,
const __global numtyp *restrict sp_lj_in,
const __global int *dev_nbor,
const __global int *dev_packed,
__global acctyp4 *restrict ans,
__global acctyp *restrict engv,
const int eflag, const int vflag, const int inum,
const int nbor_pitch,
const __global numtyp *restrict q_,
const numtyp cut_coulsq, const numtyp qqrd2e,
const numtyp g_ewald, const int t_per_atom) {
int tid, ii, offset;
atom_info(t_per_atom,ii,tid,offset);
__local numtyp sp_lj[8];
sp_lj[0]=sp_lj_in[0];
sp_lj[1]=sp_lj_in[1];
sp_lj[2]=sp_lj_in[2];
sp_lj[3]=sp_lj_in[3];
sp_lj[4]=sp_lj_in[4];
sp_lj[5]=sp_lj_in[5];
sp_lj[6]=sp_lj_in[6];
sp_lj[7]=sp_lj_in[7];
acctyp energy=(acctyp)0;
acctyp e_coul=(acctyp)0;
acctyp4 f;
f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0;
acctyp virial[6];
for (int i=0; i<6; i++)
virial[i]=(acctyp)0;
if (ii<inum) {
int nbor, nbor_end;
int i, numj;
__local int n_stride;
nbor_info(dev_nbor,dev_packed,nbor_pitch,t_per_atom,ii,offset,i,numj,
n_stride,nbor_end,nbor);
numtyp4 ix; fetch4(ix,i,pos_tex); //x_[i];
numtyp qtmp; fetch(qtmp,i,q_tex);
int itype=ix.w;
for ( ; nbor<nbor_end; nbor+=n_stride) {
int j=dev_packed[nbor];
numtyp factor_lj, factor_coul;
factor_lj = sp_lj[sbmask(j)];
factor_coul = (numtyp)1.0-sp_lj[sbmask(j)+4];
j &= NEIGHMASK;
numtyp4 jx; fetch4(jx,j,pos_tex); //x_[j];
int jtype=jx.w;
// Compute r12
numtyp delx = ix.x-jx.x;
numtyp dely = ix.y-jx.y;
numtyp delz = ix.z-jx.z;
numtyp rsq = delx*delx+dely*dely+delz*delz;
int mtype=itype*lj_types+jtype;
if (rsq<lj1[mtype].z) {
numtyp r2inv=ucl_recip(rsq);
numtyp forcecoul, force_lj, force, r6inv, prefactor, _erfc;
if (rsq < lj1[mtype].w) {
r6inv = r2inv*r2inv*r2inv;
force_lj = factor_lj*r6inv*(lj1[mtype].x*r6inv-lj1[mtype].y);
} else
force_lj = (numtyp)0.0;
if (rsq < cut_coulsq) {
numtyp r = ucl_rsqrt(r2inv);
numtyp grij = g_ewald * r;
numtyp expm2 = ucl_exp(-grij*grij);
numtyp t = ucl_recip((numtyp)1.0 + EWALD_P*grij);
_erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
fetch(prefactor,j,q_tex);
prefactor *= qqrd2e * qtmp/r;
forcecoul = prefactor * (_erfc + EWALD_F*grij*expm2-factor_coul);
} else
forcecoul = (numtyp)0.0;
force = (force_lj + forcecoul) * r2inv;
f.x+=delx*force;
f.y+=dely*force;
f.z+=delz*force;
if (eflag>0) {
if (rsq < cut_coulsq)
e_coul += prefactor*(_erfc-factor_coul);
if (rsq < lj1[mtype].w) {
numtyp e=r6inv*(lj3[mtype].x*r6inv-lj3[mtype].y);
energy+=factor_lj*(e-lj3[mtype].z);
}
}
if (vflag>0) {
virial[0] += delx*delx*force;
virial[1] += dely*dely*force;
virial[2] += delz*delz*force;
virial[3] += delx*dely*force;
virial[4] += delx*delz*force;
virial[5] += dely*delz*force;
}
}
} // for nbor
store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,
vflag,ans,engv);
} // if ii
}
__kernel void k_lj_expand_coul_long_fast(const __global numtyp4 *restrict x_,
const __global numtyp4 *restrict lj1_in,
const __global numtyp4 *restrict lj3_in,
const __global numtyp *restrict sp_lj_in,
const __global int *dev_nbor,
const __global int *dev_packed,
__global acctyp4 *restrict ans,
__global acctyp *restrict engv,
const int eflag, const int vflag,
const int inum, const int nbor_pitch,
const __global numtyp *restrict q_,
const numtyp cut_coulsq, const numtyp qqrd2e,
const numtyp g_ewald, const int t_per_atom) {
int tid, ii, offset;
atom_info(t_per_atom,ii,tid,offset);
__local numtyp4 lj1[MAX_SHARED_TYPES*MAX_SHARED_TYPES];
__local numtyp4 lj3[MAX_SHARED_TYPES*MAX_SHARED_TYPES];
__local numtyp sp_lj[8];
if (tid<8)
sp_lj[tid]=sp_lj_in[tid];
if (tid<MAX_SHARED_TYPES*MAX_SHARED_TYPES) {
lj1[tid]=lj1_in[tid];
if (eflag>0)
lj3[tid]=lj3_in[tid];
}
acctyp energy=(acctyp)0;
acctyp e_coul=(acctyp)0;
acctyp4 f;
f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0;
acctyp virial[6];
for (int i=0; i<6; i++)
virial[i]=(acctyp)0;
__syncthreads();
if (ii<inum) {
int nbor, nbor_end;
int i, numj;
__local int n_stride;
nbor_info(dev_nbor,dev_packed,nbor_pitch,t_per_atom,ii,offset,i,numj,
n_stride,nbor_end,nbor);
numtyp4 ix; fetch4(ix,i,pos_tex); //x_[i];
numtyp qtmp; fetch(qtmp,i,q_tex);
int iw=ix.w;
int itype=fast_mul((int)MAX_SHARED_TYPES,iw);
for ( ; nbor<nbor_end; nbor+=n_stride) {
int j=dev_packed[nbor];
numtyp factor_lj, factor_coul;
factor_lj = sp_lj[sbmask(j)];
factor_coul = (numtyp)1.0-sp_lj[sbmask(j)+4];
j &= NEIGHMASK;
numtyp4 jx; fetch4(jx,j,pos_tex); //x_[j];
int mtype=itype+jx.w;
// Compute r12
numtyp delx = ix.x-jx.x;
numtyp dely = ix.y-jx.y;
numtyp delz = ix.z-jx.z;
numtyp rsq = delx*delx+dely*dely+delz*delz;
if (rsq<lj1[mtype].z) {
numtyp r2inv=ucl_recip(rsq);
numtyp forcecoul, force_lj, force, r6inv, prefactor, _erfc;
if (rsq < lj1[mtype].w) {
r6inv = r2inv*r2inv*r2inv;
force_lj = factor_lj*r6inv*(lj1[mtype].x*r6inv-lj1[mtype].y);
} else
force_lj = (numtyp)0.0;
if (rsq < cut_coulsq) {
numtyp r = ucl_rsqrt(r2inv);
numtyp grij = g_ewald * r;
numtyp expm2 = ucl_exp(-grij*grij);
numtyp t = ucl_recip((numtyp)1.0 + EWALD_P*grij);
_erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
fetch(prefactor,j,q_tex);
prefactor *= qqrd2e * qtmp/r;
forcecoul = prefactor * (_erfc + EWALD_F*grij*expm2-factor_coul);
} else
forcecoul = (numtyp)0.0;
force = (force_lj + forcecoul) * r2inv;
f.x+=delx*force;
f.y+=dely*force;
f.z+=delz*force;
if (eflag>0) {
if (rsq < cut_coulsq)
e_coul += prefactor*(_erfc-factor_coul);
if (rsq < lj1[mtype].w) {
numtyp e=r6inv*(lj3[mtype].x*r6inv-lj3[mtype].y);
energy+=factor_lj*(e-lj3[mtype].z);
}
}
if (vflag>0) {
virial[0] += delx*delx*force;
virial[1] += dely*dely*force;
virial[2] += delz*delz*force;
virial[3] += delx*dely*force;
virial[4] += delx*delz*force;
virial[5] += dely*delz*force;
}
}
} // for nbor
store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,
vflag,ans,engv);
} // if ii
}

88
lib/gpu/lal_lj_expand_coul_long.h Normal file
View File

@ -0,0 +1,88 @@
/***************************************************************************
lj_expand_coul_long.h
-------------------
Trung Nguyen (Northwestern)
Class for acceleration of the lj/expand/coul/long pair style.
__________________________________________________________________________
This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
__________________________________________________________________________
begin :
email : trung.nguyen@northwestern.edu
***************************************************************************/
#ifndef LAL_LJ_EXPAND_COUL_LONG_H
#define LAL_LJ_EXPAND_COUL_LONG_H
#include "lal_base_charge.h"
namespace LAMMPS_AL {
template <class numtyp, class acctyp>
class LJExpandCoulLong : public BaseCharge<numtyp, acctyp> {
public:
LJExpandCoulLong();
~LJExpandCoulLong();
/// Clear any previous data and set up for a new LAMMPS run
/** \param max_nbors initial number of rows in the neighbor matrix
* \param cell_size cutoff + skin
* \param gpu_split fraction of particles handled by device
*
* Returns:
* - 0 if successfull
* - -1 if fix gpu not found
* - -3 if there is an out of memory error
* - -4 if the GPU library was not compiled for GPU
* - -5 Double precision is not supported on card **/
int init(const int ntypes, double **host_cutsq,
double **host_lj1, double **host_lj2, double **host_lj3,
double **host_lj4, double **host_offset, double **host_shift, double *host_special_lj,
const int nlocal, const int nall, const int max_nbors,
const int maxspecial, const double cell_size,
const double gpu_split, FILE *screen, double **host_cut_ljsq,
const double host_cut_coulsq, double *host_special_coul,
const double qqrd2e, const double g_ewald);
/// Send updated coeffs from host to device (to be compatible with fix adapt)
void reinit(const int ntypes, double **host_cutsq,
double **host_lj1, double **host_lj2, double **host_lj3,
double **host_lj4, double **host_offset, double **host_shift, double **host_cut_ljsq);
/// Clear all host and device data
/** \note This is called at the beginning of the init() routine **/
void clear();
/// Returns memory usage on device per atom
int bytes_per_atom(const int max_nbors) const;
/// Total host memory used by library for pair style
double host_memory_usage() const;
// --------------------------- TYPE DATA --------------------------
/// lj1.x = lj1, lj1.y = lj2, lj1.z = cutsq, lj1.w = cutsq_vdw
UCL_D_Vec<numtyp4> lj1;
/// lj3.x = lj3, lj3.y = lj4, lj3.z = offset, lj3.w = shift
UCL_D_Vec<numtyp4> lj3;
/// Special LJ values [0-3] and Special Coul values [4-7]
UCL_D_Vec<numtyp> sp_lj;
/// If atom type constants fit in shared memory, use fast kernels
bool shared_types;
/// Number of atom types
int _lj_types;
numtyp _cut_coulsq, _qqrd2e, _g_ewald;
private:
bool _allocated;
void loop(const bool _eflag, const bool _vflag);
};
}
#endif

152
lib/gpu/lal_lj_expand_coul_long_ext.cpp Normal file
View File

@ -0,0 +1,152 @@
/***************************************************************************
lj_coul_long_ext.cpp
------------------------
Trung Nguyen (Northwestern)
Functions for LAMMPS access to lj/expand/coul/long acceleration routines.
__________________________________________________________________________
This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
__________________________________________________________________________
begin :
email : trung.nguyen@northwestern.edu
***************************************************************************/
#include <iostream>
#include <cassert>
#include <cmath>
#include "lal_lj_expand_coul_long.h"
using namespace std;
using namespace LAMMPS_AL;
static LJExpandCoulLong<PRECISION,ACC_PRECISION> LJECLMF;
// ---------------------------------------------------------------------------
// Allocate memory on host and device and copy constants to device
// ---------------------------------------------------------------------------
int ljecl_gpu_init(const int ntypes, double **cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **offset, double **shift, double *special_lj, const int inum,
const int nall, const int max_nbors, const int maxspecial,
const double cell_size, int &gpu_mode, FILE *screen,
double **host_cut_ljsq, double host_cut_coulsq,
double *host_special_coul, const double qqrd2e,
const double g_ewald) {
LJECLMF.clear();
gpu_mode=LJECLMF.device->gpu_mode();
double gpu_split=LJECLMF.device->particle_split();
int first_gpu=LJECLMF.device->first_device();
int last_gpu=LJECLMF.device->last_device();
int world_me=LJECLMF.device->world_me();
int gpu_rank=LJECLMF.device->gpu_rank();
int procs_per_gpu=LJECLMF.device->procs_per_gpu();
LJECLMF.device->init_message(screen,"lj/expand/coul/long",first_gpu,last_gpu);
bool message=false;
if (LJECLMF.device->replica_me()==0 && screen)
message=true;
if (message) {
fprintf(screen,"Initializing Device and compiling on process 0...");
fflush(screen);
}
int init_ok=0;
if (world_me==0)
init_ok=LJECLMF.init(ntypes, cutsq, host_lj1, host_lj2, host_lj3, host_lj4,
offset, shift, special_lj, inum, nall, 300, maxspecial,
cell_size, gpu_split, screen, host_cut_ljsq,
host_cut_coulsq, host_special_coul, qqrd2e, g_ewald);
LJECLMF.device->world_barrier();
if (message)
fprintf(screen,"Done.\n");
for (int i=0; i<procs_per_gpu; i++) {
if (message) {
if (last_gpu-first_gpu==0)
fprintf(screen,"Initializing Device %d on core %d...",first_gpu,i);
else
fprintf(screen,"Initializing Devices %d-%d on core %d...",first_gpu,
last_gpu,i);
fflush(screen);
}
if (gpu_rank==i && world_me!=0)
init_ok=LJECLMF.init(ntypes, cutsq, host_lj1, host_lj2, host_lj3, host_lj4,
offset, shift, special_lj, inum, nall, 300, maxspecial,
cell_size, gpu_split, screen, host_cut_ljsq,
host_cut_coulsq, host_special_coul, qqrd2e, g_ewald);
LJECLMF.device->gpu_barrier();
if (message)
fprintf(screen,"Done.\n");
}
if (message)
fprintf(screen,"\n");
if (init_ok==0)
LJECLMF.estimate_gpu_overhead();
return init_ok;
}
// ---------------------------------------------------------------------------
// Copy updated coeffs from host to device
// ---------------------------------------------------------------------------
void ljecl_gpu_reinit(const int ntypes, double **cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **offset, double **shift, double **host_cut_ljsq) {
int world_me=LJECLMF.device->world_me();
int gpu_rank=LJECLMF.device->gpu_rank();
int procs_per_gpu=LJECLMF.device->procs_per_gpu();
if (world_me==0)
LJECLMF.reinit(ntypes, cutsq, host_lj1, host_lj2, host_lj3, host_lj4,
offset, shift, host_cut_ljsq);
LJECLMF.device->world_barrier();
for (int i=0; i<procs_per_gpu; i++) {
if (gpu_rank==i && world_me!=0)
LJECLMF.reinit(ntypes, cutsq, host_lj1, host_lj2, host_lj3, host_lj4,
offset, shift, host_cut_ljsq);
LJECLMF.device->gpu_barrier();
}
}
void ljecl_gpu_clear() {
LJECLMF.clear();
}
int** ljecl_gpu_compute_n(const int ago, const int inum_full,
const int nall, double **host_x, int *host_type,
double *sublo, double *subhi, tagint *tag, int **nspecial,
tagint **special, const bool eflag, const bool vflag,
const bool eatom, const bool vatom, int &host_start,
int **ilist, int **jnum, const double cpu_time,
bool &success, double *host_q, double *boxlo,
double *prd) {
return LJECLMF.compute(ago, inum_full, nall, host_x, host_type, sublo,
subhi, tag, nspecial, special, eflag, vflag, eatom,
vatom, host_start, ilist, jnum, cpu_time, success,
host_q, boxlo, prd);
}
void ljecl_gpu_compute(const int ago, const int inum_full, const int nall,
double **host_x, int *host_type, int *ilist, int *numj,
int **firstneigh, const bool eflag, const bool vflag,
const bool eatom, const bool vatom, int &host_start,
const double cpu_time, bool &success, double *host_q,
const int nlocal, double *boxlo, double *prd) {
LJECLMF.compute(ago,inum_full,nall,host_x,host_type,ilist,numj,
firstneigh,eflag,vflag,eatom,vatom,host_start,cpu_time,success,
host_q,nlocal,boxlo,prd);
}
double ljecl_gpu_bytes() {
return LJECLMF.host_memory_usage();
}

2
src/DIPOLE/pair_lj_cut_dipole_long.h
View File

@ -41,7 +41,7 @@ class PairLJCutDipoleLong : public Pair {
void write_restart_settings(FILE *);
void read_restart_settings(FILE *);
private:
protected:
double cut_lj_global;
double **cut_lj,**cut_ljsq;
double cut_coulsq;

2
src/GPU/Install.sh
View File

@ -133,6 +133,8 @@ action pppm_gpu.cpp pppm.cpp
action pppm_gpu.h pppm.cpp
action pair_ufm_gpu.cpp
action pair_ufm_gpu.h
action pair_lj_cut_dipole_long_gpu.cpp pair_lj_cut_dipole_long.cpp
action pair_lj_cut_dipole_long_gpu.h pair_lj_cut_dipole_long.cpp
# edit 2 Makefile.package files to include/exclude package info

434
src/GPU/pair_lj_cut_dipole_long_gpu.cpp Normal file
View File

@ -0,0 +1,434 @@
/* ----------------------------------------------------------------------
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: Trung Dac Nguyen (Northwestern)
------------------------------------------------------------------------- */
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include "pair_lj_cut_dipole_long_gpu.h"
#include "atom.h"
#include "atom_vec.h"
#include "comm.h"
#include "force.h"
#include "kspace.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "integrate.h"
#include "math_const.h"
#include "memory.h"
#include "error.h"
#include "neigh_request.h"
#include "universe.h"
#include "update.h"
#include "domain.h"
#include "gpu_extra.h"
#define EWALD_F 1.12837917
#define EWALD_P 0.3275911
#define A1 0.254829592
#define A2 -0.284496736
#define A3 1.421413741
#define A4 -1.453152027
#define A5 1.061405429
using namespace LAMMPS_NS;
using namespace MathConst;
// External functions from cuda library for atom decomposition
int dplj_gpu_init(const int ntypes, double **cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **offset, double *special_lj, const int nlocal,
const int nall, const int max_nbors, const int maxspecial,
const double cell_size, int &gpu_mode, FILE *screen,
double **host_cut_ljsq, const double host_cut_coulsq,
double *host_special_coul, const double qqrd2e, const double g_ewald);
void dplj_gpu_clear();
int ** dplj_gpu_compute_n(const int ago, const int inum,
const int nall, double **host_x, int *host_type,
double *sublo, double *subhi, tagint *tag,
int **nspecial, tagint **special, const bool eflag,
const bool vflag, const bool eatom, const bool vatom,
int &host_start, int **ilist, int **jnum,
const double cpu_time, bool &success,
double *host_q, double **host_mu,
double *boxlo, double *prd);
void dplj_gpu_compute(const int ago, const int inum,
const int nall, double **host_x, int *host_type,
int *ilist, int *numj, int **firstneigh,
const bool eflag, const bool vflag, const bool eatom,
const bool vatom, int &host_start, const double cpu_time,
bool &success, double *host_q, double **host_mu,
const int nlocal, double *boxlo, double *prd);
double dplj_gpu_bytes();
/* ---------------------------------------------------------------------- */
PairLJCutDipoleLongGPU::PairLJCutDipoleLongGPU(LAMMPS *lmp) : PairLJCutDipoleLong(lmp),
gpu_mode(GPU_FORCE)
{
respa_enable = 0;
reinitflag = 0;
cpu_time = 0.0;
GPU_EXTRA::gpu_ready(lmp->modify, lmp->error);
}
/* ----------------------------------------------------------------------
free all arrays
------------------------------------------------------------------------- */
PairLJCutDipoleLongGPU::~PairLJCutDipoleLongGPU()
{
dplj_gpu_clear();
}
/* ---------------------------------------------------------------------- */
void PairLJCutDipoleLongGPU::compute(int eflag, int vflag)
{
if (eflag || vflag) ev_setup(eflag,vflag);
else evflag = vflag_fdotr = 0;
int nall = atom->nlocal + atom->nghost;
int inum, host_start;
bool success = true;
int *ilist, *numneigh, **firstneigh;
if (gpu_mode != GPU_FORCE) {
inum = atom->nlocal;
firstneigh = dplj_gpu_compute_n(neighbor->ago, inum, nall, atom->x,
atom->type, domain->sublo, domain->subhi,
atom->tag, atom->nspecial, atom->special,
eflag, vflag, eflag_atom, vflag_atom,
host_start, &ilist, &numneigh, cpu_time,
success, atom->q, atom->mu, domain->boxlo,
domain->prd);
} else {
inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
dplj_gpu_compute(neighbor->ago, inum, nall, atom->x, atom->type,
ilist, numneigh, firstneigh, eflag, vflag, eflag_atom,
vflag_atom, host_start, cpu_time, success, atom->q,
atom->mu, atom->nlocal, domain->boxlo, domain->prd);
}
if (!success)
error->one(FLERR,"Insufficient memory on accelerator");
if (host_start<inum) {
cpu_time = MPI_Wtime();
cpu_compute(host_start, inum, eflag, vflag, ilist, numneigh, firstneigh);
cpu_time = MPI_Wtime() - cpu_time;
}
}
/* ----------------------------------------------------------------------
init specific to this pair style
------------------------------------------------------------------------- */
void PairLJCutDipoleLongGPU::init_style()
{
if (!atom->q_flag || !atom->mu_flag || !atom->torque_flag)
error->all(FLERR,"Pair dipole/cut/gpu requires atom attributes q, mu, torque");
if (force->newton_pair)
error->all(FLERR,"Cannot use newton pair with dipole/cut/gpu pair style");
if (strcmp(update->unit_style,"electron") == 0)
error->all(FLERR,"Cannot (yet) use 'electron' units with dipoles");
// Repeat cutsq calculation because done after call to init_style
double maxcut = -1.0;
double cut;
for (int i = 1; i <= atom->ntypes; i++) {
for (int j = i; j <= atom->ntypes; j++) {
if (setflag[i][j] != 0 || (setflag[i][i] != 0 && setflag[j][j] != 0)) {
cut = init_one(i,j);
cut *= cut;
if (cut > maxcut)
maxcut = cut;
cutsq[i][j] = cutsq[j][i] = cut;
} else
cutsq[i][j] = cutsq[j][i] = 0.0;
}
}
double cell_size = sqrt(maxcut) + neighbor->skin;
cut_coulsq = cut_coul * cut_coul;
// insure use of KSpace long-range solver, set g_ewald
if (force->kspace == NULL)
error->all(FLERR,"Pair style requires a KSpace style");
g_ewald = force->kspace->g_ewald;
// setup force tables
if (ncoultablebits) init_tables(cut_coul,NULL);
int maxspecial=0;
if (atom->molecular)
maxspecial=atom->maxspecial;
int success = dplj_gpu_init(atom->ntypes+1, cutsq, lj1, lj2, lj3, lj4,
offset, force->special_lj, atom->nlocal,
atom->nlocal+atom->nghost, 300, maxspecial,
cell_size, gpu_mode, screen, cut_ljsq, cut_coulsq,
force->special_coul, force->qqrd2e, g_ewald);
GPU_EXTRA::check_flag(success,error,world);
if (gpu_mode == GPU_FORCE) {
int irequest = neighbor->request(this,instance_me);
neighbor->requests[irequest]->half = 0;
neighbor->requests[irequest]->full = 1;
}
}
/* ---------------------------------------------------------------------- */
double PairLJCutDipoleLongGPU::memory_usage()
{
double bytes = Pair::memory_usage();
return bytes + dplj_gpu_bytes();
}
/* ---------------------------------------------------------------------- */
void PairLJCutDipoleLongGPU::cpu_compute(int start, int inum, int eflag, int vflag,
int *ilist, int *numneigh,
int **firstneigh)
{
int i,j,ii,jj,jnum,itype,jtype;
double qtmp,xtmp,ytmp,ztmp,delx,dely,delz;
double rsq,r,rinv,r2inv,r6inv;
double forcecoulx,forcecouly,forcecoulz,fforce;
double tixcoul,tiycoul,tizcoul,tjxcoul,tjycoul,tjzcoul;
double fx,fy,fz,fdx,fdy,fdz,fax,fay,faz;
double pdotp,pidotr,pjdotr,pre1,pre2,pre3;
double grij,expm2,t,erfc;
double g0,g1,g2,b0,b1,b2,b3,d0,d1,d2,d3;
double zdix,zdiy,zdiz,zdjx,zdjy,zdjz,zaix,zaiy,zaiz,zajx,zajy,zajz;
double g0b1_g1b2_g2b3,g0d1_g1d2_g2d3;
double forcelj,factor_coul,factor_lj,facm1;
double evdwl,ecoul;
int *jlist;
evdwl = ecoul = 0.0;
if (eflag || vflag) ev_setup(eflag,vflag);
else evflag = vflag_fdotr = 0;
double **x = atom->x;
double **f = atom->f;
double *q = atom->q;
double **mu = atom->mu;
double **torque = atom->torque;
int *type = atom->type;
double *special_coul = force->special_coul;
double *special_lj = force->special_lj;
double qqrd2e = force->qqrd2e;
pre1 = 2.0 * g_ewald / MY_PIS;
pre2 = 4.0 * pow(g_ewald,3.0) / MY_PIS;
pre3 = 8.0 * pow(g_ewald,5.0) / MY_PIS;
// loop over neighbors of my atoms
for (ii = start; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
itype = type[i];
jlist = firstneigh[i];
jnum = numneigh[i];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
factor_lj = special_lj[sbmask(j)];
factor_coul = special_coul[sbmask(j)];
j &= NEIGHMASK;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
jtype = type[j];
if (rsq < cutsq[itype][jtype]) {
r2inv = 1.0/rsq;
rinv = sqrt(r2inv);
if (rsq < cut_coulsq) {
r = sqrt(rsq);
grij = g_ewald * r;
expm2 = exp(-grij*grij);
t = 1.0 / (1.0 + EWALD_P*grij);
erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
pdotp = mu[i][0]*mu[j][0] + mu[i][1]*mu[j][1] + mu[i][2]*mu[j][2];
pidotr = mu[i][0]*delx + mu[i][1]*dely + mu[i][2]*delz;
pjdotr = mu[j][0]*delx + mu[j][1]*dely + mu[j][2]*delz;
g0 = qtmp*q[j];
g1 = qtmp*pjdotr - q[j]*pidotr + pdotp;
g2 = -pidotr*pjdotr;
if (factor_coul > 0.0) {
b0 = erfc * rinv;
b1 = (b0 + pre1*expm2) * r2inv;
b2 = (3.0*b1 + pre2*expm2) * r2inv;
b3 = (5.0*b2 + pre3*expm2) * r2inv;
g0b1_g1b2_g2b3 = g0*b1 + g1*b2 + g2*b3;
fdx = delx * g0b1_g1b2_g2b3 -
b1 * (qtmp*mu[j][0] - q[j]*mu[i][0]) +
b2 * (pjdotr*mu[i][0] + pidotr*mu[j][0]);
fdy = dely * g0b1_g1b2_g2b3 -
b1 * (qtmp*mu[j][1] - q[j]*mu[i][1]) +
b2 * (pjdotr*mu[i][1] + pidotr*mu[j][1]);
fdz = delz * g0b1_g1b2_g2b3 -
b1 * (qtmp*mu[j][2] - q[j]*mu[i][2]) +
b2 * (pjdotr*mu[i][2] + pidotr*mu[j][2]);
zdix = delx * (q[j]*b1 + b2*pjdotr) - b1*mu[j][0];
zdiy = dely * (q[j]*b1 + b2*pjdotr) - b1*mu[j][1];
zdiz = delz * (q[j]*b1 + b2*pjdotr) - b1*mu[j][2];
zdjx = delx * (-qtmp*b1 + b2*pidotr) - b1*mu[i][0];
zdjy = dely * (-qtmp*b1 + b2*pidotr) - b1*mu[i][1];
zdjz = delz * (-qtmp*b1 + b2*pidotr) - b1*mu[i][2];
if (factor_coul < 1.0) {
fdx *= factor_coul;
fdy *= factor_coul;
fdz *= factor_coul;
zdix *= factor_coul;
zdiy *= factor_coul;
zdiz *= factor_coul;
zdjx *= factor_coul;
zdjy *= factor_coul;
zdjz *= factor_coul;
}
} else {
fdx = fdy = fdz = 0.0;
zdix = zdiy = zdiz = 0.0;
zdjx = zdjy = zdjz = 0.0;
}
if (factor_coul < 1.0) {
d0 = (erfc - 1.0) * rinv;
d1 = (d0 + pre1*expm2) * r2inv;
d2 = (3.0*d1 + pre2*expm2) * r2inv;
d3 = (5.0*d2 + pre3*expm2) * r2inv;
g0d1_g1d2_g2d3 = g0*d1 + g1*d2 + g2*d3;
fax = delx * g0d1_g1d2_g2d3 -
d1 * (qtmp*mu[j][0] - q[j]*mu[i][0]) +
d2 * (pjdotr*mu[i][0] + pidotr*mu[j][0]);
fay = dely * g0d1_g1d2_g2d3 -
d1 * (qtmp*mu[j][1] - q[j]*mu[i][1]) +
d2 * (pjdotr*mu[i][1] + pidotr*mu[j][1]);
faz = delz * g0d1_g1d2_g2d3 -
d1 * (qtmp*mu[j][2] - q[j]*mu[i][2]) +
d2 * (pjdotr*mu[i][2] + pidotr*mu[j][2]);
zaix = delx * (q[j]*d1 + d2*pjdotr) - d1*mu[j][0];
zaiy = dely * (q[j]*d1 + d2*pjdotr) - d1*mu[j][1];
zaiz = delz * (q[j]*d1 + d2*pjdotr) - d1*mu[j][2];
zajx = delx * (-qtmp*d1 + d2*pidotr) - d1*mu[i][0];
zajy = dely * (-qtmp*d1 + d2*pidotr) - d1*mu[i][1];
zajz = delz * (-qtmp*d1 + d2*pidotr) - d1*mu[i][2];
if (factor_coul > 0.0) {
facm1 = 1.0 - factor_coul;
fax *= facm1;
fay *= facm1;
faz *= facm1;
zaix *= facm1;
zaiy *= facm1;
zaiz *= facm1;
zajx *= facm1;
zajy *= facm1;
zajz *= facm1;
}
} else {
fax = fay = faz = 0.0;
zaix = zaiy = zaiz = 0.0;
zajx = zajy = zajz = 0.0;
}
forcecoulx = fdx + fax;
forcecouly = fdy + fay;
forcecoulz = fdz + faz;
tixcoul = mu[i][1]*(zdiz + zaiz) - mu[i][2]*(zdiy + zaiy);
tiycoul = mu[i][2]*(zdix + zaix) - mu[i][0]*(zdiz + zaiz);
tizcoul = mu[i][0]*(zdiy + zaiy) - mu[i][1]*(zdix + zaix);
tjxcoul = mu[j][1]*(zdjz + zajz) - mu[j][2]*(zdjy + zajy);
tjycoul = mu[j][2]*(zdjx + zajx) - mu[j][0]*(zdjz + zajz);
tjzcoul = mu[j][0]*(zdjy + zajy) - mu[j][1]*(zdjx + zajx);
} else {
forcecoulx = forcecouly = forcecoulz = 0.0;
tixcoul = tiycoul = tizcoul = 0.0;
tjxcoul = tjycoul = tjzcoul = 0.0;
}
// LJ interaction
if (rsq < cut_ljsq[itype][jtype]) {
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
fforce = factor_lj * forcelj*r2inv;
} else fforce = 0.0;
// total force
fx = qqrd2e*forcecoulx + delx*fforce;
fy = qqrd2e*forcecouly + dely*fforce;
fz = qqrd2e*forcecoulz + delz*fforce;
// force & torque accumulation
f[i][0] += fx;
f[i][1] += fy;
f[i][2] += fz;
torque[i][0] += qqrd2e*tixcoul;
torque[i][1] += qqrd2e*tiycoul;
torque[i][2] += qqrd2e*tizcoul;
if (eflag) {
if (rsq < cut_coulsq && factor_coul > 0.0) {
ecoul = qqrd2e*(b0*g0 + b1*g1 + b2*g2);
if (factor_coul < 1.0) {
ecoul *= factor_coul;
ecoul += (1-factor_coul) * qqrd2e * (d0*g0 + d1*g1 + d2*g2);
}
} else ecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
evdwl = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype]) -
offset[itype][jtype];
evdwl *= factor_lj;
} else evdwl = 0.0;
}
if (evflag) ev_tally_xyz_full(i,evdwl,ecoul,fx,fy,fz,delx,dely,delz);
}
}
}
}

66
src/GPU/pair_lj_cut_dipole_long_gpu.h Normal file
View File

@ -0,0 +1,66 @@
/* -*- 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 PAIR_CLASS
PairStyle(lj/cut/dipole/long/gpu,PairLJCutDipoleLongGPU)
#else
#ifndef LMP_PAIR_LJ_CUT_DIPOLE_LONG_GPU_H
#define LMP_PAIR_LJ_CUT_DIPOLE_LONG_GPU_H
#include "pair_lj_cut_dipole_long.h"
namespace LAMMPS_NS {
class PairLJCutDipoleLongGPU : public PairLJCutDipoleLong {
public:
PairLJCutDipoleLongGPU(LAMMPS *lmp);
~PairLJCutDipoleLongGPU();
void cpu_compute(int, int, int, int, int *, int *, int **);
void compute(int, int);
void init_style();
double memory_usage();
enum { GPU_FORCE, GPU_NEIGH, GPU_HYB_NEIGH };
private:
int gpu_mode;
double cpu_time;
};
}
#endif
#endif
/* ERROR/WARNING messages:
E: Insufficient memory on accelerator
There is insufficient memory on one of the devices specified for the gpu
package
E: Pair lj/cutdipole/long/gpu requires atom attributes q, mu, torque
The atom style defined does not have this attribute.
E: Cannot use newton pair with lj/cutdipole/long/gpu pair style
Self-explanatory.
E: Cannot (yet) use 'electron' units with dipoles
This feature is not yet supported.
*/

325
src/GPU/pair_lj_expand_coul_long_gpu.cpp Normal file
View File

@ -0,0 +1,325 @@
/* ----------------------------------------------------------------------
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: Trung Nguyen (Northwestern)
------------------------------------------------------------------------- */
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include "pair_lj_expand_coul_long_gpu.h"
#include "atom.h"
#include "atom_vec.h"
#include "comm.h"
#include "force.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "integrate.h"
#include "memory.h"
#include "error.h"
#include "neigh_request.h"
#include "universe.h"
#include "update.h"
#include "domain.h"
#include "kspace.h"
#include "gpu_extra.h"
#define EWALD_F 1.12837917
#define EWALD_P 0.3275911
#define A1 0.254829592
#define A2 -0.284496736
#define A3 1.421413741
#define A4 -1.453152027
#define A5 1.061405429
using namespace LAMMPS_NS;
// External functions from cuda library for atom decomposition
int ljecl_gpu_init(const int ntypes, double **cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **offset, double **shift, double *special_lj, const int nlocal,
const int nall, const int max_nbors, const int maxspecial,
const double cell_size, int &gpu_mode, FILE *screen,
double **host_cut_ljsq, double host_cut_coulsq,
double *host_special_coul, const double qqrd2e,
const double g_ewald);
int ljecl_gpu_reinit(const int ntypes, double **cutsq, double **host_lj1,
double **host_lj2, double **host_lj3, double **host_lj4,
double **offset, double **shift, double **host_lj_cutsq);
void ljecl_gpu_clear();
int ** ljecl_gpu_compute_n(const int ago, const int inum,
const int nall, double **host_x, int *host_type,
double *sublo, double *subhi, tagint *tag,
int **nspecial, tagint **special, const bool eflag,
const bool vflag, const bool eatom, const bool vatom,
int &host_start, int **ilist, int **jnum,
const double cpu_time, bool &success, double *host_q,
double *boxlo, double *prd);
void ljecl_gpu_compute(const int ago, const int inum, const int nall,
double **host_x, int *host_type, int *ilist, int *numj,
int **firstneigh, const bool eflag, const bool vflag,
const bool eatom, const bool vatom, int &host_start,
const double cpu_time, bool &success, double *host_q,
const int nlocal, double *boxlo, double *prd);
double ljecl_gpu_bytes();
/* ---------------------------------------------------------------------- */
PairLJExpandCoulLongGPU::PairLJExpandCoulLongGPU(LAMMPS *lmp) :
PairLJExpandCoulLong(lmp), gpu_mode(GPU_FORCE)
{
respa_enable = 0;
cpu_time = 0.0;
GPU_EXTRA::gpu_ready(lmp->modify, lmp->error);
}
/* ----------------------------------------------------------------------
free all arrays
------------------------------------------------------------------------- */
PairLJExpandCoulLongGPU::~PairLJExpandCoulLongGPU()
{
ljecl_gpu_clear();
}
/* ---------------------------------------------------------------------- */
void PairLJExpandCoulLongGPU::compute(int eflag, int vflag)
{
if (eflag || vflag) ev_setup(eflag,vflag);
else evflag = vflag_fdotr = 0;
int nall = atom->nlocal + atom->nghost;
int inum, host_start;
bool success = true;
int *ilist, *numneigh, **firstneigh;
if (gpu_mode != GPU_FORCE) {
inum = atom->nlocal;
firstneigh = ljecl_gpu_compute_n(neighbor->ago, inum, nall, atom->x,
atom->type, domain->sublo, domain->subhi,
atom->tag, atom->nspecial, atom->special,
eflag, vflag, eflag_atom, vflag_atom,
host_start, &ilist, &numneigh, cpu_time,
success, atom->q, domain->boxlo,
domain->prd);
} else {
inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
ljecl_gpu_compute(neighbor->ago, inum, nall, atom->x, atom->type,
ilist, numneigh, firstneigh, eflag, vflag, eflag_atom,
vflag_atom, host_start, cpu_time, success, atom->q,
atom->nlocal, domain->boxlo, domain->prd);
}
if (!success)
error->one(FLERR,"Insufficient memory on accelerator");
if (host_start<inum) {
cpu_time = MPI_Wtime();
cpu_compute(host_start, inum, eflag, vflag, ilist, numneigh, firstneigh);
cpu_time = MPI_Wtime() - cpu_time;
}
}
/* ----------------------------------------------------------------------
init specific to this pair style
------------------------------------------------------------------------- */
void PairLJExpandCoulLongGPU::init_style()
{
cut_respa = NULL;
if (!atom->q_flag)
error->all(FLERR,"Pair style lj/cut/coul/long/gpu requires atom attribute q");
if (force->newton_pair)
error->all(FLERR,"Cannot use newton pair with lj/cut/coul/long/gpu pair style");
// Repeat cutsq calculation because done after call to init_style
double maxcut = -1.0;
double cut;
for (int i = 1; i <= atom->ntypes; i++) {
for (int j = i; j <= atom->ntypes; j++) {
if (setflag[i][j] != 0 || (setflag[i][i] != 0 && setflag[j][j] != 0)) {
cut = init_one(i,j);
cut *= cut;
if (cut > maxcut)
maxcut = cut;
cutsq[i][j] = cutsq[j][i] = cut;
} else
cutsq[i][j] = cutsq[j][i] = 0.0;
}
}
double cell_size = sqrt(maxcut) + neighbor->skin;
cut_coulsq = cut_coul * cut_coul;
// insure use of KSpace long-range solver, set g_ewald
if (force->kspace == NULL)
error->all(FLERR,"Pair style requires a KSpace style");
g_ewald = force->kspace->g_ewald;
// setup force tables
if (ncoultablebits) init_tables(cut_coul,cut_respa);
int maxspecial=0;
if (atom->molecular)
maxspecial=atom->maxspecial;
int success = ljecl_gpu_init(atom->ntypes+1, cutsq, lj1, lj2, lj3, lj4,
offset, shift, force->special_lj, atom->nlocal,
atom->nlocal+atom->nghost, 300, maxspecial,
cell_size, gpu_mode, screen, cut_ljsq, cut_coulsq,
force->special_coul, force->qqrd2e, g_ewald);
GPU_EXTRA::check_flag(success,error,world);
if (gpu_mode == GPU_FORCE) {
int irequest = neighbor->request(this,instance_me);
neighbor->requests[irequest]->half = 0;
neighbor->requests[irequest]->full = 1;
}
}
/* ---------------------------------------------------------------------- */
void PairLJExpandCoulLongGPU::reinit()
{
Pair::reinit();
ljecl_gpu_reinit(atom->ntypes+1, cutsq, lj1, lj2, lj3, lj4, offset, shift, cut_ljsq);
}
/* ---------------------------------------------------------------------- */
double PairLJExpandCoulLongGPU::memory_usage()
{
double bytes = Pair::memory_usage();
return bytes + ljecl_gpu_bytes();
}
/* ---------------------------------------------------------------------- */
void PairLJExpandCoulLongGPU::cpu_compute(int start, int inum, int eflag,
int vflag, int *ilist, int *numneigh,
int **firstneigh)
{
int i,j,ii,jj,jnum,itype,jtype,itable;
double qtmp,xtmp,ytmp,ztmp,delx,dely,delz,evdwl,ecoul,fpair;
double fraction,table;
double r,r2inv,r6inv,forcecoul,forcelj,factor_coul,factor_lj;
double grij,expm2,prefactor,t,erfc;
int *jlist;
double rsq;
evdwl = ecoul = 0.0;
double **x = atom->x;
double **f = atom->f;
double *q = atom->q;
int *type = atom->type;
double *special_coul = force->special_coul;
double *special_lj = force->special_lj;
double qqrd2e = force->qqrd2e;
// loop over neighbors of my atoms
for (ii = start; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
itype = type[i];
jlist = firstneigh[i];
jnum = numneigh[i];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
factor_lj = special_lj[sbmask(j)];
factor_coul = special_coul[sbmask(j)];
j &= NEIGHMASK;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
jtype = type[j];
if (rsq < cutsq[itype][jtype]) {
r2inv = 1.0/rsq;
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq) {
r = sqrt(rsq);
grij = g_ewald * r;
expm2 = exp(-grij*grij);
t = 1.0 / (1.0 + EWALD_P*grij);
erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
prefactor = qqrd2e * qtmp*q[j]/r;
forcecoul = prefactor * (erfc + EWALD_F*grij*expm2);
if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*prefactor;
} else {
union_int_float_t rsq_lookup;
rsq_lookup.f = rsq;
itable = rsq_lookup.i & ncoulmask;
itable >>= ncoulshiftbits;
fraction = (rsq_lookup.f - rtable[itable]) * drtable[itable];
table = ftable[itable] + fraction*dftable[itable];
forcecoul = qtmp*q[j] * table;
if (factor_coul < 1.0) {
table = ctable[itable] + fraction*dctable[itable];
prefactor = qtmp*q[j] * table;
forcecoul -= (1.0-factor_coul)*prefactor;
}
}
} else forcecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
} else forcelj = 0.0;
fpair = (forcecoul + factor_lj*forcelj) * r2inv;
f[i][0] += delx*fpair;
f[i][1] += dely*fpair;
f[i][2] += delz*fpair;
if (eflag) {
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq)
ecoul = prefactor*erfc;
else {
table = etable[itable] + fraction*detable[itable];
ecoul = qtmp*q[j] * table;
}
if (factor_coul < 1.0) ecoul -= (1.0-factor_coul)*prefactor;
} else ecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
evdwl = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype]) -
offset[itype][jtype];
evdwl *= factor_lj;
} else evdwl = 0.0;
}
if (evflag) ev_tally_full(i,evdwl,ecoul,fpair,delx,dely,delz);
}
}
}
}

67
src/GPU/pair_lj_expand_coul_long_gpu.h Normal file
View File

@ -0,0 +1,67 @@
/* -*- 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 PAIR_CLASS
PairStyle(lj/expand/coul/long/gpu,PairLJExpandCoulLongGPU)
#else
#ifndef LMP_PAIR_LJ_EXPAND_COUL_LONG_GPU_H
#define LMP_PAIR_LJ_EXPAND_COUL_LONG_GPU_H
#include "pair_lj_expand_coul_long.h"
namespace LAMMPS_NS {
class PairLJExpandCoulLongGPU : public PairLJExpandCoulLong {
public:
PairLJExpandCoulLongGPU(LAMMPS *lmp);
~PairLJExpandCoulLongGPU();
void cpu_compute(int, int, int, int, int *, int *, int **);
void compute(int, int);
void init_style();
void reinit();
double memory_usage();
enum { GPU_FORCE, GPU_NEIGH, GPU_HYB_NEIGH };
private:
int gpu_mode;
double cpu_time;
};
}
#endif
#endif
/* ERROR/WARNING messages:
E: Insufficient memory on accelerator
There is insufficient memory on one of the devices specified for the gpu
package
E: Pair style lj/cut/coul/long/gpu requires atom attribute q
The atom style defined does not have this attribute.
E: Cannot use newton pair with lj/cut/coul/long/gpu pair style
Self-explanatory.
E: Pair style requires a KSpace style
No kspace style is defined.
*/

990
src/USER-MISC/pair_lj_expand_coul_long.cpp Normal file
View File

@ -0,0 +1,990 @@
/* ----------------------------------------------------------------------
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: Trung Nguyen (Northwestern)
------------------------------------------------------------------------- */
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include "pair_lj_expand_coul_long.h"
#include "atom.h"
#include "comm.h"
#include "force.h"
#include "kspace.h"
#include "update.h"
#include "integrate.h"
#include "respa.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "math_const.h"
#include "memory.h"
#include "error.h"
using namespace LAMMPS_NS;
using namespace MathConst;
#define EWALD_F 1.12837917
#define EWALD_P 0.3275911
#define A1 0.254829592
#define A2 -0.284496736
#define A3 1.421413741
#define A4 -1.453152027
#define A5 1.061405429
/* ---------------------------------------------------------------------- */
PairLJExpandCoulLong::PairLJExpandCoulLong(LAMMPS *lmp) : Pair(lmp)
{
ewaldflag = pppmflag = 1;
respa_enable = 1;
writedata = 1;
ftable = NULL;
qdist = 0.0;
}
/* ---------------------------------------------------------------------- */
PairLJExpandCoulLong::~PairLJExpandCoulLong()
{
if (allocated) {
memory->destroy(setflag);
memory->destroy(cutsq);
memory->destroy(cut_lj);
memory->destroy(cut_ljsq);
memory->destroy(epsilon);
memory->destroy(sigma);
memory->destroy(lj1);
memory->destroy(lj2);
memory->destroy(lj3);
memory->destroy(lj4);
memory->destroy(offset);
memory->destroy(shift);
}
if (ftable) free_tables();
}
/* ---------------------------------------------------------------------- */
void PairLJExpandCoulLong::compute(int eflag, int vflag)
{
int i,ii,j,jj,inum,jnum,itype,jtype,itable;
double qtmp,xtmp,ytmp,ztmp,delx,dely,delz,evdwl,ecoul,fpair;
double fraction,table;
double r,r2inv,r6inv,forcecoul,forcelj,factor_coul,factor_lj;
double grij,expm2,prefactor,t,erfc;
int *ilist,*jlist,*numneigh,**firstneigh;
double rsq,rshift,rshiftsq;
evdwl = ecoul = 0.0;
if (eflag || vflag) ev_setup(eflag,vflag);
else evflag = vflag_fdotr = 0;
double **x = atom->x;
double **f = atom->f;
double *q = atom->q;
int *type = atom->type;
int nlocal = atom->nlocal;
double *special_coul = force->special_coul;
double *special_lj = force->special_lj;
int newton_pair = force->newton_pair;
double qqrd2e = force->qqrd2e;
inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// loop over neighbors of my atoms
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
itype = type[i];
jlist = firstneigh[i];
jnum = numneigh[i];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
factor_lj = special_lj[sbmask(j)];
factor_coul = special_coul[sbmask(j)];
j &= NEIGHMASK;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
jtype = type[j];
if (rsq < cutsq[itype][jtype]) {
r2inv = 1.0/rsq;
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq) {
r = sqrt(rsq);
grij = g_ewald * r;
expm2 = exp(-grij*grij);
t = 1.0 / (1.0 + EWALD_P*grij);
erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
prefactor = qqrd2e * qtmp*q[j]/r;
forcecoul = prefactor * (erfc + EWALD_F*grij*expm2);
if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*prefactor;
} else {
union_int_float_t rsq_lookup;
rsq_lookup.f = rsq;
itable = rsq_lookup.i & ncoulmask;
itable >>= ncoulshiftbits;
fraction = (rsq_lookup.f - rtable[itable]) * drtable[itable];
table = ftable[itable] + fraction*dftable[itable];
forcecoul = qtmp*q[j] * table;
if (factor_coul < 1.0) {
table = ctable[itable] + fraction*dctable[itable];
prefactor = qtmp*q[j] * table;
forcecoul -= (1.0-factor_coul)*prefactor;
}
}
} else forcecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
r = sqrt(rsq);
rshift = r - shift[itype][jtype];
rshiftsq = rshift*rshift;
r2inv = 1.0/rshiftsq;
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
forcelj = factor_lj*forcelj/rshift/r;
} else forcelj = 0.0;
fpair = forcecoul*r2inv + forcelj;
f[i][0] += delx*fpair;
f[i][1] += dely*fpair;
f[i][2] += delz*fpair;
if (newton_pair || j < nlocal) {
f[j][0] -= delx*fpair;
f[j][1] -= dely*fpair;
f[j][2] -= delz*fpair;
}
if (eflag) {
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq)
ecoul = prefactor*erfc;
else {
table = etable[itable] + fraction*detable[itable];
ecoul = qtmp*q[j] * table;
}
if (factor_coul < 1.0) ecoul -= (1.0-factor_coul)*prefactor;
} else ecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
evdwl = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype]) -
offset[itype][jtype];
evdwl *= factor_lj;
} else evdwl = 0.0;
}
if (evflag) ev_tally(i,j,nlocal,newton_pair,
evdwl,ecoul,fpair,delx,dely,delz);
}
}
}
if (vflag_fdotr) virial_fdotr_compute();
}
/* ---------------------------------------------------------------------- */
void PairLJExpandCoulLong::compute_inner()
{
int i,j,ii,jj,inum,jnum,itype,jtype;
double qtmp,xtmp,ytmp,ztmp,delx,dely,delz,fpair;
double rsq,r2inv,r6inv,forcecoul,forcelj,factor_coul,factor_lj;
double rsw,r,rshift,rshiftsq;
int *ilist,*jlist,*numneigh,**firstneigh;
double **x = atom->x;
double **f = atom->f;
double *q = atom->q;
int *type = atom->type;
int nlocal = atom->nlocal;
double *special_coul = force->special_coul;
double *special_lj = force->special_lj;
int newton_pair = force->newton_pair;
double qqrd2e = force->qqrd2e;
inum = list->inum_inner;
ilist = list->ilist_inner;
numneigh = list->numneigh_inner;
firstneigh = list->firstneigh_inner;
double cut_out_on = cut_respa[0];
double cut_out_off = cut_respa[1];
double cut_out_diff = cut_out_off - cut_out_on;
double cut_out_on_sq = cut_out_on*cut_out_on;
double cut_out_off_sq = cut_out_off*cut_out_off;
// loop over neighbors of my atoms
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
itype = type[i];
jlist = firstneigh[i];
jnum = numneigh[i];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
factor_lj = special_lj[sbmask(j)];
factor_coul = special_coul[sbmask(j)];
j &= NEIGHMASK;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
if (rsq < cut_out_off_sq) {
r2inv = 1.0/rsq;
forcecoul = qqrd2e * qtmp*q[j]*sqrt(r2inv);
if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*forcecoul;
jtype = type[j];
if (rsq < cut_ljsq[itype][jtype]) {
r = sqrt(rsq);
rshift = r - shift[itype][jtype];
rshiftsq = rshift*rshift;
r2inv = 1.0/rshiftsq;
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
forcelj = factor_lj*forcelj/rshift/r;
} else forcelj = 0.0;
fpair = forcecoul*r2inv + forcelj;
if (rsq > cut_out_on_sq) {
rsw = (sqrt(rsq) - cut_out_on)/cut_out_diff;
fpair *= 1.0 + rsw*rsw*(2.0*rsw-3.0);
}
f[i][0] += delx*fpair;
f[i][1] += dely*fpair;
f[i][2] += delz*fpair;
if (newton_pair || j < nlocal) {
f[j][0] -= delx*fpair;
f[j][1] -= dely*fpair;
f[j][2] -= delz*fpair;
}
}
}
}
}
/* ---------------------------------------------------------------------- */
void PairLJExpandCoulLong::compute_middle()
{
int i,j,ii,jj,inum,jnum,itype,jtype;
double qtmp,xtmp,ytmp,ztmp,delx,dely,delz,fpair;
double rsq,r2inv,r6inv,forcecoul,forcelj,factor_coul,factor_lj;
double rsw,r,rshift,rshiftsq;
int *ilist,*jlist,*numneigh,**firstneigh;
double **x = atom->x;
double **f = atom->f;
double *q = atom->q;
int *type = atom->type;
int nlocal = atom->nlocal;
double *special_coul = force->special_coul;
double *special_lj = force->special_lj;
int newton_pair = force->newton_pair;
double qqrd2e = force->qqrd2e;
inum = list->inum_middle;
ilist = list->ilist_middle;
numneigh = list->numneigh_middle;
firstneigh = list->firstneigh_middle;
double cut_in_off = cut_respa[0];
double cut_in_on = cut_respa[1];
double cut_out_on = cut_respa[2];
double cut_out_off = cut_respa[3];
double cut_in_diff = cut_in_on - cut_in_off;
double cut_out_diff = cut_out_off - cut_out_on;
double cut_in_off_sq = cut_in_off*cut_in_off;
double cut_in_on_sq = cut_in_on*cut_in_on;
double cut_out_on_sq = cut_out_on*cut_out_on;
double cut_out_off_sq = cut_out_off*cut_out_off;
// loop over neighbors of my atoms
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
itype = type[i];
jlist = firstneigh[i];
jnum = numneigh[i];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
factor_lj = special_lj[sbmask(j)];
factor_coul = special_coul[sbmask(j)];
j &= NEIGHMASK;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
if (rsq < cut_out_off_sq && rsq > cut_in_off_sq) {
r2inv = 1.0/rsq;
forcecoul = qqrd2e * qtmp*q[j]*sqrt(r2inv);
if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*forcecoul;
jtype = type[j];
if (rsq < cut_ljsq[itype][jtype]) {
r = sqrt(rsq);
rshift = r - shift[itype][jtype];
rshiftsq = rshift*rshift;
r2inv = 1.0/rshiftsq;
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
forcelj = factor_lj*forcelj/rshift/r;
} else forcelj = 0.0;
fpair = forcecoul*r2inv + forcelj;
if (rsq < cut_in_on_sq) {
rsw = (sqrt(rsq) - cut_in_off)/cut_in_diff;
fpair *= rsw*rsw*(3.0 - 2.0*rsw);
}
if (rsq > cut_out_on_sq) {
rsw = (sqrt(rsq) - cut_out_on)/cut_out_diff;
fpair *= 1.0 + rsw*rsw*(2.0*rsw - 3.0);
}
f[i][0] += delx*fpair;
f[i][1] += dely*fpair;
f[i][2] += delz*fpair;
if (newton_pair || j < nlocal) {
f[j][0] -= delx*fpair;
f[j][1] -= dely*fpair;
f[j][2] -= delz*fpair;
}
}
}
}
}
/* ---------------------------------------------------------------------- */
void PairLJExpandCoulLong::compute_outer(int eflag, int vflag)
{
int i,j,ii,jj,inum,jnum,itype,jtype,itable;
double qtmp,xtmp,ytmp,ztmp,delx,dely,delz,evdwl,ecoul,fpair;
double fraction,table;
double r,r2inv,r6inv,forcecoul,forcelj,factor_coul,factor_lj;
double grij,expm2,prefactor,t,erfc;
double rsw,rshift,rshiftsq;
int *ilist,*jlist,*numneigh,**firstneigh;
double rsq;
evdwl = ecoul = 0.0;
if (eflag || vflag) ev_setup(eflag,vflag);
else evflag = 0;
double **x = atom->x;
double **f = atom->f;
double *q = atom->q;
int *type = atom->type;
int nlocal = atom->nlocal;
double *special_coul = force->special_coul;
double *special_lj = force->special_lj;
int newton_pair = force->newton_pair;
double qqrd2e = force->qqrd2e;
inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
double cut_in_off = cut_respa[2];
double cut_in_on = cut_respa[3];
double cut_in_diff = cut_in_on - cut_in_off;
double cut_in_off_sq = cut_in_off*cut_in_off;
double cut_in_on_sq = cut_in_on*cut_in_on;
// loop over neighbors of my atoms
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
itype = type[i];
jlist = firstneigh[i];
jnum = numneigh[i];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
factor_lj = special_lj[sbmask(j)];
factor_coul = special_coul[sbmask(j)];
j &= NEIGHMASK;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
jtype = type[j];
if (rsq < cutsq[itype][jtype]) {
r2inv = 1.0/rsq;
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq) {
r = sqrt(rsq);
grij = g_ewald * r;
expm2 = exp(-grij*grij);
t = 1.0 / (1.0 + EWALD_P*grij);
erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
prefactor = qqrd2e * qtmp*q[j]/r;
forcecoul = prefactor * (erfc + EWALD_F*grij*expm2 - 1.0);
if (rsq > cut_in_off_sq) {
if (rsq < cut_in_on_sq) {
rsw = (r - cut_in_off)/cut_in_diff;
forcecoul += prefactor*rsw*rsw*(3.0 - 2.0*rsw);
if (factor_coul < 1.0)
forcecoul -=
(1.0-factor_coul)*prefactor*rsw*rsw*(3.0 - 2.0*rsw);
} else {
forcecoul += prefactor;
if (factor_coul < 1.0)
forcecoul -= (1.0-factor_coul)*prefactor;
}
}
} else {
union_int_float_t rsq_lookup;
rsq_lookup.f = rsq;
itable = rsq_lookup.i & ncoulmask;
itable >>= ncoulshiftbits;
fraction = (rsq_lookup.f - rtable[itable]) * drtable[itable];
table = ftable[itable] + fraction*dftable[itable];
forcecoul = qtmp*q[j] * table;
if (factor_coul < 1.0) {
table = ctable[itable] + fraction*dctable[itable];
prefactor = qtmp*q[j] * table;
forcecoul -= (1.0-factor_coul)*prefactor;
}
}
} else forcecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype] && rsq > cut_in_off_sq) {
r = sqrt(rsq);
rshift = r - shift[itype][jtype];
rshiftsq = rshift*rshift;
r2inv = 1.0/rshiftsq;
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
forcelj = factor_lj*forcelj/rshift/r;
if (rsq < cut_in_on_sq) {
rsw = (sqrt(rsq) - cut_in_off)/cut_in_diff;
forcelj *= rsw*rsw*(3.0 - 2.0*rsw);
}
} else forcelj = 0.0;
fpair = forcecoul*r2inv + forcelj;
f[i][0] += delx*fpair;
f[i][1] += dely*fpair;
f[i][2] += delz*fpair;
if (newton_pair || j < nlocal) {
f[j][0] -= delx*fpair;
f[j][1] -= dely*fpair;
f[j][2] -= delz*fpair;
}
if (eflag) {
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq) {
ecoul = prefactor*erfc;
if (factor_coul < 1.0) ecoul -= (1.0-factor_coul)*prefactor;
} else {
table = etable[itable] + fraction*detable[itable];
ecoul = qtmp*q[j] * table;
if (factor_coul < 1.0) {
table = ptable[itable] + fraction*dptable[itable];
prefactor = qtmp*q[j] * table;
ecoul -= (1.0-factor_coul)*prefactor;
}
}
} else ecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
r6inv = r2inv*r2inv*r2inv;
evdwl = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype]) -
offset[itype][jtype];
evdwl *= factor_lj;
} else evdwl = 0.0;
}
if (vflag) {
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq) {
forcecoul = prefactor * (erfc + EWALD_F*grij*expm2);
if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*prefactor;
} else {
table = vtable[itable] + fraction*dvtable[itable];
forcecoul = qtmp*q[j] * table;
if (factor_coul < 1.0) {
table = ptable[itable] + fraction*dptable[itable];
prefactor = qtmp*q[j] * table;
forcecoul -= (1.0-factor_coul)*prefactor;
}
}
} else forcecoul = 0.0;
if (rsq <= cut_in_off_sq) {
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
} else if (rsq <= cut_in_on_sq) {
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
}
fpair = (forcecoul + factor_lj*forcelj) * r2inv;
}
if (evflag) ev_tally(i,j,nlocal,newton_pair,
evdwl,ecoul,fpair,delx,dely,delz);
}
}
}
}
/* ----------------------------------------------------------------------
allocate all arrays
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::allocate()
{
allocated = 1;
int n = atom->ntypes;
memory->create(setflag,n+1,n+1,"pair:setflag");
for (int i = 1; i <= n; i++)
for (int j = i; j <= n; j++)
setflag[i][j] = 0;
memory->create(cutsq,n+1,n+1,"pair:cutsq");
memory->create(cut_lj,n+1,n+1,"pair:cut_lj");
memory->create(cut_ljsq,n+1,n+1,"pair:cut_ljsq");
memory->create(epsilon,n+1,n+1,"pair:epsilon");
memory->create(sigma,n+1,n+1,"pair:sigma");
memory->create(lj1,n+1,n+1,"pair:lj1");
memory->create(lj2,n+1,n+1,"pair:lj2");
memory->create(lj3,n+1,n+1,"pair:lj3");
memory->create(lj4,n+1,n+1,"pair:lj4");
memory->create(offset,n+1,n+1,"pair:offset");
memory->create(shift,n+1,n+1,"pair:shift");
}
/* ----------------------------------------------------------------------
global settings
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::settings(int narg, char **arg)
{
if (narg < 1 || narg > 2) error->all(FLERR,"Illegal pair_style command");
cut_lj_global = force->numeric(FLERR,arg[0]);
if (narg == 1) cut_coul = cut_lj_global;
else cut_coul = force->numeric(FLERR,arg[1]);
// reset cutoffs that have been explicitly set
if (allocated) {
int i,j;
for (i = 1; i <= atom->ntypes; i++)
for (j = i; j <= atom->ntypes; j++)
if (setflag[i][j]) cut_lj[i][j] = cut_lj_global;
}
}
/* ----------------------------------------------------------------------
set coeffs for one or more type pairs
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::coeff(int narg, char **arg)
{
if (narg < 5 || narg > 6)
error->all(FLERR,"Incorrect args for pair coefficients");
if (!allocated) allocate();
int ilo,ihi,jlo,jhi;
force->bounds(FLERR,arg[0],atom->ntypes,ilo,ihi);
force->bounds(FLERR,arg[1],atom->ntypes,jlo,jhi);
double epsilon_one = force->numeric(FLERR,arg[2]);
double sigma_one = force->numeric(FLERR,arg[3]);
double shift_one = force->numeric(FLERR,arg[4]);
double cut_lj_one = cut_lj_global;
if (narg == 6) cut_lj_one = force->numeric(FLERR,arg[5]);
int count = 0;
for (int i = ilo; i <= ihi; i++) {
for (int j = MAX(jlo,i); j <= jhi; j++) {
epsilon[i][j] = epsilon_one;
sigma[i][j] = sigma_one;
shift[i][j] = shift_one;
cut_lj[i][j] = cut_lj_one;
setflag[i][j] = 1;
count++;
}
}
if (count == 0) error->all(FLERR,"Incorrect args for pair coefficients");
}
/* ----------------------------------------------------------------------
init specific to this pair style
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::init_style()
{
if (!atom->q_flag)
error->all(FLERR,"Pair style lj/cut/coul/long requires atom attribute q");
// request regular or rRESPA neighbor lists
int irequest;
int respa = 0;
if (update->whichflag == 1 && strstr(update->integrate_style,"respa")) {
if (((Respa *) update->integrate)->level_inner >= 0) respa = 1;
if (((Respa *) update->integrate)->level_middle >= 0) respa = 2;
}
irequest = neighbor->request(this,instance_me);
if (respa >= 1) {
neighbor->requests[irequest]->respaouter = 1;
neighbor->requests[irequest]->respainner = 1;
}
if (respa == 2) neighbor->requests[irequest]->respamiddle = 1;
cut_coulsq = cut_coul * cut_coul;
// set rRESPA cutoffs
if (strstr(update->integrate_style,"respa") &&
((Respa *) update->integrate)->level_inner >= 0)
cut_respa = ((Respa *) update->integrate)->cutoff;
else cut_respa = NULL;
// insure use of KSpace long-range solver, set g_ewald
if (force->kspace == NULL)
error->all(FLERR,"Pair style requires a KSpace style");
g_ewald = force->kspace->g_ewald;
// setup force tables
if (ncoultablebits) init_tables(cut_coul,cut_respa);
}
/* ----------------------------------------------------------------------
init for one type pair i,j and corresponding j,i
------------------------------------------------------------------------- */
double PairLJExpandCoulLong::init_one(int i, int j)
{
if (setflag[i][j] == 0) {
epsilon[i][j] = mix_energy(epsilon[i][i],epsilon[j][j],
sigma[i][i],sigma[j][j]);
sigma[i][j] = mix_distance(sigma[i][i],sigma[j][j]);
cut_lj[i][j] = mix_distance(cut_lj[i][i],cut_lj[j][j]);
}
// include TIP4P qdist in full cutoff, qdist = 0.0 if not TIP4P
double cut = MAX(cut_lj[i][j],cut_coul+2.0*qdist);
cut_ljsq[i][j] = cut_lj[i][j] * cut_lj[i][j];
lj1[i][j] = 48.0 * epsilon[i][j] * pow(sigma[i][j],12.0);
lj2[i][j] = 24.0 * epsilon[i][j] * pow(sigma[i][j],6.0);
lj3[i][j] = 4.0 * epsilon[i][j] * pow(sigma[i][j],12.0);
lj4[i][j] = 4.0 * epsilon[i][j] * pow(sigma[i][j],6.0);
if (offset_flag && (cut_lj[i][j] > 0.0)) {
double ratio = sigma[i][j] / cut_lj[i][j];
offset[i][j] = 4.0 * epsilon[i][j] * (pow(ratio,12.0) - pow(ratio,6.0));
} else offset[i][j] = 0.0;
cut_ljsq[j][i] = cut_ljsq[i][j];
lj1[j][i] = lj1[i][j];
lj2[j][i] = lj2[i][j];
lj3[j][i] = lj3[i][j];
lj4[j][i] = lj4[i][j];
shift[j][i] = shift[i][j];
offset[j][i] = offset[i][j];
// check interior rRESPA cutoff
if (cut_respa && MIN(cut_lj[i][j],cut_coul) < cut_respa[3])
error->all(FLERR,"Pair cutoff < Respa interior cutoff");
// compute I,J contribution to long-range tail correction
// count total # of atoms of type I and J via Allreduce
if (tail_flag) {
int *type = atom->type;
int nlocal = atom->nlocal;
double count[2],all[2];
count[0] = count[1] = 0.0;
for (int k = 0; k < nlocal; k++) {
if (type[k] == i) count[0] += 1.0;
if (type[k] == j) count[1] += 1.0;
}
MPI_Allreduce(count,all,2,MPI_DOUBLE,MPI_SUM,world);
double sig2 = sigma[i][j]*sigma[i][j];
double sig6 = sig2*sig2*sig2;
double rc1 = cut_lj[i][j];
double rc2 = rc1*rc1;
double rc3 = rc2*rc1;
double rc6 = rc3*rc3;
double rc9 = rc3*rc6;
double shift1 = shift[i][j];
double shift2 = shift1*shift1;
double shift3 = shift2*shift1;
etail_ij = 8.0*MY_PI*all[0]*all[1]*epsilon[i][j] * sig6 *
((1.0/9.0 + 2.0*shift1/(10.0*rc1) + shift2/(11.0*rc2))*sig6/rc9 -
(1.0/3.0 + 2.0*shift1/(4.0*rc1) + shift2/(5.0*rc2))/rc3);
ptail_ij = 16.0*MY_PI*all[0]*all[1]*epsilon[i][j] * sig6 *
((1.0/9.0 + 3.0*shift1/(10.0*rc1) +
3.0*shift2/(11.0*rc2) + shift3/(12.0*rc3))*2.0*sig6/rc9 -
(1.0/3.0 + 3.0*shift1/(4.0*rc1) +
3.0*shift2/(5.0*rc2) + shift3/(6.0*rc3))/rc3);
}
return cut;
}
/* ----------------------------------------------------------------------
proc 0 writes to restart file
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::write_restart(FILE *fp)
{
write_restart_settings(fp);
int i,j;
for (i = 1; i <= atom->ntypes; i++)
for (j = i; j <= atom->ntypes; j++) {
fwrite(&setflag[i][j],sizeof(int),1,fp);
if (setflag[i][j]) {
fwrite(&epsilon[i][j],sizeof(double),1,fp);
fwrite(&sigma[i][j],sizeof(double),1,fp);
fwrite(&shift[i][j],sizeof(double),1,fp);
fwrite(&cut_lj[i][j],sizeof(double),1,fp);
}
}
}
/* ----------------------------------------------------------------------
proc 0 reads from restart file, bcasts
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::read_restart(FILE *fp)
{
read_restart_settings(fp);
allocate();
int i,j;
int me = comm->me;
for (i = 1; i <= atom->ntypes; i++)
for (j = i; j <= atom->ntypes; j++) {
if (me == 0) fread(&setflag[i][j],sizeof(int),1,fp);
MPI_Bcast(&setflag[i][j],1,MPI_INT,0,world);
if (setflag[i][j]) {
if (me == 0) {
fread(&epsilon[i][j],sizeof(double),1,fp);
fread(&sigma[i][j],sizeof(double),1,fp);
fread(&shift[i][j],sizeof(double),1,fp);
fread(&cut_lj[i][j],sizeof(double),1,fp);
}
MPI_Bcast(&epsilon[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&sigma[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&shift[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&cut_lj[i][j],1,MPI_DOUBLE,0,world);
}
}
}
/* ----------------------------------------------------------------------
proc 0 writes to restart file
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::write_restart_settings(FILE *fp)
{
fwrite(&cut_lj_global,sizeof(double),1,fp);
fwrite(&cut_coul,sizeof(double),1,fp);
fwrite(&offset_flag,sizeof(int),1,fp);
fwrite(&mix_flag,sizeof(int),1,fp);
fwrite(&tail_flag,sizeof(int),1,fp);
fwrite(&ncoultablebits,sizeof(int),1,fp);
fwrite(&tabinner,sizeof(double),1,fp);
}
/* ----------------------------------------------------------------------
proc 0 reads from restart file, bcasts
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::read_restart_settings(FILE *fp)
{
if (comm->me == 0) {
fread(&cut_lj_global,sizeof(double),1,fp);
fread(&cut_coul,sizeof(double),1,fp);
fread(&offset_flag,sizeof(int),1,fp);
fread(&mix_flag,sizeof(int),1,fp);
fread(&tail_flag,sizeof(int),1,fp);
fread(&ncoultablebits,sizeof(int),1,fp);
fread(&tabinner,sizeof(double),1,fp);
}
MPI_Bcast(&cut_lj_global,1,MPI_DOUBLE,0,world);
MPI_Bcast(&cut_coul,1,MPI_DOUBLE,0,world);
MPI_Bcast(&offset_flag,1,MPI_INT,0,world);
MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
MPI_Bcast(&tail_flag,1,MPI_INT,0,world);
MPI_Bcast(&ncoultablebits,1,MPI_INT,0,world);
MPI_Bcast(&tabinner,1,MPI_DOUBLE,0,world);
}
/* ----------------------------------------------------------------------
proc 0 writes to data file
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::write_data(FILE *fp)
{
for (int i = 1; i <= atom->ntypes; i++)
fprintf(fp,"%d %g %g %g\n",i,epsilon[i][i],sigma[i][i],shift[i][i]);
}
/* ----------------------------------------------------------------------
proc 0 writes all pairs to data file
------------------------------------------------------------------------- */
void PairLJExpandCoulLong::write_data_all(FILE *fp)
{
for (int i = 1; i <= atom->ntypes; i++)
for (int j = i; j <= atom->ntypes; j++)
fprintf(fp,"%d %d %g %g %g %g\n",i,j,epsilon[i][j],sigma[i][j],shift[i][j],cut_lj[i][j]);
}
/* ---------------------------------------------------------------------- */
double PairLJExpandCoulLong::single(int i, int j, int itype, int jtype,
double rsq,
double factor_coul, double factor_lj,
double &fforce)
{
double r2inv,r6inv,r,grij,expm2,t,erfc,prefactor;
double fraction,table,forcecoul,forcelj,phicoul,philj;
int itable;
r2inv = 1.0/rsq;
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq) {
r = sqrt(rsq);
grij = g_ewald * r;
expm2 = exp(-grij*grij);
t = 1.0 / (1.0 + EWALD_P*grij);
erfc = t * (A1+t*(A2+t*(A3+t*(A4+t*A5)))) * expm2;
prefactor = force->qqrd2e * atom->q[i]*atom->q[j]/r;
forcecoul = prefactor * (erfc + EWALD_F*grij*expm2);
if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*prefactor;
} else {
union_int_float_t rsq_lookup_single;
rsq_lookup_single.f = rsq;
itable = rsq_lookup_single.i & ncoulmask;
itable >>= ncoulshiftbits;
fraction = (rsq_lookup_single.f - rtable[itable]) * drtable[itable];
table = ftable[itable] + fraction*dftable[itable];
forcecoul = atom->q[i]*atom->q[j] * table;
if (factor_coul < 1.0) {
table = ctable[itable] + fraction*dctable[itable];
prefactor = atom->q[i]*atom->q[j] * table;
forcecoul -= (1.0-factor_coul)*prefactor;
}
}
} else forcecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
r6inv = r2inv*r2inv*r2inv;
forcelj = r6inv * (lj1[itype][jtype]*r6inv - lj2[itype][jtype]);
} else forcelj = 0.0;
fforce = (forcecoul + factor_lj*forcelj) * r2inv;
double eng = 0.0;
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq)
phicoul = prefactor*erfc;
else {
table = etable[itable] + fraction*detable[itable];
phicoul = atom->q[i]*atom->q[j] * table;
}
if (factor_coul < 1.0) phicoul -= (1.0-factor_coul)*prefactor;
eng += phicoul;
}
if (rsq < cut_ljsq[itype][jtype]) {
philj = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype]) -
offset[itype][jtype];
eng += factor_lj*philj;
}
return eng;
}
/* ---------------------------------------------------------------------- */
void *PairLJExpandCoulLong::extract(const char *str, int &dim)
{
dim = 0;
if (strcmp(str,"cut_coul") == 0) return (void *) &cut_coul;
dim = 2;
if (strcmp(str,"epsilon") == 0) return (void *) epsilon;
if (strcmp(str,"sigma") == 0) return (void *) sigma;
if (strcmp(str,"delta") == 0) return (void *) shift;
return NULL;
}

93
src/USER-MISC/pair_lj_expand_coul_long.h Normal file
View File

@ -0,0 +1,93 @@
/* -*- 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 PAIR_CLASS
PairStyle(lj/expand/coul/long,PairLJExpandCoulLong)
#else
#ifndef LMP_PAIR_LJ_EXPAND_COUL_LONG_H
#define LMP_PAIR_LJ_EXPAND_COUL_LONG_H
#include "pair.h"
namespace LAMMPS_NS {
class PairLJExpandCoulLong : public Pair {
public:
PairLJExpandCoulLong(class LAMMPS *);
virtual ~PairLJExpandCoulLong();
virtual void compute(int, int);
virtual void settings(int, char **);
void coeff(int, char **);
virtual void init_style();
virtual double init_one(int, int);
void write_restart(FILE *);
void read_restart(FILE *);
virtual void write_restart_settings(FILE *);
virtual void read_restart_settings(FILE *);
void write_data(FILE *);
void write_data_all(FILE *);
virtual double single(int, int, int, int, double, double, double, double &);
void compute_inner();
void compute_middle();
virtual void compute_outer(int, int);
virtual void *extract(const char *, int &);
protected:
double cut_lj_global;
double **cut_lj,**cut_ljsq;
double cut_coul,cut_coulsq;
double **epsilon,**sigma;
double **lj1,**lj2,**lj3,**lj4,**offset,**shift;
double *cut_respa;
double qdist; // TIP4P distance from O site to negative charge
double g_ewald;
virtual void allocate();
};
}
#endif
#endif
/* ERROR/WARNING messages:
E: Illegal ... command
Self-explanatory. Check the input script syntax and compare to the
documentation for the command. You can use -echo screen as a
command-line option when running LAMMPS to see the offending line.
E: Incorrect args for pair coefficients
Self-explanatory. Check the input script or data file.
E: Pair style lj/cut/coul/long requires atom attribute q
The atom style defined does not have this attribute.
E: Pair style requires a KSpace style
No kspace style is defined.
E: Pair cutoff < Respa interior cutoff
One or more pairwise cutoffs are too short to use with the specified
rRESPA cutoffs.
*/