1 // clang-format off
2 /* ----------------------------------------------------------------------
3 LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
4 https://www.lammps.org/, Sandia National Laboratories
5 Steve Plimpton, sjplimp@sandia.gov
6
7 Copyright (2003) Sandia Corporation. Under the terms of Contract
8 DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
9 certain rights in this software. This software is distributed under
10 the GNU General Public License.
11
12 See the README file in the top-level LAMMPS directory.
13 ------------------------------------------------------------------------- */
14
15 #include "pair_lj_sdk_kokkos.h"
16 #include <cmath>
17 #include <cstdio>
18 #include <cstring>
19 #include "kokkos.h"
20 #include "atom_kokkos.h"
21 #include "comm.h"
22 #include "force.h"
23 #include "neighbor.h"
24 #include "neigh_request.h"
25 #include "update.h"
26 #include "respa.h"
27 #include "memory_kokkos.h"
28 #include "error.h"
29 #include "atom_masks.h"
30
31 #include "lj_sdk_common.h"
32
33 using namespace LAMMPS_NS;
34 using namespace LJSDKParms;
35
36 #define KOKKOS_CUDA_MAX_THREADS 256
37 #define KOKKOS_CUDA_MIN_BLOCKS 8
38
39 /* ---------------------------------------------------------------------- */
40
41 template<class DeviceType>
PairLJSDKKokkos(LAMMPS * lmp)42 PairLJSDKKokkos<DeviceType>::PairLJSDKKokkos(LAMMPS *lmp) : PairLJSDK(lmp)
43 {
44 respa_enable = 0;
45
46 kokkosable = 1;
47 atomKK = (AtomKokkos *) atom;
48 execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
49 datamask_read = X_MASK | F_MASK | TYPE_MASK | ENERGY_MASK | VIRIAL_MASK;
50 datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
51 }
52
53 /* ---------------------------------------------------------------------- */
54
55 template<class DeviceType>
~PairLJSDKKokkos()56 PairLJSDKKokkos<DeviceType>::~PairLJSDKKokkos()
57 {
58 if (copymode) return;
59
60 if (allocated) {
61 memoryKK->destroy_kokkos(k_eatom,eatom);
62 memoryKK->destroy_kokkos(k_vatom,vatom);
63 memoryKK->destroy_kokkos(k_cutsq,cutsq);
64 }
65 }
66
67 /* ---------------------------------------------------------------------- */
68
69 template<class DeviceType>
compute(int eflag_in,int vflag_in)70 void PairLJSDKKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
71 {
72 eflag = eflag_in;
73 vflag = vflag_in;
74
75
76 if (neighflag == FULL) no_virial_fdotr_compute = 1;
77
78 ev_init(eflag,vflag,0);
79
80 // reallocate per-atom arrays if necessary
81
82 if (eflag_atom) {
83 memoryKK->destroy_kokkos(k_eatom,eatom);
84 memoryKK->create_kokkos(k_eatom,eatom,maxeatom,"pair:eatom");
85 d_eatom = k_eatom.view<DeviceType>();
86 }
87 if (vflag_atom) {
88 memoryKK->destroy_kokkos(k_vatom,vatom);
89 memoryKK->create_kokkos(k_vatom,vatom,maxvatom,"pair:vatom");
90 d_vatom = k_vatom.view<DeviceType>();
91 }
92
93 atomKK->sync(execution_space,datamask_read);
94 k_cutsq.template sync<DeviceType>();
95 k_params.template sync<DeviceType>();
96 if (eflag || vflag) atomKK->modified(execution_space,datamask_modify);
97 else atomKK->modified(execution_space,F_MASK);
98
99 x = atomKK->k_x.view<DeviceType>();
100 c_x = atomKK->k_x.view<DeviceType>();
101 f = atomKK->k_f.view<DeviceType>();
102 type = atomKK->k_type.view<DeviceType>();
103 tag = atomKK->k_tag.view<DeviceType>();
104 nlocal = atom->nlocal;
105 nall = atom->nlocal + atom->nghost;
106 newton_pair = force->newton_pair;
107 special_lj[0] = force->special_lj[0];
108 special_lj[1] = force->special_lj[1];
109 special_lj[2] = force->special_lj[2];
110 special_lj[3] = force->special_lj[3];
111
112 // loop over neighbors of my atoms
113
114 EV_FLOAT ev = pair_compute<PairLJSDKKokkos<DeviceType>,void >(this,(NeighListKokkos<DeviceType>*)list);
115
116 if (eflag) eng_vdwl += ev.evdwl;
117 if (vflag_global) {
118 virial[0] += ev.v[0];
119 virial[1] += ev.v[1];
120 virial[2] += ev.v[2];
121 virial[3] += ev.v[3];
122 virial[4] += ev.v[4];
123 virial[5] += ev.v[5];
124 }
125
126 if (eflag_atom) {
127 k_eatom.template modify<DeviceType>();
128 k_eatom.template sync<LMPHostType>();
129 }
130
131 if (vflag_atom) {
132 k_vatom.template modify<DeviceType>();
133 k_vatom.template sync<LMPHostType>();
134 }
135
136 if (vflag_fdotr) pair_virial_fdotr_compute(this);
137
138 }
139
140 template<class DeviceType>
141 template<bool STACKPARAMS, class Specialisation>
142 KOKKOS_INLINE_FUNCTION
143 F_FLOAT PairLJSDKKokkos<DeviceType>::
compute_fpair(const F_FLOAT & rsq,const int & i,const int & j,const int & itype,const int & jtype) const144 compute_fpair(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
145 (void) i;
146 (void) j;
147 const F_FLOAT r2inv = 1.0/rsq;
148 const int ljt = (STACKPARAMS?m_params[itype][jtype].lj_type:params(itype,jtype).lj_type);
149
150 const F_FLOAT lj_1 = (STACKPARAMS?m_params[itype][jtype].lj1:params(itype,jtype).lj1);
151 const F_FLOAT lj_2 = (STACKPARAMS?m_params[itype][jtype].lj2:params(itype,jtype).lj2);
152
153 /*if (ljt == LJ12_4) {
154
155 const F_FLOAT r4inv=r2inv*r2inv;
156 return r4inv*(lj_1*r4inv*r4inv - lj_2) * r2inv;
157
158 } else if (ljt == LJ9_6) {
159
160 const F_FLOAT r3inv = r2inv*sqrt(r2inv);
161 const F_FLOAT r6inv = r3inv*r3inv;
162 return r6inv*(lj_1*r3inv - lj_2) * r2inv;
163
164 } else if (ljt == LJ12_6) {
165
166 const double r6inv = r2inv*r2inv*r2inv;
167 return r6inv*(lj_1*r6inv - lj_2) * r2inv;
168
169 }
170 if (ljt!=LJ12_4 && ljt!=LJ9_6 && ljt!=LJ12_6) return 0.0;*/
171 const F_FLOAT r4inv=r2inv*r2inv;
172 const F_FLOAT r6inv=r2inv*r4inv;
173 const F_FLOAT a = ljt==LJ12_4?r4inv:r6inv;
174 const F_FLOAT b = ljt==LJ12_4?r4inv:(ljt==LJ9_6?1.0/sqrt(r2inv):r2inv);
175 return a* ( lj_1*r6inv*b - lj_2 * r2inv);
176 }
177
178 template<class DeviceType>
179 template<bool STACKPARAMS, class Specialisation>
180 KOKKOS_INLINE_FUNCTION
181 F_FLOAT PairLJSDKKokkos<DeviceType>::
compute_evdwl(const F_FLOAT & rsq,const int & i,const int & j,const int & itype,const int & jtype) const182 compute_evdwl(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
183 (void) i;
184 (void) j;
185 const F_FLOAT r2inv = 1.0/rsq;
186 const int ljt = (STACKPARAMS?m_params[itype][jtype].lj_type:params(itype,jtype).lj_type);
187
188 const F_FLOAT lj_3 = (STACKPARAMS?m_params[itype][jtype].lj3:params(itype,jtype).lj3);
189 const F_FLOAT lj_4 = (STACKPARAMS?m_params[itype][jtype].lj4:params(itype,jtype).lj4);
190 const F_FLOAT offset = (STACKPARAMS?m_params[itype][jtype].offset:params(itype,jtype).offset);
191
192 if (ljt == LJ12_4) {
193 const F_FLOAT r4inv=r2inv*r2inv;
194
195 return r4inv*(lj_3*r4inv*r4inv - lj_4) - offset;
196
197 } else if (ljt == LJ9_6) {
198 const F_FLOAT r3inv = r2inv*sqrt(r2inv);
199 const F_FLOAT r6inv = r3inv*r3inv;
200 return r6inv*(lj_3*r3inv - lj_4) - offset;
201
202 } else if (ljt == LJ12_6) {
203 const double r6inv = r2inv*r2inv*r2inv;
204 return r6inv*(lj_3*r6inv - lj_4) - offset;
205 } else
206 return 0.0;
207 }
208
209 /* ----------------------------------------------------------------------
210 allocate all arrays
211 ------------------------------------------------------------------------- */
212
213 template<class DeviceType>
allocate()214 void PairLJSDKKokkos<DeviceType>::allocate()
215 {
216 PairLJSDK::allocate();
217
218 int n = atom->ntypes;
219 memory->destroy(cutsq);
220 memoryKK->create_kokkos(k_cutsq,cutsq,n+1,n+1,"pair:cutsq");
221 d_cutsq = k_cutsq.template view<DeviceType>();
222 k_params = Kokkos::DualView<params_lj**,Kokkos::LayoutRight,DeviceType>("PairLJSDK::params",n+1,n+1);
223 params = k_params.template view<DeviceType>();
224 }
225
226 /* ----------------------------------------------------------------------
227 global settings
228 ------------------------------------------------------------------------- */
229
230 template<class DeviceType>
settings(int narg,char ** arg)231 void PairLJSDKKokkos<DeviceType>::settings(int narg, char **arg)
232 {
233 if (narg > 2) error->all(FLERR,"Illegal pair_style command");
234
235 PairLJSDK::settings(1,arg);
236 }
237
238 /* ----------------------------------------------------------------------
239 init specific to this pair style
240 ------------------------------------------------------------------------- */
241
242 template<class DeviceType>
init_style()243 void PairLJSDKKokkos<DeviceType>::init_style()
244 {
245 PairLJSDK::init_style();
246
247 // error if rRESPA with inner levels
248
249 if (update->whichflag == 1 && utils::strmatch(update->integrate_style,"^respa")) {
250 int respa = 0;
251 if (((Respa *) update->integrate)->level_inner >= 0) respa = 1;
252 if (((Respa *) update->integrate)->level_middle >= 0) respa = 2;
253 if (respa)
254 error->all(FLERR,"Cannot use Kokkos pair style with rRESPA inner/middle");
255 }
256
257 // irequest = neigh request made by parent class
258
259 neighflag = lmp->kokkos->neighflag;
260 int irequest = neighbor->nrequest - 1;
261
262 neighbor->requests[irequest]->
263 kokkos_host = std::is_same<DeviceType,LMPHostType>::value &&
264 !std::is_same<DeviceType,LMPDeviceType>::value;
265 neighbor->requests[irequest]->
266 kokkos_device = std::is_same<DeviceType,LMPDeviceType>::value;
267
268 if (neighflag == FULL) {
269 neighbor->requests[irequest]->full = 1;
270 neighbor->requests[irequest]->half = 0;
271 } else if (neighflag == HALF || neighflag == HALFTHREAD) {
272 neighbor->requests[irequest]->full = 0;
273 neighbor->requests[irequest]->half = 1;
274 } else {
275 error->all(FLERR,"Cannot use chosen neighbor list style with lj/sdk/kk");
276 }
277 }
278
279 /* ----------------------------------------------------------------------
280 init for one type pair i,j and corresponding j,i
281 ------------------------------------------------------------------------- */
282
283 template<class DeviceType>
init_one(int i,int j)284 double PairLJSDKKokkos<DeviceType>::init_one(int i, int j)
285 {
286 double cutone = PairLJSDK::init_one(i,j);
287
288 k_params.h_view(i,j).lj1 = lj1[i][j];
289 k_params.h_view(i,j).lj2 = lj2[i][j];
290 k_params.h_view(i,j).lj3 = lj3[i][j];
291 k_params.h_view(i,j).lj4 = lj4[i][j];
292 k_params.h_view(i,j).offset = offset[i][j];
293 k_params.h_view(i,j).cutsq = cutone*cutone;
294 k_params.h_view(i,j).lj_type = lj_type[i][j];
295 k_params.h_view(j,i) = k_params.h_view(i,j);
296 if (i<MAX_TYPES_STACKPARAMS+1 && j<MAX_TYPES_STACKPARAMS+1) {
297 m_params[i][j] = m_params[j][i] = k_params.h_view(i,j);
298 m_cutsq[j][i] = m_cutsq[i][j] = cutone*cutone;
299 }
300
301 k_cutsq.h_view(i,j) = k_cutsq.h_view(j,i) = cutone*cutone;
302 k_cutsq.template modify<LMPHostType>();
303 k_params.template modify<LMPHostType>();
304
305 return cutone;
306 }
307
308
309
310 namespace LAMMPS_NS {
311 template class PairLJSDKKokkos<LMPDeviceType>;
312 #ifdef LMP_KOKKOS_GPU
313 template class PairLJSDKKokkos<LMPHostType>;
314 #endif
315 }
316
317