From eeff0b8633236526b57c1811f25781cfec5b10e7 Mon Sep 17 00:00:00 2001 From: Anders Hafreager Date: Wed, 14 Jun 2017 10:24:16 +0200 Subject: [PATCH] Added vashishta GPU package for NVidia --- lib/gpu/Nvidia.makefile | 14 + lib/gpu/lal_vashishta.cpp | 283 +++++++++++++ lib/gpu/lal_vashishta.cu | 744 +++++++++++++++++++++++++++++++++ lib/gpu/lal_vashishta.h | 97 +++++ lib/gpu/lal_vashishta_ext.cpp | 139 ++++++ src/GPU/Install.sh | 2 + src/GPU/pair_vashishta_gpu.cpp | 258 ++++++++++++ src/GPU/pair_vashishta_gpu.h | 69 +++ 8 files changed, 1606 insertions(+) create mode 100644 lib/gpu/lal_vashishta.cpp create mode 100644 lib/gpu/lal_vashishta.cu create mode 100644 lib/gpu/lal_vashishta.h create mode 100644 lib/gpu/lal_vashishta_ext.cpp create mode 100644 src/GPU/pair_vashishta_gpu.cpp create mode 100644 src/GPU/pair_vashishta_gpu.h diff --git a/lib/gpu/Nvidia.makefile b/lib/gpu/Nvidia.makefile index 660544cfaa..ee2eb72632 100644 --- a/lib/gpu/Nvidia.makefile +++ b/lib/gpu/Nvidia.makefile @@ -63,6 +63,7 @@ OBJS = $(OBJ_DIR)/lal_atom.o $(OBJ_DIR)/lal_ans.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 \ @@ -117,6 +118,7 @@ CBNS = $(OBJ_DIR)/device.cubin $(OBJ_DIR)/device_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 \ @@ -613,6 +615,18 @@ $(OBJ_DIR)/lal_coul_dsf.o: $(ALL_H) lal_coul_dsf.h lal_coul_dsf.cpp $(OBJ_DIR)/c $(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) --cubin -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) --cubin -DNV_KERNEL -o $@ lal_sw.cu diff --git a/lib/gpu/lal_vashishta.cpp b/lib/gpu/lal_vashishta.cpp new file mode 100644 index 0000000000..96537e65d3 --- /dev/null +++ b/lib/gpu/lal_vashishta.cpp @@ -0,0 +1,283 @@ +/*************************************************************************** + vashishta.cpp + ------------------- + Anders Hafreager (UiO) + + Class for acceleration of the vashishta pair style. + + __________________________________________________________________________ + This file is part of the LAMMPS Accelerator Library (LAMMPS_AL) + __________________________________________________________________________ + + begin : Mon June 12, 2017 + email : andershaf@gmail.com + ***************************************************************************/ + +#if defined(USE_OPENCL) +#include "vashishta_cl.h" +#elif defined(USE_CUDART) +const char *vashishta=0; +#else +#include "vashishta_cubin.h" +#endif + +#include "lal_vashishta.h" +#include +using namespace LAMMPS_AL; +#define VashishtaT Vashishta + +extern Device device; + +template +VashishtaT::Vashishta() : BaseThree(), _allocated(false) { +} + +template +VashishtaT::~Vashishta() { + clear(); +} + +template +int VashishtaT::bytes_per_atom(const int max_nbors) const { + return this->bytes_per_atom_atomic(max_nbors); +} + +template +int VashishtaT::init(const int ntypes, const int nlocal, const int nall, const int max_nbors, + const double cell_size, const double gpu_split, FILE *_screen, + int* host_map, const int nelements, int*** host_elem2param, const int nparams, + const double* cutsq, const double* r0, + const double* gamma, const double* eta, + const double* lam1inv, const double* lam4inv, + const double* zizj, const double* mbigd, + const double* dvrc, const double* big6w, + const double* heta, const double* bigh, + const double* bigw, const double* c0, + const double* costheta, const double* bigb, + const double* big2b, const double* bigc) +{ + int success; + success=this->init_three(nlocal,nall,max_nbors,0,cell_size,gpu_split, + _screen,vashishta,"k_vashishta","k_vashishta_three_center", + "k_vashishta_three_end"); + 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; + + _nparams = nparams; + _nelements = nelements; + + UCL_H_Vec dview(nparams,*(this->ucl_device), + UCL_WRITE_ONLY); + + for (int i=0; iucl_device),UCL_READ_ONLY); + + for (int i=0; i(eta[i]); + dview[i].y=static_cast(lam1inv[i]); + dview[i].z=static_cast(lam4inv[i]); + dview[i].w=static_cast(zizj[i]); + } + + ucl_copy(param1,dview,false); + param1_tex.get_texture(*(this->pair_program),"param1_tex"); + param1_tex.bind_float(param1,4); + + param2.alloc(nparams,*(this->ucl_device),UCL_READ_ONLY); + + for (int i=0; i(mbigd[i]); + dview[i].y=static_cast(dvrc[i]); + dview[i].z=static_cast(big6w[i]); + dview[i].w=static_cast(heta[i]); + } + + ucl_copy(param2,dview,false); + param2_tex.get_texture(*(this->pair_program),"param2_tex"); + param2_tex.bind_float(param2,4); + + param3.alloc(nparams,*(this->ucl_device),UCL_READ_ONLY); + + for (int i=0; i(bigh[i]); + dview[i].y=static_cast(bigw[i]); + dview[i].z=static_cast(dvrc[i]); + dview[i].w=static_cast(c0[i]); + } + + ucl_copy(param3,dview,false); + param3_tex.get_texture(*(this->pair_program),"param3_tex"); + param3_tex.bind_float(param3,4); + + param4.alloc(nparams,*(this->ucl_device),UCL_READ_ONLY); + + for (int i=0; i(r0sq); + dview[i].y=static_cast(gamma[i]); + dview[i].z=static_cast(cutsq[i]); + dview[i].w=static_cast(r0[i]); + } + + ucl_copy(param4,dview,false); + param4_tex.get_texture(*(this->pair_program),"param4_tex"); + param4_tex.bind_float(param4,4); + + param5.alloc(nparams,*(this->ucl_device),UCL_READ_ONLY); + + for (int i=0; i(bigc[i]); + dview[i].y=static_cast(costheta[i]); + dview[i].z=static_cast(bigb[i]); + dview[i].w=static_cast(big2b[i]); + } + + ucl_copy(param5,dview,false); + param5_tex.get_texture(*(this->pair_program),"param5_tex"); + param5_tex.bind_float(param5,4); + + UCL_H_Vec dview_elem2param(nelements*nelements*nelements, + *(this->ucl_device), UCL_WRITE_ONLY); + + elem2param.alloc(nelements*nelements*nelements,*(this->ucl_device), + UCL_READ_ONLY); + + for (int i = 0; i < nelements; i++) + for (int j = 0; j < nelements; j++) + for (int k = 0; k < nelements; k++) { + int idx = i*nelements*nelements+j*nelements+k; + dview_elem2param[idx] = host_elem2param[i][j][k]; + } + + ucl_copy(elem2param,dview_elem2param,false); + + UCL_H_Vec dview_map(lj_types, *(this->ucl_device), UCL_WRITE_ONLY); + for (int i = 0; i < ntypes; i++) + dview_map[i] = host_map[i]; + + map.alloc(lj_types,*(this->ucl_device), UCL_READ_ONLY); + ucl_copy(map,dview_map,false); + + _allocated=true; + this->_max_bytes=param1.row_bytes()+param2.row_bytes()+param3.row_bytes()+param4.row_bytes()+param5.row_bytes()+ + map.row_bytes()+elem2param.row_bytes(); + return 0; +} + +template +void VashishtaT::clear() { + if (!_allocated) + return; + _allocated=false; + + param1.clear(); + param2.clear(); + param3.clear(); + param4.clear(); + param5.clear(); + map.clear(); + elem2param.clear(); + this->clear_atomic(); +} + +template +double VashishtaT::host_memory_usage() const { + return this->host_memory_usage_atomic()+sizeof(Vashishta); +} + +#define KTHREADS this->_threads_per_atom +#define JTHREADS this->_threads_per_atom +// --------------------------------------------------------------------------- +// Calculate energies, forces, and torques +// --------------------------------------------------------------------------- +template +void VashishtaT::loop(const bool _eflag, const bool _vflag, const int evatom) { + // Compute the block size and grid size to keep all cores busy + int BX=this->block_pair(); + int eflag, vflag; + if (_eflag) + eflag=1; + else + eflag=0; + + if (_vflag) + vflag=1; + else + vflag=0; + + int GX=static_cast(ceil(static_cast(this->ans->inum())/ + (BX/this->_threads_per_atom))); + + // this->_nbor_data == nbor->dev_packed for gpu_nbor == 0 and tpa > 1 + // this->_nbor_data == nbor->dev_nbor for gpu_nbor == 1 or tpa == 1 + int ainum=this->ans->inum(); + int nbor_pitch=this->nbor->nbor_pitch(); + this->time_pair.start(); + + this->k_pair.set_size(GX,BX); + this->k_pair.run(&this->atom->x, ¶m1, ¶m2, ¶m3, ¶m4, ¶m5, + &map, &elem2param, &_nelements, + &this->nbor->dev_nbor, &this->_nbor_data->begin(), + &this->ans->force, &this->ans->engv, + &eflag, &vflag, &ainum, &nbor_pitch, + &this->_threads_per_atom); + + BX=this->block_size(); + GX=static_cast(ceil(static_cast(this->ans->inum())/ + (BX/(KTHREADS*JTHREADS)))); + + this->k_three_center.set_size(GX,BX); + this->k_three_center.run(&this->atom->x, ¶m1, ¶m2, ¶m3, ¶m4, ¶m5, + &map, &elem2param, &_nelements, + &this->nbor->dev_nbor, &this->_nbor_data->begin(), + &this->ans->force, &this->ans->engv, &eflag, &vflag, &ainum, + &nbor_pitch, &this->_threads_per_atom, &evatom); + Answer *end_ans; + #ifdef THREE_CONCURRENT + end_ans=this->ans2; + #else + end_ans=this->ans; + #endif + if (evatom!=0) { + + this->k_three_end_vatom.set_size(GX,BX); + this->k_three_end_vatom.run(&this->atom->x, ¶m1, ¶m2, ¶m3, ¶m4, ¶m5, + &map, &elem2param, &_nelements, + &this->nbor->dev_nbor, &this->_nbor_data->begin(), + &this->nbor->dev_acc, + &end_ans->force, &end_ans->engv, &eflag, &vflag, &ainum, + &nbor_pitch, &this->_threads_per_atom, &this->_gpu_nbor); + } else { + + this->k_three_end.set_size(GX,BX); + this->k_three_end.run(&this->atom->x, ¶m1, ¶m2, ¶m3, ¶m4, ¶m5, + &map, &elem2param, &_nelements, + &this->nbor->dev_nbor, &this->_nbor_data->begin(), + &this->nbor->dev_acc, + &end_ans->force, &end_ans->engv, &eflag, &vflag, &ainum, + &nbor_pitch, &this->_threads_per_atom, &this->_gpu_nbor); + } + + this->time_pair.stop(); +} + +template class Vashishta; + diff --git a/lib/gpu/lal_vashishta.cu b/lib/gpu/lal_vashishta.cu new file mode 100644 index 0000000000..9abbbd10e5 --- /dev/null +++ b/lib/gpu/lal_vashishta.cu @@ -0,0 +1,744 @@ +// ************************************************************************** +// vashishta.cu +// ------------------- +// Anders Hafreager (UiO) +// +// Device code for acceleration of the vashishta pair style +// +// __________________________________________________________________________ +// This file is part of the LAMMPS Accelerator Library (LAMMPS_AL) +// __________________________________________________________________________ +// +// begin : Mon June 12, 2017 +// email : andershaf@gmail.com +// ***************************************************************************/ + +#ifdef NV_KERNEL +#include "lal_aux_fun1.h" + +#ifndef _DOUBLE_DOUBLE +texture pos_tex; +texture param1_tex; +texture param2_tex; +texture param3_tex; +texture param4_tex; +texture param5_tex; +#else +texture pos_tex; +texture param1_tex; +texture param2_tex; +texture param3_tex; +texture param4_tex; +texture param5_tex; +#endif + +#else +#define pos_tex x_ +#define param1_tex param1 +#define param2_tex param2 +#define param3_tex param3 +#define param3_tex param4 +#define param3_tex param5 +#endif + +#define THIRD (numtyp)0.66666666666666666667 + +//#define THREE_CONCURRENT + +#if (ARCH < 300) + +#define store_answers_p(f, energy, virial, ii, inum, tid, t_per_atom, offset, \ + eflag, vflag, ans, engv) \ + if (t_per_atom>1) { \ + __local acctyp red_acc[6][BLOCK_ELLIPSE]; \ + red_acc[0][tid]=f.x; \ + red_acc[1][tid]=f.y; \ + red_acc[2][tid]=f.z; \ + red_acc[3][tid]=energy; \ + for (unsigned int s=t_per_atom/2; s>0; s>>=1) { \ + if (offset < s) { \ + for (int r=0; r<4; 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]; \ + energy=red_acc[3][tid]; \ + if (vflag>0) { \ + for (int r=0; r<6; r++) \ + red_acc[r][tid]=virial[r]; \ + 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]; \ + } \ + } \ + for (int r=0; r<6; r++) \ + virial[r]=red_acc[r][tid]; \ + } \ + } \ + if (offset==0) { \ + int ei=ii; \ + if (eflag>0) { \ + engv[ei]+=energy*(acctyp)0.5; \ + ei+=inum; \ + } \ + if (vflag>0) { \ + for (int i=0; i<6; i++) { \ + engv[ei]+=virial[i]*(acctyp)0.5; \ + ei+=inum; \ + } \ + } \ + acctyp4 old=ans[ii]; \ + old.x+=f.x; \ + old.y+=f.y; \ + old.z+=f.z; \ + ans[ii]=old; \ + } + +#else + +#define store_answers_p(f, energy, 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); \ + energy += shfl_xor(energy, 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; \ + } \ + if (vflag>0) { \ + for (int i=0; i<6; i++) { \ + engv[ei]+=virial[i]*(acctyp)0.5; \ + ei+=inum; \ + } \ + } \ + acctyp4 old=ans[ii]; \ + old.x+=f.x; \ + old.y+=f.y; \ + old.z+=f.z; \ + ans[ii]=old; \ + } + +#endif + + +__kernel void k_vashishta(const __global numtyp4 *restrict x_, + const __global numtyp4 *restrict param1, + const __global numtyp4 *restrict param2, + const __global numtyp4 *restrict param3, + const __global numtyp4 *restrict param4, + const __global numtyp4 *restrict param5, + const __global int *restrict map, + const __global int *restrict elem2param, + const int nelements, + 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 int t_per_atom) { + __local int n_stride; + int tid, ii, offset; + atom_info(t_per_atom,ii,tid,offset); + + acctyp energy=(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 (ii0) + energy += (param3_bigh*reta+vc2-vc3-param3_bigw*r6inv-r*param3_dvrc+param3_c0); + + 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(f,energy,virial,ii,inum,tid,t_per_atom,offset,eflag,vflag, + ans,engv); + } // if ii + +} + +#define threebody(delr1x, delr1y, delr1z, eflag, energy) \ +{ \ + numtyp r1 = ucl_sqrt(rsq1); \ + numtyp rinvsq1 = ucl_recip(rsq1); \ + numtyp rainv1 = ucl_recip(r1 - param_r0_ij); \ + numtyp gsrainv1 = param_gamma_ij * rainv1; \ + numtyp gsrainvsq1 = gsrainv1*rainv1/r1; \ + numtyp expgsrainv1 = ucl_exp(gsrainv1); \ + \ + numtyp r2 = ucl_sqrt(rsq2); \ + numtyp rinvsq2 = ucl_recip(rsq2); \ + numtyp rainv2 = ucl_recip(r2 - param_r0_ik); \ + numtyp gsrainv2 = param_gamma_ik * rainv2; \ + numtyp gsrainvsq2 = gsrainv2*rainv2/r2; \ + numtyp expgsrainv2 = ucl_exp(gsrainv2); \ + \ + numtyp rinv12 = ucl_recip(r1*r2); \ + numtyp cs = (delr1x*delr2x + delr1y*delr2y + delr1z*delr2z) * rinv12; \ + numtyp delcs = cs - param_costheta_ijk; \ + numtyp delcssq = delcs*delcs; \ + numtyp pcsinv = param_bigc_ijk*delcssq+1.0; \ + numtyp pcsinvsq = pcsinv*pcsinv; \ + numtyp pcs = delcssq/pcsinv; \ + \ + numtyp facexp = expgsrainv1*expgsrainv2; \ + \ + numtyp facrad = param_bigb_ijk * facexp*pcs; \ + numtyp frad1 = facrad*gsrainvsq1; \ + numtyp frad2 = facrad*gsrainvsq2; \ + numtyp facang = param_big2b_ijk * facexp*delcs/pcsinvsq; \ + numtyp facang12 = rinv12*facang; \ + numtyp csfacang = cs*facang; \ + numtyp csfac1 = rinvsq1*csfacang; \ + \ + fjx = delr1x*(frad1+csfac1)-delr2x*facang12; \ + fjy = delr1y*(frad1+csfac1)-delr2y*facang12; \ + fjz = delr1z*(frad1+csfac1)-delr2z*facang12; \ + \ + numtyp csfac2 = rinvsq2*csfacang; \ + \ + fkx = delr2x*(frad2+csfac2)-delr1x*facang12; \ + fky = delr2y*(frad2+csfac2)-delr1y*facang12; \ + fkz = delr2z*(frad2+csfac2)-delr1z*facang12; \ + \ + if (eflag>0) \ + energy+=facrad; \ + if (vflag>0) { \ + virial[0] += delr1x*fjx + delr2x*fkx; \ + virial[1] += delr1y*fjy + delr2y*fky; \ + virial[2] += delr1z*fjz + delr2z*fkz; \ + virial[3] += delr1x*fjy + delr2x*fky; \ + virial[4] += delr1x*fjz + delr2x*fkz; \ + virial[5] += delr1y*fjz + delr2y*fkz; \ + } \ +} + +#define threebody_half(delr1x, delr1y, delr1z) \ +{ \ + numtyp r1 = ucl_sqrt(rsq1); \ + numtyp rinvsq1 = ucl_recip(rsq1); \ + numtyp rainv1 = ucl_recip(r1 - param_r0_ij); \ + numtyp gsrainv1 = param_gamma_ij * rainv1; \ + numtyp gsrainvsq1 = gsrainv1*rainv1/r1; \ + numtyp expgsrainv1 = ucl_exp(gsrainv1); \ + \ + numtyp r2 = ucl_sqrt(rsq2); \ + numtyp rainv2 = ucl_recip(r2 - param_r0_ik); \ + numtyp gsrainv2 = param_gamma_ik * rainv2; \ + numtyp expgsrainv2 = ucl_exp(gsrainv2); \ + \ + numtyp rinv12 = ucl_recip(r1*r2); \ + numtyp cs = (delr1x*delr2x + delr1y*delr2y + delr1z*delr2z) * rinv12; \ + numtyp delcs = cs - param_costheta_ijk; \ + numtyp delcssq = delcs*delcs; \ + numtyp pcsinv = param_bigc_ijk*delcssq+1.0; \ + numtyp pcsinvsq = pcsinv*pcsinv; \ + numtyp pcs = delcssq/pcsinv; \ + \ + numtyp facexp = expgsrainv1*expgsrainv2; \ + \ + numtyp facrad = param_bigb_ijk * facexp*pcs; \ + numtyp frad1 = facrad*gsrainvsq1; \ + numtyp facang = param_big2b_ijk * facexp*delcs/pcsinvsq; \ + numtyp facang12 = rinv12*facang; \ + numtyp csfacang = cs*facang; \ + numtyp csfac1 = rinvsq1*csfacang; \ + \ + fjx = delr1x*(frad1+csfac1)-delr2x*facang12; \ + fjy = delr1y*(frad1+csfac1)-delr2y*facang12; \ + fjz = delr1z*(frad1+csfac1)-delr2z*facang12; \ +} + +__kernel void k_vashishta_three_center(const __global numtyp4 *restrict x_, + const __global numtyp4 *restrict param1, + const __global numtyp4 *restrict param2, + const __global numtyp4 *restrict param3, + const __global numtyp4 *restrict param4, + const __global numtyp4 *restrict param5, + const __global int *restrict map, + const __global int *restrict elem2param, + const int nelements, + 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 int t_per_atom, const int evatom) { + __local int tpa_sq, n_stride; + tpa_sq=fast_mul(t_per_atom,t_per_atom); + numtyp param_gamma_ij, param_r0sq_ij, param_r0_ij, param_gamma_ik, param_r0sq_ik, param_r0_ik; + numtyp param_costheta_ijk, param_bigc_ijk, param_bigb_ijk, param_big2b_ijk; + + int tid, ii, offset; + atom_info(tpa_sq,ii,tid,offset); + + acctyp energy=(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 param_r0sq_ij) continue; + param_gamma_ij=param4_ijparam.y; + param_r0_ij=param4_ijparam.w; + + int nbor_k=nbor_j-offset_j+offset_k; + if (nbor_k<=nbor_j) + nbor_k+=n_stride; + + for ( ; nbor_k param_r0sq_ij) continue; + + param_gamma_ij=param4_ijparam.y; + param_r0_ij = param4_ijparam.w; + + int nbor_k,numk; + if (dev_nbor==dev_packed) { + if (gpu_nbor) nbor_k=j+nbor_pitch; + else nbor_k=dev_acc[j]+nbor_pitch; + numk=dev_nbor[nbor_k]; + nbor_k+=nbor_pitch+fast_mul(j,t_per_atom-1); + k_end=nbor_k+fast_mul(numk/t_per_atom,n_stride)+(numk & (t_per_atom-1)); + nbor_k+=offset_k; + } else { + nbor_k=dev_acc[j]+nbor_pitch; + numk=dev_nbor[nbor_k]; + nbor_k+=nbor_pitch; + nbor_k=dev_nbor[nbor_k]; + k_end=nbor_k+numk; + nbor_k+=offset_k; + } + + for ( ; nbor_k param_r0sq_ij) continue; + + param_gamma_ij=param4_ijparam.y; + param_r0_ij=param4_ijparam.w; + + int nbor_k,numk; + if (dev_nbor==dev_packed) { + if (gpu_nbor) nbor_k=j+nbor_pitch; + else nbor_k=dev_acc[j]+nbor_pitch; + numk=dev_nbor[nbor_k]; + nbor_k+=nbor_pitch+fast_mul(j,t_per_atom-1); + k_end=nbor_k+fast_mul(numk/t_per_atom,n_stride)+(numk & (t_per_atom-1)); + nbor_k+=offset_k; + } else { + nbor_k=dev_acc[j]+nbor_pitch; + numk=dev_nbor[nbor_k]; + nbor_k+=nbor_pitch; + nbor_k=dev_nbor[nbor_k]; + k_end=nbor_k+numk; + nbor_k+=offset_k; + } + + for ( ; nbor_k +class Vashishta : public BaseThree { + public: + Vashishta(); + ~Vashishta(); + + /// 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, const int nlocal, const int nall, const int max_nbors, + const double cell_size, const double gpu_split, FILE *screen, + int* host_map, const int nelements, int*** host_elem2param, const int nparams, + const double* cutsq, const double* r0, + const double* gamma, const double* eta, + const double* lam1inv, const double* lam4inv, + const double* zizj, const double* mbigd, + const double* dvrc, const double* big6w, + const double* heta, const double* bigh, + const double* bigw, const double* c0, + const double* costheta, const double* bigb, + const double* big2b, const double* bigc); + + /// 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 -------------------------- + + /// If atom type constants fit in shared memory, use fast kernels + bool shared_types; + + /// Number of atom types + int _lj_types; + + /// param1.x = eta, param1.y = lam1inv, param1.z = lam4inv, param1.w = zizj + UCL_D_Vec param1; + /// param2.x = mbigd, param2.y = dvrc, param2.z = big6w, param2.w = heta + UCL_D_Vec param2; + /// param3.x = bigh, param3.y = bigw, param3.z = dvrc, param3.w = c0 + UCL_D_Vec param3; + /// param4.x = r0sq, param4.y = gamma, param4.z = cutsq, param4.w = r0 + UCL_D_Vec param4; + /// param5.x = bigc, param5.y = costheta, param5.z = bigb, param5.w = big2b + UCL_D_Vec param5; + + UCL_D_Vec elem2param; + UCL_D_Vec map; + int _nparams,_nelements; + + UCL_Texture param1_tex, param2_tex, param3_tex, param4_tex, param5_tex; + + private: + bool _allocated; + void loop(const bool _eflag, const bool _vflag, const int evatom); + +}; + +} + +#endif + diff --git a/lib/gpu/lal_vashishta_ext.cpp b/lib/gpu/lal_vashishta_ext.cpp new file mode 100644 index 0000000000..fd4ae02525 --- /dev/null +++ b/lib/gpu/lal_vashishta_ext.cpp @@ -0,0 +1,139 @@ +/*************************************************************************** + vashishta_ext.cpp + ------------------- + Anders Hafreager (UiO) + + Class for acceleration of the vashishta pair style. + + __________________________________________________________________________ + This file is part of the LAMMPS Accelerator Library (LAMMPS_AL) + __________________________________________________________________________ + + begin : Mon June 12, 2017 + email : andershaf@gmail.com + ***************************************************************************/ + +#include +#include +#include + +#include "lal_vashishta.h" +#include +using namespace std; + +using namespace std; +using namespace LAMMPS_AL; + +static Vashishta VashishtaMF; + +// --------------------------------------------------------------------------- +// Allocate memory on host and device and copy constants to device +// --------------------------------------------------------------------------- +int vashishta_gpu_init(const int ntypes, const int inum, const int nall, const int max_nbors, + const double cell_size, int &gpu_mode, FILE *screen, + int* host_map, const int nelements, int*** host_elem2param, const int nparams, + const double* cutsq, const double* r0, + const double* gamma, const double* eta, + const double* lam1inv, const double* lam4inv, + const double* zizj, const double* mbigd, + const double* dvrc, const double* big6w, + const double* heta, const double* bigh, + const double* bigw, const double* c0, + const double* costheta, const double* bigb, + const double* big2b, const double* bigc) { + VashishtaMF.clear(); + gpu_mode=VashishtaMF.device->gpu_mode(); + double gpu_split=VashishtaMF.device->particle_split(); + int first_gpu=VashishtaMF.device->first_device(); + int last_gpu=VashishtaMF.device->last_device(); + int world_me=VashishtaMF.device->world_me(); + int gpu_rank=VashishtaMF.device->gpu_rank(); + int procs_per_gpu=VashishtaMF.device->procs_per_gpu(); + + // disable host/device split for now + if (gpu_split != 1.0) + return -8; + + VashishtaMF.device->init_message(screen,"vashishta/gpu",first_gpu,last_gpu); + + bool message=false; + if (VashishtaMF.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=VashishtaMF.init(ntypes, inum, nall, 500, cell_size, gpu_split, screen, + host_map, nelements, host_elem2param, nparams, + cutsq, r0, gamma, eta, lam1inv, + lam4inv, zizj, mbigd, dvrc, big6w, heta, bigh, bigw, + c0, costheta, bigb, big2b, bigc); + + VashishtaMF.device->world_barrier(); + if (message) + fprintf(screen,"Done.\n"); + + for (int i=0; igpu_barrier(); + if (message) + fprintf(screen,"Done.\n"); + } + if (message) + fprintf(screen,"\n"); + + if (init_ok==0) + VashishtaMF.estimate_gpu_overhead(); + + cout << "Seems like this was ok!" << endl; + return init_ok; +} + +void vashishta_gpu_clear() { + VashishtaMF.clear(); +} + +int ** vashishta_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) { + return VashishtaMF.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); +} + +void vashishta_gpu_compute(const int ago, const int nlocal, const int nall, + const int nlist, 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) { + VashishtaMF.compute(ago,nlocal,nall,nlist,host_x,host_type,ilist,numj, + firstneigh,eflag,vflag,eatom,vatom,host_start,cpu_time,success); +} + +double vashishta_gpu_bytes() { + return VashishtaMF.host_memory_usage(); +} + + diff --git a/src/GPU/Install.sh b/src/GPU/Install.sh index ca73adbf82..f4aeaa2706 100644 --- a/src/GPU/Install.sh +++ b/src/GPU/Install.sh @@ -113,6 +113,8 @@ action pair_soft_gpu.cpp action pair_soft_gpu.h action pair_sw_gpu.cpp pair_sw.cpp action pair_sw_gpu.h pair_sw.h +action pair_vashishta_gpu.cpp pair_vashishta.cpp +action pair_vashishta_gpu.h pair_vashishta.h action pair_table_gpu.cpp pair_table.cpp action pair_table_gpu.h pair_table.cpp action pair_tersoff_gpu.cpp pair_tersoff.cpp diff --git a/src/GPU/pair_vashishta_gpu.cpp b/src/GPU/pair_vashishta_gpu.cpp new file mode 100644 index 0000000000..19e0799671 --- /dev/null +++ b/src/GPU/pair_vashishta_gpu.cpp @@ -0,0 +1,258 @@ +/* ---------------------------------------------------------------------- + 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: Anders Hafreager (UiO) +------------------------------------------------------------------------- */ +#include +#include +#include +#include +#include +#include "pair_vashishta_gpu.h" +#include "atom.h" +#include "neighbor.h" +#include "neigh_request.h" +#include "force.h" +#include "comm.h" +#include "memory.h" +#include "neighbor.h" +#include "neigh_list.h" +#include "memory.h" +#include "error.h" +#include "domain.h" +#include "gpu_extra.h" + +using namespace LAMMPS_NS; + +// External functions from cuda library for atom decomposition + +int vashishta_gpu_init(const int ntypes, const int inum, const int nall, const int max_nbors, + const double cell_size, int &gpu_mode, FILE *screen, + int* host_map, const int nelements, int*** host_elem2param, const int nparams, + const double* cutsq, const double* r0, + const double* gamma, const double* eta, + const double* lam1inv, const double* lam4inv, + const double* zizj, const double* mbigd, + const double* dvrc, const double* big6w, + const double* heta, const double* bigh, + const double* bigw, const double* c0, + const double* costheta, const double* bigb, + const double* big2b, const double* bigc); +void vashishta_gpu_clear(); +int ** vashishta_gpu_compute_n(const int ago, const int inum, + const int nall, double **host_x, int *host_type, + double *sublo, double *subhi, tagint *tag, int **nspecial, + tagint **special, const bool eflag, const bool vflag, + const bool eatom, const bool vatom, int &host_start, + int **ilist, int **jnum, + const double cpu_time, bool &success); +void vashishta_gpu_compute(const int ago, const int nloc, const int nall, const int ln, + 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 vashishta_gpu_bytes(); +extern double lmp_gpu_forces(double **f, double **tor, double *eatom, + double **vatom, double *virial, double &ecoul); + +/* ---------------------------------------------------------------------- */ + +PairVashishtaGPU::PairVashishtaGPU(LAMMPS *lmp) : PairVashishta(lmp), gpu_mode(GPU_FORCE) +{ + cpu_time = 0.0; + reinitflag = 0; + gpu_allocated = false; + GPU_EXTRA::gpu_ready(lmp->modify, lmp->error); + + cutghost = NULL; + ghostneigh = 1; +} + +/* ---------------------------------------------------------------------- + check if allocated, since class can be destructed when incomplete +------------------------------------------------------------------------- */ + +PairVashishtaGPU::~PairVashishtaGPU() +{ + vashishta_gpu_clear(); + if (allocated) + memory->destroy(cutghost); +} + +/* ---------------------------------------------------------------------- */ + +void PairVashishtaGPU::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 = vashishta_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); + } else { + inum = list->inum; + ilist = list->ilist; + numneigh = list->numneigh; + firstneigh = list->firstneigh; + + vashishta_gpu_compute(neighbor->ago, inum, nall, inum+list->gnum, + atom->x, atom->type, ilist, numneigh, firstneigh, eflag, + vflag, eflag_atom, vflag_atom, host_start, cpu_time, + success); + } + if (!success) + error->one(FLERR,"Insufficient memory on accelerator"); +} + +/* ---------------------------------------------------------------------- */ + +void PairVashishtaGPU::allocate() +{ + if(!allocated) { + PairVashishta::allocate(); + } + int n = atom->ntypes; + + memory->create(cutghost,n+1,n+1,"pair:cutghost"); + gpu_allocated = true; +} + +/* ---------------------------------------------------------------------- + init specific to this pair style +------------------------------------------------------------------------- */ + +void PairVashishtaGPU::init_style() +{ + double cell_size = cutmax + neighbor->skin; + + if (atom->tag_enable == 0) + error->all(FLERR,"Pair style vashishta/gpu requires atom IDs"); + if (force->newton_pair != 0) + error->all(FLERR,"Pair style vashishta/gpu requires newton pair off"); + + double *cutsq, *r0, *r0eps, *gamma, *eta; + double *lam1inv, *lam4inv, *zizj, *mbigd; + double *dvrc, *big6w, *heta, *bigh; + double *bigw, *c0, *costheta, *bigb; + double *big2b, *bigc; + + cutsq = r0 = gamma = eta = NULL; + lam1inv = lam4inv = zizj = mbigd = NULL; + dvrc = big6w = heta = bigh = NULL; + bigw = c0 = costheta = bigb = NULL; + big2b = bigc = NULL; + + memory->create(cutsq,nparams,"pair:cutsq"); + memory->create(r0,nparams,"pair:r0"); + memory->create(gamma,nparams,"pair:gamma"); + memory->create(eta,nparams,"pair:eta"); + memory->create(lam1inv,nparams,"pair:lam1inv"); + memory->create(lam4inv,nparams,"pair:lam4inv"); + memory->create(zizj,nparams,"pair:zizj"); + memory->create(mbigd,nparams,"pair:mbigd"); + memory->create(dvrc,nparams,"pair:dvrc"); + memory->create(big6w,nparams,"pair:big6w"); + memory->create(heta,nparams,"pair:heta"); + memory->create(bigh,nparams,"pair:bigh"); + memory->create(bigw,nparams,"pair:bigw"); + memory->create(c0,nparams,"pair:c0"); + memory->create(costheta,nparams,"pair:costheta"); + memory->create(bigb,nparams,"pair:bigb"); + memory->create(big2b,nparams,"pair:big2b"); + memory->create(bigc,nparams,"pair:bigc"); + + for (int i = 0; i < nparams; i++) { + cutsq[i] = params[i].cutsq; + r0[i] = params[i].r0; + gamma[i] = params[i].gamma; + eta[i] = params[i].eta; + lam1inv[i] = params[i].lam1inv; + lam4inv[i] = params[i].lam4inv; + zizj[i] = params[i].zizj; + mbigd[i] = params[i].mbigd; + dvrc[i] = params[i].dvrc; + big6w[i] = params[i].big6w; + heta[i] = params[i].heta; + bigh[i] = params[i].bigh; + bigw[i] = params[i].bigw; + c0[i] = params[i].c0; + costheta[i] = params[i].costheta; + bigb[i] = params[i].bigb; + big2b[i] = params[i].big2b; + bigc[i] = params[i].bigc; + } + int success = vashishta_gpu_init(atom->ntypes+1, atom->nlocal, atom->nlocal+atom->nghost, 500, + cell_size, gpu_mode, screen, map, nelements, + elem2param, nparams, cutsq, r0, gamma, eta, lam1inv, + lam4inv, zizj, mbigd, dvrc, big6w, heta, bigh, bigw, + c0, costheta, bigb, big2b, bigc); + memory->destroy(cutsq); + memory->destroy(r0); + memory->destroy(gamma); + memory->destroy(eta); + memory->destroy(lam1inv); + memory->destroy(lam4inv); + memory->destroy(zizj); + memory->destroy(mbigd); + memory->destroy(dvrc); + memory->destroy(big6w); + memory->destroy(heta); + memory->destroy(bigh); + memory->destroy(bigw); + memory->destroy(c0); + memory->destroy(costheta); + memory->destroy(bigb); + memory->destroy(big2b); + memory->destroy(bigc); + + 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; + neighbor->requests[irequest]->ghost = 1; + } + + if (comm->cutghostuser < (2.0*cutmax + neighbor->skin) ) + comm->cutghostuser=2.0*cutmax + neighbor->skin; + +} + +/* ---------------------------------------------------------------------- + init for one type pair i,j and corresponding j,i +------------------------------------------------------------------------- */ + +double PairVashishtaGPU::init_one(int i, int j) +{ + if(!gpu_allocated) { + allocate(); + } + if (setflag[i][j] == 0) error->all(FLERR,"All pair coeffs are not set"); + cutghost[i][j] = cutmax; + cutghost[j][i] = cutmax; + + return cutmax; +} + diff --git a/src/GPU/pair_vashishta_gpu.h b/src/GPU/pair_vashishta_gpu.h new file mode 100644 index 0000000000..d54ede505b --- /dev/null +++ b/src/GPU/pair_vashishta_gpu.h @@ -0,0 +1,69 @@ +/* -*- 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(vashishta/gpu,PairVashishtaGPU) + +#else + +#ifndef LMP_PAIR_VASHISHTA_GPU_H +#define LMP_PAIR_VASHISHTA_GPU_H + +#include "pair_vashishta.h" + +namespace LAMMPS_NS { + +class PairVashishtaGPU : public PairVashishta { + public: + PairVashishtaGPU(class LAMMPS *); + ~PairVashishtaGPU(); + void compute(int, int); + double init_one(int, int); + void init_style(); + + enum { GPU_FORCE, GPU_NEIGH, GPU_HYB_NEIGH }; + + protected: + void allocate(); + int gpu_allocated; + 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 vashishta/gpu requires atom IDs + +This is a requirement to use this potential. + +E: Pair style vashishta/gpu requires newton pair off + +See the newton command. This is a restriction to use this potential. + +E: All pair coeffs are not set + +All pair coefficients must be set in the data file or by the +pair_coeff command before running a simulation. + +*/