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 /* ----------------------------------------------------------------------
16 Contributing authors: Stefan Paquay (Eindhoven University of Technology)
17 ------------------------------------------------------------------------- */
18
19 #include "pair_morse_kokkos.h"
20 #include <cmath>
21 #include <cstring>
22 #include "kokkos.h"
23 #include "atom_kokkos.h"
24 #include "force.h"
25 #include "neighbor.h"
26 #include "neigh_list.h"
27 #include "neigh_request.h"
28 #include "update.h"
29 #include "respa.h"
30 #include "math_const.h"
31 #include "memory_kokkos.h"
32 #include "error.h"
33 #include "atom_masks.h"
34
35 using namespace LAMMPS_NS;
36 using namespace MathConst;
37
38 #define KOKKOS_CUDA_MAX_THREADS 256
39 #define KOKKOS_CUDA_MIN_BLOCKS 8
40
41 /* ---------------------------------------------------------------------- */
42
43 template<class DeviceType>
PairMorseKokkos(LAMMPS * lmp)44 PairMorseKokkos<DeviceType>::PairMorseKokkos(LAMMPS *lmp) : PairMorse(lmp)
45 {
46 respa_enable = 0;
47
48 kokkosable = 1;
49 atomKK = (AtomKokkos *) atom;
50 execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
51 datamask_read = X_MASK | F_MASK | TYPE_MASK | ENERGY_MASK | VIRIAL_MASK;
52 datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
53 }
54
55 /* ---------------------------------------------------------------------- */
56
57 template<class DeviceType>
~PairMorseKokkos()58 PairMorseKokkos<DeviceType>::~PairMorseKokkos()
59 {
60 if (copymode) return;
61
62 if (allocated) {
63 memoryKK->destroy_kokkos(k_eatom,eatom);
64 memoryKK->destroy_kokkos(k_vatom,vatom);
65 memoryKK->destroy_kokkos(k_cutsq,cutsq);
66 }
67 }
68
69 /* ---------------------------------------------------------------------- */
70
71 template<class DeviceType>
compute(int eflag_in,int vflag_in)72 void PairMorseKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
73 {
74 eflag = eflag_in;
75 vflag = vflag_in;
76
77
78 if (neighflag == FULL) no_virial_fdotr_compute = 1;
79
80 ev_init(eflag,vflag,0);
81
82 // reallocate per-atom arrays if necessary
83
84 if (eflag_atom) {
85 memoryKK->destroy_kokkos(k_eatom,eatom);
86 memoryKK->create_kokkos(k_eatom,eatom,maxeatom,"pair:eatom");
87 d_eatom = k_eatom.view<DeviceType>();
88 }
89 if (vflag_atom) {
90 memoryKK->destroy_kokkos(k_vatom,vatom);
91 memoryKK->create_kokkos(k_vatom,vatom,maxvatom,"pair:vatom");
92 d_vatom = k_vatom.view<DeviceType>();
93 }
94
95 atomKK->sync(execution_space,datamask_read);
96 k_cutsq.template sync<DeviceType>();
97 k_params.template sync<DeviceType>();
98 if (eflag || vflag) atomKK->modified(execution_space,datamask_modify);
99 else atomKK->modified(execution_space,F_MASK);
100
101 x = atomKK->k_x.view<DeviceType>();
102 c_x = atomKK->k_x.view<DeviceType>();
103 f = atomKK->k_f.view<DeviceType>();
104 type = atomKK->k_type.view<DeviceType>();
105 tag = atomKK->k_tag.view<DeviceType>();
106 nlocal = atom->nlocal;
107 nall = atom->nlocal + atom->nghost;
108 newton_pair = force->newton_pair;
109 special_lj[0] = force->special_lj[0];
110 special_lj[1] = force->special_lj[1];
111 special_lj[2] = force->special_lj[2];
112 special_lj[3] = force->special_lj[3];
113
114 // loop over neighbors of my atoms
115
116 EV_FLOAT ev = pair_compute<PairMorseKokkos<DeviceType>,void >(this,(NeighListKokkos<DeviceType>*)list);
117
118 if (eflag_global) eng_vdwl += ev.evdwl;
119 if (vflag_global) {
120 virial[0] += ev.v[0];
121 virial[1] += ev.v[1];
122 virial[2] += ev.v[2];
123 virial[3] += ev.v[3];
124 virial[4] += ev.v[4];
125 virial[5] += ev.v[5];
126 }
127
128 if (vflag_fdotr) pair_virial_fdotr_compute(this);
129
130 if (eflag_atom) {
131 k_eatom.template modify<DeviceType>();
132 k_eatom.template sync<LMPHostType>();
133 }
134
135 if (vflag_atom) {
136 k_vatom.template modify<DeviceType>();
137 k_vatom.template sync<LMPHostType>();
138 }
139 }
140
141 template<class DeviceType>
142 template<bool STACKPARAMS, class Specialisation>
143 KOKKOS_INLINE_FUNCTION
144 F_FLOAT PairMorseKokkos<DeviceType>::
compute_fpair(const F_FLOAT & rsq,const int & i,const int & j,const int & itype,const int & jtype) const145 compute_fpair(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
146 (void) i;
147 (void) j;
148 const F_FLOAT rr = sqrt(rsq);
149 const F_FLOAT r0 = STACKPARAMS ? m_params[itype][jtype].r0 : params(itype,jtype).r0;
150 const F_FLOAT d0 = STACKPARAMS ? m_params[itype][jtype].d0 : params(itype,jtype).d0;
151 const F_FLOAT aa = STACKPARAMS ? m_params[itype][jtype].alpha : params(itype,jtype).alpha;
152 const F_FLOAT dr = rr - r0;
153
154 // U = d0 * [ exp( -2*a*(x-r0)) - 2*exp(-a*(x-r0)) ]
155 // f = -2*a*d0*[ -exp( -2*a*(x-r0) ) + exp( -a*(x-r0) ) ] * grad(r)
156 // = +2*a*d0*[ exp( -2*a*(x-r0) ) - exp( -a*(x-r0) ) ] * grad(r)
157 const F_FLOAT dexp = exp( -aa*dr );
158 const F_FLOAT forcelj = 2*aa*d0*dexp*(dexp-1.0);
159
160 return forcelj / rr;
161 }
162
163 template<class DeviceType>
164 template<bool STACKPARAMS, class Specialisation>
165 KOKKOS_INLINE_FUNCTION
166 F_FLOAT PairMorseKokkos<DeviceType>::
compute_evdwl(const F_FLOAT & rsq,const int & i,const int & j,const int & itype,const int & jtype) const167 compute_evdwl(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
168 (void) i;
169 (void) j;
170 const F_FLOAT rr = sqrt(rsq);
171 const F_FLOAT r0 = STACKPARAMS ? m_params[itype][jtype].r0 : params(itype,jtype).r0;
172 const F_FLOAT d0 = STACKPARAMS ? m_params[itype][jtype].d0 : params(itype,jtype).d0;
173 const F_FLOAT aa = STACKPARAMS ? m_params[itype][jtype].alpha : params(itype,jtype).alpha;
174 const F_FLOAT dr = rr - r0;
175
176 // U = d0 * [ exp( -2*a*(x-r0)) - 2*exp(-a*(x-r0)) ]
177 // f = -2*a*d0*[ -exp( -2*a*(x-r0) ) + exp( -a*(x-r0) ) ] * grad(r)
178 // = +2*a*d0*[ exp( -2*a*(x-r0) ) - exp( -a*(x-r0) ) ] * grad(r)
179 const F_FLOAT dexp = exp( -aa*dr );
180
181 return d0 * dexp * ( dexp - 2.0 );
182 }
183
184 /* ----------------------------------------------------------------------
185 allocate all arrays
186 ------------------------------------------------------------------------- */
187
188 template<class DeviceType>
allocate()189 void PairMorseKokkos<DeviceType>::allocate()
190 {
191 PairMorse::allocate();
192
193 int n = atom->ntypes;
194 memory->destroy(cutsq);
195 memoryKK->create_kokkos(k_cutsq,cutsq,n+1,n+1,"pair:cutsq");
196 d_cutsq = k_cutsq.template view<DeviceType>();
197 k_params = Kokkos::DualView<params_morse**,Kokkos::LayoutRight,DeviceType>("PairMorse::params",n+1,n+1);
198 params = k_params.template view<DeviceType>();
199 }
200
201 /* ----------------------------------------------------------------------
202 global settings
203 ------------------------------------------------------------------------- */
204
205 template<class DeviceType>
settings(int narg,char ** arg)206 void PairMorseKokkos<DeviceType>::settings(int narg, char **arg)
207 {
208 if (narg > 2) error->all(FLERR,"Illegal pair_style command");
209
210 PairMorse::settings(1,arg);
211 }
212
213 /* ----------------------------------------------------------------------
214 init specific to this pair style
215 ------------------------------------------------------------------------- */
216
217 template<class DeviceType>
init_style()218 void PairMorseKokkos<DeviceType>::init_style()
219 {
220 PairMorse::init_style();
221
222 // error if rRESPA with inner levels
223
224 if (update->whichflag == 1 && utils::strmatch(update->integrate_style,"^respa")) {
225 int respa = 0;
226 if (((Respa *) update->integrate)->level_inner >= 0) respa = 1;
227 if (((Respa *) update->integrate)->level_middle >= 0) respa = 2;
228 if (respa)
229 error->all(FLERR,"Cannot use Kokkos pair style with rRESPA inner/middle");
230 }
231
232 // irequest = neigh request made by parent class
233
234 neighflag = lmp->kokkos->neighflag;
235 int irequest = neighbor->nrequest - 1;
236
237 neighbor->requests[irequest]->
238 kokkos_host = std::is_same<DeviceType,LMPHostType>::value &&
239 !std::is_same<DeviceType,LMPDeviceType>::value;
240 neighbor->requests[irequest]->
241 kokkos_device = std::is_same<DeviceType,LMPDeviceType>::value;
242
243 if (neighflag == FULL) {
244 neighbor->requests[irequest]->full = 1;
245 neighbor->requests[irequest]->half = 0;
246 } else if (neighflag == HALF || neighflag == HALFTHREAD) {
247 neighbor->requests[irequest]->full = 0;
248 neighbor->requests[irequest]->half = 1;
249 } else {
250 error->all(FLERR,"Cannot use chosen neighbor list style with morse/kk");
251 }
252 }
253
254 /* ----------------------------------------------------------------------
255 init for one type pair i,j and corresponding j,i
256 ------------------------------------------------------------------------- */
257 // Rewrite this.
258 template<class DeviceType>
init_one(int i,int j)259 double PairMorseKokkos<DeviceType>::init_one(int i, int j)
260 {
261 double cutone = PairMorse::init_one(i,j);
262
263 k_params.h_view(i,j).d0 = d0[i][j];
264 k_params.h_view(i,j).alpha = alpha[i][j];
265 k_params.h_view(i,j).r0 = r0[i][j];
266 k_params.h_view(i,j).offset = offset[i][j];
267 k_params.h_view(i,j).cutsq = cutone*cutone;
268 k_params.h_view(j,i) = k_params.h_view(i,j);
269
270 if (i<MAX_TYPES_STACKPARAMS+1 && j<MAX_TYPES_STACKPARAMS+1) {
271 m_params[i][j] = m_params[j][i] = k_params.h_view(i,j);
272 m_cutsq[j][i] = m_cutsq[i][j] = cutone*cutone;
273 }
274
275 k_cutsq.h_view(i,j) = k_cutsq.h_view(j,i) = cutone*cutone;
276 k_cutsq.template modify<LMPHostType>();
277 k_params.template modify<LMPHostType>();
278
279 return cutone;
280 }
281
282
283
284 namespace LAMMPS_NS {
285 template class PairMorseKokkos<LMPDeviceType>;
286 #ifdef LMP_KOKKOS_GPU
287 template class PairMorseKokkos<LMPHostType>;
288 #endif
289 }
290
291