← Back to branch summary

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

« back to all changes in this revision

Viewing changes to src/KOKKOS/fix_langevin_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:
src/KOKKOS/fix_langevin_kokkos.cpp
 
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 "string.h"
 
17
#include "fix_langevin_kokkos.h"
 
18
#include "atom_masks.h"
 
19
#include "atom_kokkos.h"
 
20
#include "force.h"
 
21
#include "update.h"
 
22
#include "respa.h"
 
23
#include "error.h"
 
24
#include "memory.h"
 
25
#include "group.h"
 
26
#include "random_mars.h"
 
27
#include "compute.h"
 
28
#include "comm.h"
 
29
#include "modify.h"
 
30
#include "input.h"
 
31
#include "variable.h"
 
32
 
 
33
using namespace LAMMPS_NS;
 
34
using namespace FixConst;
 
35
 
 
36
enum{NOBIAS,BIAS};
 
37
enum{CONSTANT,EQUAL,ATOM};
 
38
#define SINERTIA 0.4          // moment of inertia prefactor for sphere
 
39
#define EINERTIA 0.2          // moment of inertia prefactor for ellipsoid
 
40
 
 
41
/* ---------------------------------------------------------------------- */
 
42
 
 
43
template<class DeviceType>
 
44
FixLangevinKokkos<DeviceType>::FixLangevinKokkos(LAMMPS *lmp, int narg, char **arg) :
 
45
  FixLangevin(lmp, narg, arg),rand_pool(seed + comm->me)
 
46
{
 
47
  atomKK = (AtomKokkos *) atom;
 
48
  int ntypes = atomKK->ntypes;
 
49
 
 
50
  // allocate per-type arrays for force prefactors
 
51
  memory->create_kokkos(k_gfactor1,gfactor1,ntypes+1,"langevin:gfactor1");
 
52
  memory->create_kokkos(k_gfactor2,gfactor2,ntypes+1,"langevin:gfactor2");
 
53
  memory->create_kokkos(k_ratio,ratio,ntypes+1,"langevin:ratio");
 
54
  d_gfactor1 = k_gfactor1.template view<DeviceType>();
 
55
  h_gfactor1 = k_gfactor1.template view<LMPHostType>();
 
56
  d_gfactor2 = k_gfactor2.template view<DeviceType>();
 
57
  h_gfactor2 = k_gfactor2.template view<LMPHostType>();
 
58
  d_ratio = k_ratio.template view<DeviceType>();
 
59
  h_ratio = k_ratio.template view<LMPHostType>();
 
60
 
 
61
  // optional args
 
62
  for (int i = 1; i <= ntypes; i++) ratio[i] = 1.0;
 
63
  k_ratio.template modify<LMPHostType>();
 
64
 
 
65
  if(gjfflag){
 
66
    nvalues = 3;
 
67
    grow_arrays(atomKK->nmax);
 
68
    atom->add_callback(0);
 
69
    // initialize franprev to zero
 
70
    for (int i = 0; i < atomKK->nlocal; i++) {
 
71
      franprev[i][0] = 0.0;
 
72
      franprev[i][1] = 0.0;
 
73
      franprev[i][2] = 0.0;
 
74
    }
 
75
    k_franprev.template modify<LMPHostType>();
 
76
  }
 
77
  if(zeroflag){
 
78
    k_fsumall = tdual_double_1d_3n("langevin:fsumall");
 
79
    h_fsumall = k_fsumall.template view<LMPHostType>();
 
80
    d_fsumall = k_fsumall.template view<DeviceType>();
 
81
  }
 
82
 
 
83
  execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
 
84
  datamask_read =  V_MASK | F_MASK | MASK_MASK | RMASS_MASK | TYPE_MASK;
 
85
  datamask_modify = F_MASK;
 
86
 
 
87
}
 
88
 
 
89
/* ---------------------------------------------------------------------- */
 
90
 
 
91
template<class DeviceType>
 
92
FixLangevinKokkos<DeviceType>::~FixLangevinKokkos()
 
93
{
 
94
  memory->destroy_kokkos(k_gfactor1,gfactor1);
 
95
  memory->destroy_kokkos(k_gfactor2,gfactor2);
 
96
  memory->destroy_kokkos(k_ratio,ratio);
 
97
  memory->destroy_kokkos(k_flangevin,flangevin);
 
98
  if(gjfflag) memory->destroy_kokkos(k_franprev,franprev);
 
99
  memory->destroy_kokkos(k_tforce,tforce);
 
100
}
 
101
 
 
102
/* ---------------------------------------------------------------------- */
 
103
 
 
104
template<class DeviceType>
 
105
void FixLangevinKokkos<DeviceType>::init()
 
106
{
 
107
  FixLangevin::init();
 
108
  if(oflag)
 
109
    error->all(FLERR,"Fix langevin omega is not yet implemented with kokkos");
 
110
  if(ascale)
 
111
    error->all(FLERR,"Fix langevin angmom is not yet implemented with kokkos");
 
112
 
 
113
  // prefactors are modified in the init
 
114
  k_gfactor1.template modify<LMPHostType>();
 
115
  k_gfactor2.template modify<LMPHostType>();
 
116
}
 
117
 
 
118
/* ---------------------------------------------------------------------- */
 
119
 
 
120
template<class DeviceType>
 
121
void FixLangevinKokkos<DeviceType>::grow_arrays(int nmax)
 
122
{
 
123
  memory->grow_kokkos(k_franprev,franprev,nmax,3,"langevin:franprev");
 
124
  d_franprev = k_franprev.template view<DeviceType>();
 
125
  h_franprev = k_franprev.template view<LMPHostType>();
 
126
}
 
127
 
 
128
/* ---------------------------------------------------------------------- */
 
129
 
 
130
template<class DeviceType>
 
131
void FixLangevinKokkos<DeviceType>::post_force(int vflag)
 
132
{
 
133
  // sync the device views which might have been modified on host
 
134
  atomKK->sync(execution_space,datamask_read);
 
135
  rmass = atomKK->rmass;
 
136
  f = atomKK->k_f.template view<DeviceType>();
 
137
  v = atomKK->k_v.template view<DeviceType>();
 
138
  type = atomKK->k_type.template view<DeviceType>();
 
139
  mask = atomKK->k_mask.template view<DeviceType>();
 
140
 
 
141
  k_gfactor1.template sync<DeviceType>();
 
142
  k_gfactor2.template sync<DeviceType>();
 
143
  k_ratio.template sync<DeviceType>();
 
144
  if(gjfflag) k_franprev.template sync<DeviceType>();
 
145
 
 
146
  boltz = force->boltz;
 
147
  dt = update->dt;
 
148
  mvv2e = force->mvv2e;
 
149
  ftm2v = force->ftm2v;
 
150
  fran_prop_const = sqrt(24.0*boltz/t_period/dt/mvv2e);
 
151
 
 
152
  compute_target(); // modifies tforce vector, hence sync here
 
153
  k_tforce.template sync<DeviceType>();
 
154
 
 
155
  double fsum[3],fsumall[3];
 
156
  bigint count;
 
157
  int nlocal = atomKK->nlocal;
 
158
 
 
159
  if (zeroflag) {
 
160
    fsum[0] = fsum[1] = fsum[2] = 0.0;
 
161
    count = group->count(igroup);
 
162
    if (count == 0)
 
163
      error->all(FLERR,"Cannot zero Langevin force of 0 atoms");
 
164
  }
 
165
 
 
166
  // reallocate flangevin if necessary
 
167
  if (tallyflag) {
 
168
    if (nlocal > maxatom1) {
 
169
      memory->destroy_kokkos(k_flangevin,flangevin);
 
170
      maxatom1 = atomKK->nmax;
 
171
      memory->create_kokkos(k_flangevin,flangevin,maxatom1,3,"langevin:flangevin");
 
172
      d_flangevin = k_flangevin.template view<DeviceType>();
 
173
      h_flangevin = k_flangevin.template view<LMPHostType>();
 
174
    }
 
175
  }
 
176
 
 
177
  // account for bias velocity
 
178
  if(tbiasflag == BIAS){
 
179
    temperature->compute_scalar();
 
180
    temperature->remove_bias_all(); // modifies velocities
 
181
    // if temeprature compute is kokkosized host-devcie comm won't be needed
 
182
    atomKK->modified(Host,V_MASK);
 
183
    atomKK->sync(execution_space,V_MASK);
 
184
  }
 
185
 
 
186
  // compute langevin force in parallel on the device
 
187
  FSUM s_fsum;
 
188
  if (tstyle == ATOM)
 
189
    if (gjfflag)
 
190
      if (tallyflag)
 
191
        if (tbiasflag == BIAS)
 
192
          if (rmass)
 
193
            if (zeroflag) {
 
194
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,1,1,1> post_functor(this);
 
195
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
196
            }
 
197
            else{
 
198
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,1,1,0> post_functor(this);
 
199
              Kokkos::parallel_for(nlocal,post_functor);
 
200
            }
 
201
          else
 
202
            if (zeroflag) {
 
203
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,1,0,1> post_functor(this);
 
204
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
205
            }
 
206
            else {
 
207
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,1,0,0> post_functor(this);
 
208
              Kokkos::parallel_for(nlocal,post_functor);
 
209
            }
 
210
        else
 
211
          if (rmass)
 
212
            if (zeroflag) {
 
213
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,0,1,1> post_functor(this);
 
214
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
215
            }
 
216
            else {
 
217
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,0,1,0> post_functor(this);
 
218
              Kokkos::parallel_for(nlocal,post_functor);
 
219
            }
 
220
          else
 
221
            if (zeroflag) {
 
222
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,0,0,1> post_functor(this);
 
223
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
224
            }
 
225
            else{
 
226
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,1,0,0,0> post_functor(this);
 
227
              Kokkos::parallel_for(nlocal,post_functor);
 
228
            }
 
229
      else
 
230
        if (tbiasflag == BIAS)
 
231
          if (rmass)
 
232
            if (zeroflag) {
 
233
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,1,1,1> post_functor(this);
 
234
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
235
            }
 
236
            else {
 
237
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,1,1,0> post_functor(this);
 
238
              Kokkos::parallel_for(nlocal,post_functor);
 
239
            }
 
240
          else
 
241
            if (zeroflag) {
 
242
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,1,0,1> post_functor(this);
 
243
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
244
            }
 
245
            else {
 
246
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,1,0,0> post_functor(this);
 
247
              Kokkos::parallel_for(nlocal,post_functor);
 
248
            }
 
249
        else
 
250
          if (rmass)
 
251
            if (zeroflag) {
 
252
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,0,1,1> post_functor(this);
 
253
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
254
            }
 
255
            else {
 
256
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,0,1,0> post_functor(this);
 
257
              Kokkos::parallel_for(nlocal,post_functor);
 
258
            }
 
259
          else
 
260
            if (zeroflag) {
 
261
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,0,0,1> post_functor(this);
 
262
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
263
            }
 
264
            else {
 
265
              FixLangevinKokkosPostForceFunctor<DeviceType,1,1,0,0,0,0> post_functor(this);
 
266
              Kokkos::parallel_for(nlocal,post_functor);
 
267
            }
 
268
    else
 
269
      if (tallyflag)
 
270
        if (tbiasflag == BIAS)
 
271
          if (rmass)
 
272
            if (zeroflag) {
 
273
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,1,1,1> post_functor(this);
 
274
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
275
            }
 
276
            else {
 
277
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,1,1,0> post_functor(this);
 
278
              Kokkos::parallel_for(nlocal,post_functor);
 
279
            }
 
280
          else
 
281
            if (zeroflag) {
 
282
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,1,0,1> post_functor(this);
 
283
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
284
            }
 
285
            else {
 
286
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,1,0,0> post_functor(this);
 
287
              Kokkos::parallel_for(nlocal,post_functor);
 
288
            }
 
289
        else
 
290
          if (rmass)
 
291
            if (zeroflag) {
 
292
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,0,1,1> post_functor(this);
 
293
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
294
            }
 
295
            else {
 
296
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,0,1,0> post_functor(this);
 
297
              Kokkos::parallel_for(nlocal,post_functor);
 
298
            }
 
299
          else
 
300
            if (zeroflag) {
 
301
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,0,0,1> post_functor(this);
 
302
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
303
            }
 
304
            else {
 
305
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,1,0,0,0> post_functor(this);
 
306
              Kokkos::parallel_for(nlocal,post_functor);
 
307
            }
 
308
      else
 
309
        if (tbiasflag == BIAS)
 
310
          if (rmass)
 
311
            if (zeroflag) {
 
312
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,1,1,1> post_functor(this);
 
313
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
314
            }
 
315
            else {
 
316
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,1,1,0> post_functor(this);
 
317
              Kokkos::parallel_for(nlocal,post_functor);
 
318
            }
 
319
          else
 
320
            if (zeroflag) {
 
321
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,1,0,1> post_functor(this);
 
322
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
323
            }
 
324
            else {
 
325
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,1,0,0> post_functor(this);
 
326
              Kokkos::parallel_for(nlocal,post_functor);
 
327
            }
 
328
        else
 
329
          if (rmass)
 
330
            if (zeroflag) {
 
331
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,0,1,1> post_functor(this);
 
332
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
333
            }
 
334
            else {
 
335
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,0,1,0> post_functor(this);
 
336
              Kokkos::parallel_for(nlocal,post_functor);
 
337
            }
 
338
          else
 
339
            if (zeroflag) {
 
340
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,0,0,1> post_functor(this);
 
341
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
342
            }
 
343
            else {
 
344
              FixLangevinKokkosPostForceFunctor<DeviceType,1,0,0,0,0,0> post_functor(this);
 
345
              Kokkos::parallel_for(nlocal,post_functor);
 
346
            }
 
347
  else
 
348
    if (gjfflag)
 
349
      if (tallyflag)
 
350
        if (tbiasflag == BIAS)
 
351
          if (rmass)
 
352
            if (zeroflag) {
 
353
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,1,1,1> post_functor(this);
 
354
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
355
            }
 
356
            else {
 
357
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,1,1,0> post_functor(this);
 
358
              Kokkos::parallel_for(nlocal,post_functor);
 
359
            }
 
360
          else
 
361
            if (zeroflag) {
 
362
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,1,0,1> post_functor(this);
 
363
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
364
            }
 
365
            else {
 
366
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,1,0,0> post_functor(this);
 
367
              Kokkos::parallel_for(nlocal,post_functor);
 
368
            }
 
369
        else
 
370
          if (rmass)
 
371
            if (zeroflag) {
 
372
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,0,1,1> post_functor(this);
 
373
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
374
            }
 
375
            else {
 
376
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,0,1,0> post_functor(this);
 
377
              Kokkos::parallel_for(nlocal,post_functor);
 
378
            }
 
379
          else
 
380
            if (zeroflag) {
 
381
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,0,0,1> post_functor(this);
 
382
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
383
            }
 
384
            else {
 
385
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,1,0,0,0> post_functor(this);
 
386
              Kokkos::parallel_for(nlocal,post_functor);
 
387
            }
 
388
      else
 
389
        if (tbiasflag == BIAS)
 
390
          if (rmass)
 
391
            if (zeroflag) {
 
392
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,1,1,1> post_functor(this);
 
393
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
394
            }
 
395
            else {
 
396
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,1,1,0> post_functor(this);
 
397
              Kokkos::parallel_for(nlocal,post_functor);
 
398
            }
 
399
          else
 
400
            if (zeroflag) {
 
401
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,1,0,1> post_functor(this);
 
402
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
403
            }
 
404
            else {
 
405
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,1,0,0> post_functor(this);
 
406
              Kokkos::parallel_for(nlocal,post_functor);
 
407
            }
 
408
        else
 
409
          if (rmass)
 
410
            if (zeroflag) {
 
411
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,0,1,1> post_functor(this);
 
412
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
413
            }
 
414
            else {
 
415
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,0,1,0> post_functor(this);
 
416
              Kokkos::parallel_for(nlocal,post_functor);
 
417
            }
 
418
          else
 
419
            if (zeroflag) {
 
420
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,0,0,1> post_functor(this);
 
421
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
422
            }
 
423
            else {
 
424
              FixLangevinKokkosPostForceFunctor<DeviceType,0,1,0,0,0,0> post_functor(this);
 
425
              Kokkos::parallel_for(nlocal,post_functor);
 
426
            }
 
427
    else
 
428
      if (tallyflag)
 
429
        if (tbiasflag == BIAS)
 
430
          if (rmass)
 
431
            if (zeroflag) {
 
432
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,1,1,1> post_functor(this);
 
433
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
434
            }
 
435
            else {
 
436
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,1,1,0> post_functor(this);
 
437
              Kokkos::parallel_for(nlocal,post_functor);
 
438
            }
 
439
          else
 
440
            if (zeroflag) {
 
441
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,1,0,1> post_functor(this);
 
442
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
443
            }
 
444
            else {
 
445
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,1,0,0> post_functor(this);
 
446
              Kokkos::parallel_for(nlocal,post_functor);
 
447
            }
 
448
        else
 
449
          if (rmass)
 
450
            if (zeroflag) {
 
451
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,0,1,1> post_functor(this);
 
452
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
453
            }
 
454
            else {
 
455
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,0,1,0> post_functor(this);
 
456
              Kokkos::parallel_for(nlocal,post_functor);
 
457
            }
 
458
          else
 
459
            if (zeroflag) {
 
460
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,0,0,1> post_functor(this);
 
461
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
462
            }
 
463
            else {
 
464
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,1,0,0,0> post_functor(this);
 
465
              Kokkos::parallel_for(nlocal,post_functor);
 
466
            }
 
467
      else
 
468
        if (tbiasflag == BIAS)
 
469
          if (rmass)
 
470
            if (zeroflag) {
 
471
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,1,1,1> post_functor(this);
 
472
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
473
            }
 
474
            else {
 
475
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,1,1,0> post_functor(this);
 
476
              Kokkos::parallel_for(nlocal,post_functor);
 
477
            }
 
478
          else
 
479
            if (zeroflag) {
 
480
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,1,0,1> post_functor(this);
 
481
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
482
            }
 
483
            else {
 
484
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,1,0,0> post_functor(this);
 
485
              Kokkos::parallel_for(nlocal,post_functor);
 
486
            }
 
487
        else
 
488
          if (rmass)
 
489
            if (zeroflag) {
 
490
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,0,1,1> post_functor(this);
 
491
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
492
            }
 
493
            else {
 
494
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,0,1,0> post_functor(this);
 
495
              Kokkos::parallel_for(nlocal,post_functor);
 
496
            }
 
497
          else
 
498
            if (zeroflag) {
 
499
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,0,0,1> post_functor(this);
 
500
              Kokkos::parallel_reduce(nlocal,post_functor,s_fsum);
 
501
            }
 
502
            else {
 
503
              FixLangevinKokkosPostForceFunctor<DeviceType,0,0,0,0,0,0> post_functor(this);
 
504
              Kokkos::parallel_for(nlocal,post_functor);
 
505
            }
 
506
 
 
507
  DeviceType::fence();
 
508
 
 
509
  if(tbiasflag == BIAS){
 
510
    temperature->restore_bias_all(); // modifies velocities
 
511
    atomKK->modified(Host,V_MASK);
 
512
  }
 
513
 
 
514
  // set modify flags for the views modified in post_force functor
 
515
  if (gjfflag) k_franprev.template modify<DeviceType>();
 
516
  if (tallyflag) k_flangevin.template modify<DeviceType>();
 
517
 
 
518
  // set total force to zero
 
519
  if (zeroflag) {
 
520
    fsum[0] = s_fsum.fx; fsum[1] = s_fsum.fy; fsum[2] = s_fsum.fz;
 
521
    MPI_Allreduce(fsum,fsumall,3,MPI_DOUBLE,MPI_SUM,world);
 
522
    h_fsumall(0) = fsumall[0]/count;
 
523
    h_fsumall(1) = fsumall[1]/count;
 
524
    h_fsumall(2) = fsumall[2]/count;
 
525
    k_fsumall.template modify<LMPHostType>();
 
526
    k_fsumall.template sync<DeviceType>();
 
527
    // set total force zero in parallel on the device
 
528
    FixLangevinKokkosZeroForceFunctor<DeviceType> zero_functor(this);
 
529
    Kokkos::parallel_for(nlocal,zero_functor);
 
530
    DeviceType::fence();
 
531
  }
 
532
  // f is modified by both post_force and zero_force functors
 
533
  atomKK->modified(execution_space,datamask_modify);
 
534
 
 
535
  // thermostat omega and angmom
 
536
  //  if (oflag) omega_thermostat();
 
537
  //  if (ascale) angmom_thermostat();
 
538
 
 
539
}
 
540
 
 
541
/* ---------------------------------------------------------------------- */
 
542
 
 
543
template<class DeviceType>
 
544
template<int Tp_TSTYLEATOM, int Tp_GJF, int Tp_TALLY,
 
545
         int Tp_BIAS, int Tp_RMASS, int Tp_ZERO>
 
546
KOKKOS_INLINE_FUNCTION
 
547
FSUM FixLangevinKokkos<DeviceType>::post_force_item(int i) const
 
548
{
 
549
  FSUM fsum;
 
550
  double fdrag[3],fran[3];
 
551
  double gamma1,gamma2;
 
552
  double fswap;
 
553
  double tsqrt_t = tsqrt;
 
554
 
 
555
  if (mask[i] & groupbit) {
 
556
    rand_type rand_gen = rand_pool.get_state();
 
557
    if(Tp_TSTYLEATOM) tsqrt_t = sqrt(d_tforce[i]);
 
558
    if(Tp_RMASS){
 
559
      gamma1 = -rmass[i] / t_period / ftm2v;
 
560
      gamma2 = sqrt(rmass[i]) * fran_prop_const / ftm2v;
 
561
      gamma1 *= 1.0/d_ratio[type[i]];
 
562
      gamma2 *= 1.0/sqrt(d_ratio[type[i]]) * tsqrt_t;
 
563
    } else {
 
564
      gamma1 = d_gfactor1[type[i]];
 
565
      gamma2 = d_gfactor2[type[i]] * tsqrt_t;
 
566
    }
 
567
 
 
568
    fran[0] = gamma2 * (rand_gen.drand() - 0.5); //(random->uniform()-0.5);
 
569
    fran[1] = gamma2 * (rand_gen.drand() - 0.5); //(random->uniform()-0.5);
 
570
    fran[2] = gamma2 * (rand_gen.drand() - 0.5); //(random->uniform()-0.5);
 
571
 
 
572
    if(Tp_BIAS){
 
573
      fdrag[0] = gamma1*v(i,0);
 
574
      fdrag[1] = gamma1*v(i,1);
 
575
      fdrag[2] = gamma1*v(i,2);
 
576
      if (v(i,0) == 0.0) fran[0] = 0.0;
 
577
      if (v(i,1) == 0.0) fran[1] = 0.0;
 
578
      if (v(i,2) == 0.0) fran[2] = 0.0;
 
579
    }else{
 
580
      fdrag[0] = gamma1*v(i,0);
 
581
      fdrag[1] = gamma1*v(i,1);
 
582
      fdrag[2] = gamma1*v(i,2);
 
583
    }
 
584
 
 
585
    if (Tp_GJF) {
 
586
      fswap = 0.5*(fran[0]+d_franprev(i,0));
 
587
      d_franprev(i,0) = fran[0];
 
588
      fran[0] = fswap;
 
589
      fswap = 0.5*(fran[1]+d_franprev(i,1));
 
590
      d_franprev(i,1) = fran[1];
 
591
      fran[1] = fswap;
 
592
      fswap = 0.5*(fran[2]+d_franprev(i,2));
 
593
      d_franprev(i,2) = fran[2];
 
594
      fran[2] = fswap;
 
595
 
 
596
      fdrag[0] *= gjffac;
 
597
      fdrag[1] *= gjffac;
 
598
      fdrag[2] *= gjffac;
 
599
      fran[0] *= gjffac;
 
600
      fran[1] *= gjffac;
 
601
      fran[2] *= gjffac;
 
602
      f(i,0) *= gjffac;
 
603
      f(i,1) *= gjffac;
 
604
      f(i,2) *= gjffac;
 
605
    }
 
606
 
 
607
    f(i,0) += fdrag[0] + fran[0];
 
608
    f(i,1) += fdrag[1] + fran[1];
 
609
    f(i,2) += fdrag[2] + fran[2];
 
610
 
 
611
    if (Tp_TALLY) {
 
612
      d_flangevin(i,0) = fdrag[0] + fran[0];
 
613
      d_flangevin(i,1) = fdrag[1] + fran[1];
 
614
      d_flangevin(i,2) = fdrag[2] + fran[2];
 
615
    }
 
616
 
 
617
    if (Tp_ZERO) {
 
618
      fsum.fx = fran[0];
 
619
      fsum.fy = fran[1];
 
620
      fsum.fz = fran[2];
 
621
    }
 
622
    rand_pool.free_state(rand_gen);
 
623
  }
 
624
 
 
625
  return fsum;
 
626
}
 
627
 
 
628
/* ---------------------------------------------------------------------- */
 
629
 
 
630
template<class DeviceType>
 
631
KOKKOS_INLINE_FUNCTION
 
632
void FixLangevinKokkos<DeviceType>::zero_force_item(int i) const
 
633
{
 
634
  if (mask[i] & groupbit) {
 
635
    f(i,0) -= d_fsumall[0];
 
636
    f(i,1) -= d_fsumall[1];
 
637
    f(i,2) -= d_fsumall[2];
 
638
  }
 
639
 
 
640
}
 
641
 
 
642
/* ----------------------------------------------------------------------
 
643
   set current t_target and t_sqrt
 
644
   ------------------------------------------------------------------------- */
 
645
 
 
646
template<class DeviceType>
 
647
void FixLangevinKokkos<DeviceType>::compute_target()
 
648
{
 
649
  atomKK->sync(Host, MASK_MASK);
 
650
  mask = atomKK->k_mask.template view<DeviceType>();
 
651
  int nlocal = atomKK->nlocal;
 
652
 
 
653
  double delta = update->ntimestep - update->beginstep;
 
654
  if (delta != 0.0) delta /= update->endstep - update->beginstep;
 
655
 
 
656
  // if variable temp, evaluate variable, wrap with clear/add
 
657
  // reallocate tforce array if necessary
 
658
 
 
659
  if (tstyle == CONSTANT) {
 
660
    t_target = t_start + delta * (t_stop-t_start);
 
661
    tsqrt = sqrt(t_target);
 
662
  } else {
 
663
    modify->clearstep_compute();
 
664
    if (tstyle == EQUAL) {
 
665
      t_target = input->variable->compute_equal(tvar);
 
666
      if (t_target < 0.0)
 
667
        error->one(FLERR,"Fix langevin variable returned negative temperature");
 
668
      tsqrt = sqrt(t_target);
 
669
    } else {
 
670
      if (nlocal > maxatom2) {
 
671
        maxatom2 = atom->nmax;
 
672
        memory->destroy_kokkos(k_tforce,tforce);
 
673
        memory->create_kokkos(k_tforce,tforce,maxatom2,"langevin:tforce");
 
674
        d_tforce = k_tforce.template view<DeviceType>();
 
675
        h_tforce = k_tforce.template view<LMPHostType>();
 
676
      }
 
677
      input->variable->compute_atom(tvar,igroup,tforce,1,0); // tforce is modified on host
 
678
      k_tforce.template modify<LMPHostType>();
 
679
      for (int i = 0; i < nlocal; i++)
 
680
        if (mask[i] & groupbit)
 
681
          if (h_tforce[i] < 0.0)
 
682
            error->one(FLERR,
 
683
                       "Fix langevin variable returned negative temperature");
 
684
    }
 
685
    modify->addstep_compute(update->ntimestep + 1);
 
686
  }
 
687
}
 
688
 
 
689
/* ---------------------------------------------------------------------- */
 
690
 
 
691
template<class DeviceType>
 
692
void FixLangevinKokkos<DeviceType>::reset_dt()
 
693
{
 
694
  if (atomKK->mass) {
 
695
    for (int i = 1; i <= atomKK->ntypes; i++) {
 
696
      h_gfactor2[i] = sqrt(atomKK->mass[i]) *
 
697
        sqrt(24.0*force->boltz/t_period/update->dt/force->mvv2e) /
 
698
        force->ftm2v;
 
699
      h_gfactor2[i] *= 1.0/sqrt(h_ratio[i]);
 
700
    }
 
701
    k_gfactor2.template modify<LMPHostType>();
 
702
  }
 
703
 
 
704
}
 
705
 
 
706
/* ---------------------------------------------------------------------- */
 
707
 
 
708
template<class DeviceType>
 
709
double FixLangevinKokkos<DeviceType>::compute_scalar()
 
710
{
 
711
  if (!tallyflag || flangevin == NULL) return 0.0;
 
712
 
 
713
  v = atomKK->k_v.template view<DeviceType>();
 
714
  mask = atomKK->k_mask.template view<DeviceType>();
 
715
 
 
716
  // capture the very first energy transfer to thermal reservoir
 
717
 
 
718
  if (update->ntimestep == update->beginstep) {
 
719
    energy_onestep = 0.0;
 
720
    atomKK->sync(execution_space,V_MASK | MASK_MASK);
 
721
    int nlocal = atomKK->nlocal;
 
722
    k_flangevin.template sync<DeviceType>();
 
723
    FixLangevinKokkosTallyEnergyFunctor<DeviceType> scalar_functor(this);
 
724
    Kokkos::parallel_reduce(nlocal,scalar_functor,energy_onestep);
 
725
    DeviceType::fence();
 
726
    energy = 0.5*energy_onestep*update->dt;
 
727
  }
 
728
 
 
729
  // convert midstep energy back to previous fullstep energy
 
730
  double energy_me = energy - 0.5*energy_onestep*update->dt;
 
731
  double energy_all;
 
732
  MPI_Allreduce(&energy_me,&energy_all,1,MPI_DOUBLE,MPI_SUM,world);
 
733
  return -energy_all;
 
734
}
 
735
 
 
736
/* ---------------------------------------------------------------------- */
 
737
 
 
738
template<class DeviceType>
 
739
KOKKOS_INLINE_FUNCTION
 
740
double FixLangevinKokkos<DeviceType>::compute_energy_item(int i) const
 
741
{
 
742
  double energy;
 
743
  if (mask[i] & groupbit)
 
744
    energy = d_flangevin(i,0)*v(i,0) + d_flangevin(i,1)*v(i,1) +
 
745
      d_flangevin(i,2)*v(i,2);
 
746
  return energy;
 
747
}
 
748
 
 
749
/* ----------------------------------------------------------------------
 
750
   tally energy transfer to thermal reservoir
 
751
   ------------------------------------------------------------------------- */
 
752
 
 
753
template<class DeviceType>
 
754
void FixLangevinKokkos<DeviceType>::end_of_step()
 
755
{
 
756
  if (!tallyflag) return;
 
757
 
 
758
  v = atomKK->k_v.template view<DeviceType>();
 
759
  mask = atomKK->k_mask.template view<DeviceType>();
 
760
 
 
761
  atomKK->sync(execution_space,V_MASK | MASK_MASK);
 
762
  int nlocal = atomKK->nlocal;
 
763
 
 
764
  energy_onestep = 0.0;
 
765
 
 
766
  k_flangevin.template sync<DeviceType>();
 
767
  FixLangevinKokkosTallyEnergyFunctor<DeviceType> tally_functor(this);
 
768
  Kokkos::parallel_reduce(nlocal,tally_functor,energy_onestep);
 
769
  DeviceType::fence();
 
770
 
 
771
  energy += energy_onestep*update->dt;
 
772
}
 
773
 
 
774
/* ----------------------------------------------------------------------
 
775
   copy values within local atom-based array
 
776
   ------------------------------------------------------------------------- */
 
777
 
 
778
template<class DeviceType>
 
779
void FixLangevinKokkos<DeviceType>::copy_arrays(int i, int j, int delflag)
 
780
{
 
781
  for (int m = 0; m < nvalues; m++)
 
782
    h_franprev(j,m) = h_franprev(i,m);
 
783
 
 
784
  k_franprev.template modify<LMPHostType>();
 
785
 
 
786
}
 
787
 
 
788
/* ---------------------------------------------------------------------- */
 
789
 
 
790
template<class DeviceType>
 
791
void FixLangevinKokkos<DeviceType>::cleanup_copy()
 
792
{
 
793
  random = NULL;
 
794
  tstr = NULL;
 
795
  gfactor1 = NULL;
 
796
  gfactor2 = NULL;
 
797
  ratio = NULL;
 
798
  id_temp = NULL;
 
799
  flangevin = NULL;
 
800
  tforce = NULL;
 
801
  gjfflag = 0;
 
802
  franprev = NULL;
 
803
  id = style = NULL;
 
804
  vatom = NULL;
 
805
}
 
806
 
 
807
template class FixLangevinKokkos<LMPDeviceType>;
 
808
#ifdef KOKKOS_HAVE_CUDA
 
809
template class FixLangevinKokkos<LMPHostType>;
 
810
#endif