~ubuntu-branches/debian/sid/lammps/sid

« back to all changes in this revision

Viewing changes to src/KOKKOS/pair_lj_cut_kokkos.cpp

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2015-04-29 23:44:49 UTC
  • mfrom: (5.1.3 experimental)
  • Revision ID: package-import@ubuntu.com-20150429234449-mbhy9utku6hp6oq8
Tags: 0~20150313.gitfa668e1-1
Upload into unstable.

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/* ----------------------------------------------------------------------
 
2
   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
 
3
   http://lammps.sandia.gov, Sandia National Laboratories
 
4
   Steve Plimpton, sjplimp@sandia.gov
 
5
 
 
6
   Copyright (2003) Sandia Corporation.  Under the terms of Contract
 
7
   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
 
8
   certain rights in this software.  This software is distributed under
 
9
   the GNU General Public License.
 
10
 
 
11
   See the README file in the top-level LAMMPS directory.
 
12
------------------------------------------------------------------------- */
 
13
 
 
14
#include "math.h"
 
15
#include "stdio.h"
 
16
#include "stdlib.h"
 
17
#include "string.h"
 
18
#include "pair_lj_cut_kokkos.h"
 
19
#include "kokkos.h"
 
20
#include "atom_kokkos.h"
 
21
#include "comm.h"
 
22
#include "force.h"
 
23
#include "neighbor.h"
 
24
#include "neigh_list.h"
 
25
#include "neigh_request.h"
 
26
#include "update.h"
 
27
#include "integrate.h"
 
28
#include "respa.h"
 
29
#include "math_const.h"
 
30
#include "memory.h"
 
31
#include "error.h"
 
32
#include "atom_masks.h"
 
33
 
 
34
using namespace LAMMPS_NS;
 
35
using namespace MathConst;
 
36
 
 
37
#define KOKKOS_CUDA_MAX_THREADS 256
 
38
#define KOKKOS_CUDA_MIN_BLOCKS 8
 
39
 
 
40
/* ---------------------------------------------------------------------- */
 
41
 
 
42
template<class DeviceType>
 
43
PairLJCutKokkos<DeviceType>::PairLJCutKokkos(LAMMPS *lmp) : PairLJCut(lmp)
 
44
{
 
45
  respa_enable = 0;
 
46
 
 
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
  cutsq = NULL;
 
52
}
 
53
 
 
54
/* ---------------------------------------------------------------------- */
 
55
 
 
56
template<class DeviceType>
 
57
PairLJCutKokkos<DeviceType>::~PairLJCutKokkos()
 
58
{
 
59
  if (allocated) {
 
60
    k_cutsq = DAT::tdual_ffloat_2d();
 
61
    memory->sfree(cutsq);
 
62
    cutsq = NULL;
 
63
  }
 
64
}
 
65
 
 
66
/* ---------------------------------------------------------------------- */
 
67
 
 
68
template<class DeviceType>
 
69
void PairLJCutKokkos<DeviceType>::cleanup_copy() {
 
70
  // WHY needed: this prevents parent copy from deallocating any arrays
 
71
  allocated = 0;
 
72
  cutsq = NULL;
 
73
  eatom = NULL;
 
74
  vatom = NULL;
 
75
}
 
76
 
 
77
/* ---------------------------------------------------------------------- */
 
78
 
 
79
template<class DeviceType>
 
80
void PairLJCutKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
 
81
{
 
82
  eflag = eflag_in;
 
83
  vflag = vflag_in;
 
84
 
 
85
 
 
86
  if (neighflag == FULL || neighflag == FULLCLUSTER) no_virial_fdotr_compute = 1;
 
87
 
 
88
  if (eflag || vflag) ev_setup(eflag,vflag);
 
89
  else evflag = vflag_fdotr = 0;
 
90
 
 
91
  atomKK->sync(execution_space,datamask_read);
 
92
  k_cutsq.template sync<DeviceType>();
 
93
  k_params.template sync<DeviceType>();
 
94
  if (eflag || vflag) atomKK->modified(execution_space,datamask_modify);
 
95
  else atomKK->modified(execution_space,F_MASK);
 
96
 
 
97
  x = atomKK->k_x.view<DeviceType>();
 
98
  c_x = atomKK->k_x.view<DeviceType>();
 
99
  f = atomKK->k_f.view<DeviceType>();
 
100
  type = atomKK->k_type.view<DeviceType>();
 
101
  tag = atomKK->k_tag.view<DeviceType>();
 
102
  nlocal = atom->nlocal;
 
103
  nall = atom->nlocal + atom->nghost;
 
104
  newton_pair = force->newton_pair;
 
105
  special_lj[0] = force->special_lj[0];
 
106
  special_lj[1] = force->special_lj[1];
 
107
  special_lj[2] = force->special_lj[2];
 
108
  special_lj[3] = force->special_lj[3];
 
109
 
 
110
  // loop over neighbors of my atoms
 
111
 
 
112
  EV_FLOAT ev = pair_compute<PairLJCutKokkos<DeviceType>,void >(this,(NeighListKokkos<DeviceType>*)list);
 
113
  DeviceType::fence();
 
114
 
 
115
  if (eflag) eng_vdwl += ev.evdwl;
 
116
  if (vflag_global) {
 
117
    virial[0] += ev.v[0];
 
118
    virial[1] += ev.v[1];
 
119
    virial[2] += ev.v[2];
 
120
    virial[3] += ev.v[3];
 
121
    virial[4] += ev.v[4];
 
122
    virial[5] += ev.v[5];
 
123
  }
 
124
 
 
125
  if (vflag_fdotr) pair_virial_fdotr_compute(this);
 
126
}
 
127
 
 
128
template<class DeviceType>
 
129
template<bool STACKPARAMS, class Specialisation>
 
130
KOKKOS_INLINE_FUNCTION
 
131
F_FLOAT PairLJCutKokkos<DeviceType>::
 
132
compute_fpair(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
 
133
  (void) i;
 
134
  (void) j;
 
135
  const F_FLOAT r2inv = 1.0/rsq;
 
136
  const F_FLOAT r6inv = r2inv*r2inv*r2inv;
 
137
 
 
138
  const F_FLOAT forcelj = r6inv *
 
139
    ((STACKPARAMS?m_params[itype][jtype].lj1:params(itype,jtype).lj1)*r6inv -
 
140
     (STACKPARAMS?m_params[itype][jtype].lj2:params(itype,jtype).lj2));
 
141
 
 
142
  return forcelj*r2inv;
 
143
}
 
144
 
 
145
template<class DeviceType>
 
146
template<bool STACKPARAMS, class Specialisation>
 
147
KOKKOS_INLINE_FUNCTION
 
148
F_FLOAT PairLJCutKokkos<DeviceType>::
 
149
compute_evdwl(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
 
150
  (void) i;
 
151
  (void) j;
 
152
  const F_FLOAT r2inv = 1.0/rsq;
 
153
  const F_FLOAT r6inv = r2inv*r2inv*r2inv;
 
154
 
 
155
  return r6inv*((STACKPARAMS?m_params[itype][jtype].lj3:params(itype,jtype).lj3)*r6inv -
 
156
                (STACKPARAMS?m_params[itype][jtype].lj4:params(itype,jtype).lj4)) -
 
157
                (STACKPARAMS?m_params[itype][jtype].offset:params(itype,jtype).offset);
 
158
}
 
159
 
 
160
/* ----------------------------------------------------------------------
 
161
   allocate all arrays
 
162
------------------------------------------------------------------------- */
 
163
 
 
164
template<class DeviceType>
 
165
void PairLJCutKokkos<DeviceType>::allocate()
 
166
{
 
167
  PairLJCut::allocate();
 
168
 
 
169
  int n = atom->ntypes;
 
170
  memory->destroy(cutsq);
 
171
  memory->create_kokkos(k_cutsq,cutsq,n+1,n+1,"pair:cutsq");
 
172
  d_cutsq = k_cutsq.template view<DeviceType>();
 
173
  k_params = Kokkos::DualView<params_lj**,Kokkos::LayoutRight,DeviceType>("PairLJCut::params",n+1,n+1);
 
174
  params = k_params.d_view;
 
175
}
 
176
 
 
177
/* ----------------------------------------------------------------------
 
178
   global settings
 
179
------------------------------------------------------------------------- */
 
180
 
 
181
template<class DeviceType>
 
182
void PairLJCutKokkos<DeviceType>::settings(int narg, char **arg)
 
183
{
 
184
  if (narg > 2) error->all(FLERR,"Illegal pair_style command");
 
185
 
 
186
  PairLJCut::settings(1,arg);
 
187
}
 
188
 
 
189
/* ----------------------------------------------------------------------
 
190
   init specific to this pair style
 
191
------------------------------------------------------------------------- */
 
192
 
 
193
template<class DeviceType>
 
194
void PairLJCutKokkos<DeviceType>::init_style()
 
195
{
 
196
  PairLJCut::init_style();
 
197
 
 
198
  // error if rRESPA with inner levels
 
199
 
 
200
  if (update->whichflag == 1 && strstr(update->integrate_style,"respa")) {
 
201
    int respa = 0;
 
202
    if (((Respa *) update->integrate)->level_inner >= 0) respa = 1;
 
203
    if (((Respa *) update->integrate)->level_middle >= 0) respa = 2;
 
204
    if (respa) 
 
205
      error->all(FLERR,"Cannot use Kokkos pair style with rRESPA inner/middle");
 
206
  }
 
207
 
 
208
  // irequest = neigh request made by parent class
 
209
 
 
210
  neighflag = lmp->kokkos->neighflag;
 
211
  int irequest = neighbor->nrequest - 1;
 
212
 
 
213
  neighbor->requests[irequest]->
 
214
    kokkos_host = Kokkos::Impl::is_same<DeviceType,LMPHostType>::value && 
 
215
    !Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
 
216
  neighbor->requests[irequest]->
 
217
    kokkos_device = Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
 
218
 
 
219
  if (neighflag == FULL) {
 
220
    neighbor->requests[irequest]->full = 1;
 
221
    neighbor->requests[irequest]->half = 0;
 
222
    neighbor->requests[irequest]->full_cluster = 0;
 
223
  } else if (neighflag == HALF || neighflag == HALFTHREAD) {
 
224
    neighbor->requests[irequest]->full = 0;
 
225
    neighbor->requests[irequest]->half = 1;
 
226
    neighbor->requests[irequest]->full_cluster = 0;
 
227
  } else if (neighflag == N2) {
 
228
    neighbor->requests[irequest]->full = 0;
 
229
    neighbor->requests[irequest]->half = 0;
 
230
    neighbor->requests[irequest]->full_cluster = 0;
 
231
  } else if (neighflag == FULLCLUSTER) {
 
232
    neighbor->requests[irequest]->full_cluster = 1;
 
233
    neighbor->requests[irequest]->full = 1;
 
234
    neighbor->requests[irequest]->half = 0;
 
235
  } else {
 
236
    error->all(FLERR,"Cannot use chosen neighbor list style with lj/cut/kk");
 
237
  }
 
238
}
 
239
 
 
240
/* ----------------------------------------------------------------------
 
241
   init for one type pair i,j and corresponding j,i
 
242
------------------------------------------------------------------------- */
 
243
 
 
244
template<class DeviceType>
 
245
double PairLJCutKokkos<DeviceType>::init_one(int i, int j)
 
246
{
 
247
  double cutone = PairLJCut::init_one(i,j);
 
248
 
 
249
  k_params.h_view(i,j).lj1 = lj1[i][j];
 
250
  k_params.h_view(i,j).lj2 = lj2[i][j];
 
251
  k_params.h_view(i,j).lj3 = lj3[i][j];
 
252
  k_params.h_view(i,j).lj4 = lj4[i][j];
 
253
  k_params.h_view(i,j).offset = offset[i][j];
 
254
  k_params.h_view(i,j).cutsq = cutone*cutone;
 
255
  k_params.h_view(j,i) = k_params.h_view(i,j);
 
256
  if(i<MAX_TYPES_STACKPARAMS+1 && j<MAX_TYPES_STACKPARAMS+1) {
 
257
    m_params[i][j] = m_params[j][i] = k_params.h_view(i,j);
 
258
    m_cutsq[j][i] = m_cutsq[i][j] = cutone*cutone;
 
259
  }
 
260
  k_cutsq.h_view(i,j) = cutone*cutone;
 
261
  k_cutsq.template modify<LMPHostType>();
 
262
  k_params.template modify<LMPHostType>();
 
263
 
 
264
  return cutone;
 
265
}
 
266
 
 
267
 
 
268
 
 
269
template class PairLJCutKokkos<LMPDeviceType>;
 
270
#ifdef KOKKOS_HAVE_CUDA
 
271
template class PairLJCutKokkos<LMPHostType>;
 
272
#endif