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 "verlet_kokkos.h"
16 #include "neighbor.h"
17 #include "domain.h"
18 #include "comm.h"
19 #include "atom_kokkos.h"
20 #include "atom_masks.h"
21 #include "force.h"
22 #include "pair.h"
23 #include "bond.h"
24 #include "angle.h"
25 #include "dihedral.h"
26 #include "improper.h"
27 #include "kspace.h"
28 #include "output.h"
29 #include "update.h"
30 #include "modify.h"
31 #include "timer.h"
32 #include "memory_kokkos.h"
33 #include "kokkos.h"
34
35 using namespace LAMMPS_NS;
36
37 template<class ViewA, class ViewB>
38 struct ForceAdder {
39 ViewA a;
40 ViewB b;
ForceAdderForceAdder41 ForceAdder(const ViewA& a_, const ViewB& b_):a(a_),b(b_) {}
42 KOKKOS_INLINE_FUNCTION
operator ()ForceAdder43 void operator() (const int& i) const {
44 a(i,0) += b(i,0);
45 a(i,1) += b(i,1);
46 a(i,2) += b(i,2);
47 }
48 };
49
50 /* ---------------------------------------------------------------------- */
51
52 template<class View>
53 struct Zero {
54 View v;
ZeroZero55 Zero(const View &v_):v(v_) {}
56 KOKKOS_INLINE_FUNCTION
operator ()Zero57 void operator()(const int &i) const {
58 v(i,0) = 0;
59 v(i,1) = 0;
60 v(i,2) = 0;
61 }
62 };
63
64 /* ---------------------------------------------------------------------- */
65
VerletKokkos(LAMMPS * lmp,int narg,char ** arg)66 VerletKokkos::VerletKokkos(LAMMPS *lmp, int narg, char **arg) :
67 Verlet(lmp, narg, arg)
68 {
69 atomKK = (AtomKokkos *) atom;
70 }
71
72 /* ----------------------------------------------------------------------
73 setup before run
74 ------------------------------------------------------------------------- */
75
setup(int flag)76 void VerletKokkos::setup(int flag)
77 {
78 if (comm->me == 0 && screen) {
79 fputs("Setting up Verlet run ...\n",screen);
80 if (flag) {
81 fmt::print(screen," Unit style : {}\n"
82 " Current step : {}\n"
83 " Time step : {}\n",
84 update->unit_style,update->ntimestep,update->dt);
85 timer->print_timeout(screen);
86 }
87 }
88
89 update->setupflag = 1;
90
91 // setup domain, communication and neighboring
92 // acquire ghosts
93 // build neighbor lists
94
95 lmp->kokkos->auto_sync = 1;
96
97 atom->setup();
98 modify->setup_pre_exchange();
99 if (triclinic) domain->x2lamda(atom->nlocal);
100 domain->pbc();
101 domain->reset_box();
102 comm->setup();
103 if (neighbor->style) neighbor->setup_bins();
104 comm->exchange();
105 if (atom->sortfreq > 0) atom->sort();
106 comm->borders();
107 if (triclinic) domain->lamda2x(atom->nlocal+atom->nghost);
108 domain->image_check();
109 domain->box_too_small_check();
110 modify->setup_pre_neighbor();
111 neighbor->build(1);
112 modify->setup_post_neighbor();
113 neighbor->ncalls = 0;
114
115 // compute all forces
116
117 force->setup();
118 ev_set(update->ntimestep);
119 force_clear();
120 modify->setup_pre_force(vflag);
121
122 if (pair_compute_flag) {
123 atomKK->sync(force->pair->execution_space,force->pair->datamask_read);
124 force->pair->compute(eflag,vflag);
125 atomKK->modified(force->pair->execution_space,force->pair->datamask_modify);
126 }
127 else if (force->pair) force->pair->compute_dummy(eflag,vflag);
128
129 if (atom->molecular != Atom::ATOMIC) {
130 if (force->bond) {
131 atomKK->sync(force->bond->execution_space,force->bond->datamask_read);
132 force->bond->compute(eflag,vflag);
133 atomKK->modified(force->bond->execution_space,force->bond->datamask_modify);
134 }
135 if (force->angle) {
136 atomKK->sync(force->angle->execution_space,force->angle->datamask_read);
137 force->angle->compute(eflag,vflag);
138 atomKK->modified(force->angle->execution_space,force->angle->datamask_modify);
139 }
140 if (force->dihedral) {
141 atomKK->sync(force->dihedral->execution_space,force->dihedral->datamask_read);
142 force->dihedral->compute(eflag,vflag);
143 atomKK->modified(force->dihedral->execution_space,force->dihedral->datamask_modify);
144 }
145 if (force->improper) {
146 atomKK->sync(force->improper->execution_space,force->improper->datamask_read);
147 force->improper->compute(eflag,vflag);
148 atomKK->modified(force->improper->execution_space,force->improper->datamask_modify);
149 }
150 }
151
152 if (force->kspace) {
153 force->kspace->setup();
154 if (kspace_compute_flag) {
155 atomKK->sync(force->kspace->execution_space,force->kspace->datamask_read);
156 force->kspace->compute(eflag,vflag);
157 atomKK->modified(force->kspace->execution_space,force->kspace->datamask_modify);
158 } else force->kspace->compute_dummy(eflag,vflag);
159 }
160
161 modify->setup_pre_reverse(eflag,vflag);
162 if (force->newton) comm->reverse_comm();
163
164 lmp->kokkos->auto_sync = 0;
165 modify->setup(vflag);
166 output->setup(flag);
167 lmp->kokkos->auto_sync = 1;
168 update->setupflag = 0;
169 }
170
171 /* ----------------------------------------------------------------------
172 setup without output
173 flag = 0 = just force calculation
174 flag = 1 = reneighbor and force calculation
175 ------------------------------------------------------------------------- */
176
setup_minimal(int flag)177 void VerletKokkos::setup_minimal(int flag)
178 {
179 update->setupflag = 1;
180
181 // setup domain, communication and neighboring
182 // acquire ghosts
183 // build neighbor lists
184
185 lmp->kokkos->auto_sync = 1;
186
187 if (flag) {
188 modify->setup_pre_exchange();
189 if (triclinic) domain->x2lamda(atom->nlocal);
190 domain->pbc();
191 domain->reset_box();
192 comm->setup();
193 if (neighbor->style) neighbor->setup_bins();
194 comm->exchange();
195 comm->borders();
196 if (triclinic) domain->lamda2x(atom->nlocal+atom->nghost);
197 domain->image_check();
198 domain->box_too_small_check();
199 modify->setup_pre_neighbor();
200 neighbor->build(1);
201 modify->setup_post_neighbor();
202 neighbor->ncalls = 0;
203 }
204
205 // compute all forces
206
207 ev_set(update->ntimestep);
208 force_clear();
209 modify->setup_pre_force(vflag);
210
211 if (pair_compute_flag) {
212 atomKK->sync(force->pair->execution_space,force->pair->datamask_read);
213 force->pair->compute(eflag,vflag);
214 atomKK->modified(force->pair->execution_space,force->pair->datamask_modify);
215 }
216 else if (force->pair) force->pair->compute_dummy(eflag,vflag);
217
218 if (atom->molecular != Atom::ATOMIC) {
219 if (force->bond) {
220 atomKK->sync(force->bond->execution_space,force->bond->datamask_read);
221 force->bond->compute(eflag,vflag);
222 atomKK->modified(force->bond->execution_space,force->bond->datamask_modify);
223 }
224 if (force->angle) {
225 atomKK->sync(force->angle->execution_space,force->angle->datamask_read);
226 force->angle->compute(eflag,vflag);
227 atomKK->modified(force->angle->execution_space,force->angle->datamask_modify);
228 }
229 if (force->dihedral) {
230 atomKK->sync(force->dihedral->execution_space,force->dihedral->datamask_read);
231 force->dihedral->compute(eflag,vflag);
232 atomKK->modified(force->dihedral->execution_space,force->dihedral->datamask_modify);
233 }
234 if (force->improper) {
235 atomKK->sync(force->improper->execution_space,force->improper->datamask_read);
236 force->improper->compute(eflag,vflag);
237 atomKK->modified(force->improper->execution_space,force->improper->datamask_modify);
238 }
239 }
240
241 if (force->kspace) {
242 force->kspace->setup();
243 if (kspace_compute_flag) {
244 atomKK->sync(force->kspace->execution_space,force->kspace->datamask_read);
245 force->kspace->compute(eflag,vflag);
246 atomKK->modified(force->kspace->execution_space,force->kspace->datamask_modify);
247 } else force->kspace->compute_dummy(eflag,vflag);
248 }
249
250 modify->setup_pre_reverse(eflag,vflag);
251 if (force->newton) comm->reverse_comm();
252
253 lmp->kokkos->auto_sync = 0;
254 modify->setup(vflag);
255 lmp->kokkos->auto_sync = 1;
256 update->setupflag = 0;
257 }
258
259 /* ----------------------------------------------------------------------
260 run for N steps
261 ------------------------------------------------------------------------- */
262
run(int n)263 void VerletKokkos::run(int n)
264 {
265 bigint ntimestep;
266 int nflag,sortflag;
267
268 int n_post_integrate = modify->n_post_integrate;
269 int n_pre_exchange = modify->n_pre_exchange;
270 int n_pre_neighbor = modify->n_pre_neighbor;
271 int n_post_neighbor = modify->n_post_neighbor;
272 int n_pre_force = modify->n_pre_force;
273 int n_pre_reverse = modify->n_pre_reverse;
274 int n_post_force = modify->n_post_force;
275 int n_end_of_step = modify->n_end_of_step;
276
277 lmp->kokkos->auto_sync = 0;
278
279 if (atomKK->sortfreq > 0) sortflag = 1;
280 else sortflag = 0;
281
282 f_merge_copy = DAT::t_f_array("VerletKokkos::f_merge_copy",atomKK->k_f.extent(0));
283
284 atomKK->sync(Device,ALL_MASK);
285 //static double time = 0.0;
286 //Kokkos::Impl::Timer ktimer;
287
288 timer->init_timeout();
289 for (int i = 0; i < n; i++) {
290 if (timer->check_timeout(i)) {
291 update->nsteps = i;
292 break;
293 }
294
295 ntimestep = ++update->ntimestep;
296 ev_set(ntimestep);
297
298 // initial time integration
299
300 //ktimer.reset();
301 timer->stamp();
302 modify->initial_integrate(vflag);
303 //time += ktimer.seconds();
304 if (n_post_integrate) modify->post_integrate();
305 timer->stamp(Timer::MODIFY);
306
307 // regular communication vs neighbor list rebuild
308
309 nflag = neighbor->decide();
310
311 if (nflag == 0) {
312 timer->stamp();
313 comm->forward_comm();
314 timer->stamp(Timer::COMM);
315 } else {
316 // added debug
317 //atomKK->sync(Host,ALL_MASK);
318 //atomKK->modified(Host,ALL_MASK);
319
320 if (n_pre_exchange) {
321 timer->stamp();
322 modify->pre_exchange();
323 timer->stamp(Timer::MODIFY);
324 }
325 // debug
326 //atomKK->sync(Host,ALL_MASK);
327 //atomKK->modified(Host,ALL_MASK);
328 if (triclinic) domain->x2lamda(atomKK->nlocal);
329 domain->pbc();
330 if (domain->box_change) {
331 domain->reset_box();
332 comm->setup();
333 if (neighbor->style) neighbor->setup_bins();
334 }
335 timer->stamp();
336
337 // added debug
338 //atomKK->sync(Device,ALL_MASK);
339 //atomKK->modified(Device,ALL_MASK);
340
341 comm->exchange();
342 if (sortflag && ntimestep >= atomKK->nextsort) atomKK->sort();
343 comm->borders();
344
345 // added debug
346 //atomKK->sync(Host,ALL_MASK);
347 //atomKK->modified(Host,ALL_MASK);
348
349 if (triclinic) domain->lamda2x(atomKK->nlocal+atomKK->nghost);
350
351 timer->stamp(Timer::COMM);
352 if (n_pre_neighbor) {
353 modify->pre_neighbor();
354 timer->stamp(Timer::MODIFY);
355 }
356 neighbor->build(1);
357 timer->stamp(Timer::NEIGH);
358 if (n_post_neighbor) {
359 modify->post_neighbor();
360 timer->stamp(Timer::MODIFY);
361 }
362 }
363
364 // force computations
365 // important for pair to come before bonded contributions
366 // since some bonded potentials tally pairwise energy/virial
367 // and Pair:ev_tally() needs to be called before any tallying
368
369 force_clear();
370
371 timer->stamp();
372
373 if (n_pre_force) {
374 modify->pre_force(vflag);
375 timer->stamp(Timer::MODIFY);
376 }
377
378 bool execute_on_host = false;
379 unsigned int datamask_read_device = 0;
380 unsigned int datamask_modify_device = 0;
381 unsigned int datamask_read_host = 0;
382
383 if (pair_compute_flag) {
384 if (force->pair->execution_space==Host) {
385 execute_on_host = true;
386 datamask_read_host |= force->pair->datamask_read;
387 datamask_modify_device |= force->pair->datamask_modify;
388 } else {
389 datamask_read_device |= force->pair->datamask_read;
390 datamask_modify_device |= force->pair->datamask_modify;
391 }
392 }
393 if (atomKK->molecular && force->bond) {
394 if (force->bond->execution_space==Host) {
395 execute_on_host = true;
396 datamask_read_host |= force->bond->datamask_read;
397 datamask_modify_device |= force->bond->datamask_modify;
398 } else {
399 datamask_read_device |= force->bond->datamask_read;
400 datamask_modify_device |= force->bond->datamask_modify;
401 }
402 }
403 if (atomKK->molecular && force->angle) {
404 if (force->angle->execution_space==Host) {
405 execute_on_host = true;
406 datamask_read_host |= force->angle->datamask_read;
407 datamask_modify_device |= force->angle->datamask_modify;
408 } else {
409 datamask_read_device |= force->angle->datamask_read;
410 datamask_modify_device |= force->angle->datamask_modify;
411 }
412 }
413 if (atomKK->molecular && force->dihedral) {
414 if (force->dihedral->execution_space==Host) {
415 execute_on_host = true;
416 datamask_read_host |= force->dihedral->datamask_read;
417 datamask_modify_device |= force->dihedral->datamask_modify;
418 } else {
419 datamask_read_device |= force->dihedral->datamask_read;
420 datamask_modify_device |= force->dihedral->datamask_modify;
421 }
422 }
423 if (atomKK->molecular && force->improper) {
424 if (force->improper->execution_space==Host) {
425 execute_on_host = true;
426 datamask_read_host |= force->improper->datamask_read;
427 datamask_modify_device |= force->improper->datamask_modify;
428 } else {
429 datamask_read_device |= force->improper->datamask_read;
430 datamask_modify_device |= force->improper->datamask_modify;
431 }
432 }
433 if (kspace_compute_flag) {
434 if (force->kspace->execution_space==Host) {
435 execute_on_host = true;
436 datamask_read_host |= force->kspace->datamask_read;
437 datamask_modify_device |= force->kspace->datamask_modify;
438 } else {
439 datamask_read_device |= force->kspace->datamask_read;
440 datamask_modify_device |= force->kspace->datamask_modify;
441 }
442 }
443
444
445 if (pair_compute_flag) {
446 atomKK->sync(force->pair->execution_space,force->pair->datamask_read);
447 atomKK->sync(force->pair->execution_space,~(~force->pair->datamask_read|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
448 Kokkos::Impl::Timer ktimer;
449 force->pair->compute(eflag,vflag);
450 atomKK->modified(force->pair->execution_space,force->pair->datamask_modify);
451 atomKK->modified(force->pair->execution_space,~(~force->pair->datamask_modify|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
452 timer->stamp(Timer::PAIR);
453 }
454
455 if (execute_on_host) {
456 if (pair_compute_flag && force->pair->datamask_modify!=(F_MASK | ENERGY_MASK | VIRIAL_MASK))
457 Kokkos::fence();
458 atomKK->sync_overlapping_device(Host,~(~datamask_read_host|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
459 if (pair_compute_flag && force->pair->execution_space!=Host) {
460 Kokkos::deep_copy(LMPHostType(),atomKK->k_f.h_view,0.0);
461 }
462 }
463
464 if (atomKK->molecular) {
465 if (force->bond) {
466 atomKK->sync(force->bond->execution_space,~(~force->bond->datamask_read|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
467 force->bond->compute(eflag,vflag);
468 atomKK->modified(force->bond->execution_space,~(~force->bond->datamask_modify|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
469 }
470 if (force->angle) {
471 atomKK->sync(force->angle->execution_space,~(~force->angle->datamask_read|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
472 force->angle->compute(eflag,vflag);
473 atomKK->modified(force->angle->execution_space,~(~force->angle->datamask_modify|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
474 }
475 if (force->dihedral) {
476 atomKK->sync(force->dihedral->execution_space,~(~force->dihedral->datamask_read|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
477 force->dihedral->compute(eflag,vflag);
478 atomKK->modified(force->dihedral->execution_space,~(~force->dihedral->datamask_modify|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
479 }
480 if (force->improper) {
481 atomKK->sync(force->improper->execution_space,~(~force->improper->datamask_read|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
482 force->improper->compute(eflag,vflag);
483 atomKK->modified(force->improper->execution_space,~(~force->improper->datamask_modify|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
484 }
485 timer->stamp(Timer::BOND);
486 }
487
488 if (kspace_compute_flag) {
489 atomKK->sync(force->kspace->execution_space,~(~force->kspace->datamask_read|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
490 force->kspace->compute(eflag,vflag);
491 atomKK->modified(force->kspace->execution_space,~(~force->kspace->datamask_modify|(F_MASK | ENERGY_MASK | VIRIAL_MASK)));
492 timer->stamp(Timer::KSPACE);
493 }
494
495 if (execute_on_host && !std::is_same<LMPHostType,LMPDeviceType>::value) {
496 if (f_merge_copy.extent(0)<atomKK->k_f.extent(0)) {
497 f_merge_copy = DAT::t_f_array("VerletKokkos::f_merge_copy",atomKK->k_f.extent(0));
498 }
499 f = atomKK->k_f.d_view;
500 Kokkos::deep_copy(LMPHostType(),f_merge_copy,atomKK->k_f.h_view);
501 Kokkos::parallel_for(atomKK->k_f.extent(0),
502 ForceAdder<DAT::t_f_array,DAT::t_f_array>(atomKK->k_f.d_view,f_merge_copy));
503 atomKK->k_f.clear_sync_state(); // special case
504 atomKK->k_f.modify<LMPDeviceType>();
505 }
506
507 if (n_pre_reverse) {
508 modify->pre_reverse(eflag,vflag);
509 timer->stamp(Timer::MODIFY);
510 }
511
512 // reverse communication of forces
513
514 if (force->newton) {
515 Kokkos::fence();
516 comm->reverse_comm();
517 timer->stamp(Timer::COMM);
518 }
519
520 // force modifications, final time integration, diagnostics
521
522 if (n_post_force) modify->post_force(vflag);
523 modify->final_integrate();
524 if (n_end_of_step) modify->end_of_step();
525 timer->stamp(Timer::MODIFY);
526
527 // all output
528
529 if (ntimestep == output->next) {
530 atomKK->sync(Host,ALL_MASK);
531
532 timer->stamp();
533 output->write(ntimestep);
534 timer->stamp(Timer::OUTPUT);
535 }
536 }
537
538 atomKK->sync(Host,ALL_MASK);
539 lmp->kokkos->auto_sync = 1;
540 }
541
542 /* ----------------------------------------------------------------------
543 clear force on own & ghost atoms
544 clear other arrays as needed
545 ------------------------------------------------------------------------- */
546
force_clear()547 void VerletKokkos::force_clear()
548 {
549 if (external_force_clear) return;
550
551 atomKK->k_f.clear_sync_state(); // ignore host forces/torques since device views
552 atomKK->k_torque.clear_sync_state(); // will be cleared below
553
554 // clear force on all particles
555 // if either newton flag is set, also include ghosts
556 // when using threads always clear all forces.
557
558 if (neighbor->includegroup == 0) {
559 int nall = atomKK->nlocal;
560 if (force->newton) nall += atomKK->nghost;
561
562 Kokkos::parallel_for(nall, Zero<typename ArrayTypes<LMPDeviceType>::t_f_array>(atomKK->k_f.view<LMPDeviceType>()));
563 atomKK->modified(Device,F_MASK);
564
565 if (torqueflag) {
566 Kokkos::parallel_for(nall, Zero<typename ArrayTypes<LMPDeviceType>::t_f_array>(atomKK->k_torque.view<LMPDeviceType>()));
567 atomKK->modified(Device,TORQUE_MASK);
568 }
569
570 // reset SPIN forces
571
572 if (extraflag) {
573 Kokkos::parallel_for(nall, Zero<typename ArrayTypes<LMPDeviceType>::t_fm_array>(atomKK->k_fm.view<LMPDeviceType>()));
574 atomKK->modified(Device,FM_MASK);
575 Kokkos::parallel_for(nall, Zero<typename ArrayTypes<LMPDeviceType>::t_fm_array>(atomKK->k_fm_long.view<LMPDeviceType>()));
576 atomKK->modified(Device,FML_MASK);
577 }
578
579 // neighbor includegroup flag is set
580 // clear force only on initial nfirst particles
581 // if either newton flag is set, also include ghosts
582
583 } else {
584 Kokkos::parallel_for(atomKK->nfirst, Zero<typename ArrayTypes<LMPDeviceType>::t_f_array>(atomKK->k_f.view<LMPDeviceType>()));
585 atomKK->modified(Device,F_MASK);
586
587 if (torqueflag) {
588 Kokkos::parallel_for(atomKK->nfirst, Zero<typename ArrayTypes<LMPDeviceType>::t_f_array>(atomKK->k_torque.view<LMPDeviceType>()));
589 atomKK->modified(Device,TORQUE_MASK);
590 }
591
592 // reset SPIN forces
593
594 if (extraflag) {
595 Kokkos::parallel_for(atomKK->nfirst, Zero<typename ArrayTypes<LMPDeviceType>::t_fm_array>(atomKK->k_fm.view<LMPDeviceType>()));
596 atomKK->modified(Device,FM_MASK);
597 Kokkos::parallel_for(atomKK->nfirst, Zero<typename ArrayTypes<LMPDeviceType>::t_fm_array>(atomKK->k_fm_long.view<LMPDeviceType>()));
598 atomKK->modified(Device,FML_MASK);
599 }
600
601 if (force->newton) {
602 auto range = Kokkos::RangePolicy<LMPDeviceType>(atomKK->nlocal, atomKK->nlocal + atomKK->nghost);
603 Kokkos::parallel_for(range, Zero<typename ArrayTypes<LMPDeviceType>::t_f_array>(atomKK->k_f.view<LMPDeviceType>()));
604 atomKK->modified(Device,F_MASK);
605
606 if (torqueflag) {
607 Kokkos::parallel_for(range, Zero<typename ArrayTypes<LMPDeviceType>::t_f_array>(atomKK->k_torque.view<LMPDeviceType>()));
608 atomKK->modified(Device,TORQUE_MASK);
609 }
610
611 // reset SPIN forces
612
613 if (extraflag) {
614 Kokkos::parallel_for(range, Zero<typename ArrayTypes<LMPDeviceType>::t_fm_array>(atomKK->k_fm.view<LMPDeviceType>()));
615 atomKK->modified(Device,FM_MASK);
616 Kokkos::parallel_for(range, Zero<typename ArrayTypes<LMPDeviceType>::t_fm_array>(atomKK->k_fm_long.view<LMPDeviceType>()));
617 atomKK->modified(Device,FML_MASK);
618 }
619 }
620 }
621 }
622