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 author: W. Michael Brown (Intel)
17 ------------------------------------------------------------------------- */
18
19 #include "npair_half_bin_newton_intel.h"
20
21 #include "atom.h"
22 #include "comm.h"
23 #include "error.h"
24 #include "neigh_list.h"
25 #include "neighbor.h"
26
27 using namespace LAMMPS_NS;
28
29 /* ---------------------------------------------------------------------- */
30
NPairHalfBinNewtonIntel(LAMMPS * lmp)31 NPairHalfBinNewtonIntel::NPairHalfBinNewtonIntel(LAMMPS *lmp) :
32 NPairIntel(lmp) {}
33
34 /* ----------------------------------------------------------------------
35 binned neighbor list construction with full Newton's 3rd law
36 each owned atom i checks its own bin and other bins in Newton stencil
37 every pair stored exactly once by some processor
38 ------------------------------------------------------------------------- */
39
build(NeighList * list)40 void NPairHalfBinNewtonIntel::build(NeighList *list)
41 {
42 if (nstencil / 2 > INTEL_MAX_STENCIL_CHECK)
43 error->all(FLERR, "Too many neighbor bins for INTEL package.");
44
45 #ifdef _LMP_INTEL_OFFLOAD
46 if (exclude)
47 error->all(FLERR, "Exclusion lists not yet supported for Intel offload");
48 #endif
49
50 if (_fix->precision() == FixIntel::PREC_MODE_MIXED)
51 hbni(list, _fix->get_mixed_buffers());
52 else if (_fix->precision() == FixIntel::PREC_MODE_DOUBLE)
53 hbni(list, _fix->get_double_buffers());
54 else
55 hbni(list, _fix->get_single_buffers());
56
57 _fix->stop_watch(TIME_HOST_NEIGHBOR);
58 }
59
60 template <class flt_t, class acc_t>
61 void NPairHalfBinNewtonIntel::
hbni(NeighList * list,IntelBuffers<flt_t,acc_t> * buffers)62 hbni(NeighList *list, IntelBuffers<flt_t,acc_t> *buffers) {
63 const int nlocal = (includegroup) ? atom->nfirst : atom->nlocal;
64 list->inum = nlocal;
65
66 int host_start = _fix->host_start_neighbor();
67 const int off_end = _fix->offload_end_neighbor();
68
69 #ifdef _LMP_INTEL_OFFLOAD
70 if (off_end) grow_stencil();
71 if (_fix->full_host_list()) host_start = 0;
72 int offload_noghost = _fix->offload_noghost();
73 #endif
74
75 buffers->grow_list(list, atom->nlocal, comm->nthreads, 0, off_end);
76
77 int need_ic = 0;
78 if (atom->molecular != Atom::ATOMIC)
79 dminimum_image_check(need_ic, neighbor->cutneighmax, neighbor->cutneighmax,
80 neighbor->cutneighmax);
81
82 #ifdef _LMP_INTEL_OFFLOAD
83 if (need_ic) {
84 if (offload_noghost) {
85 bin_newton<flt_t,acc_t,1,1,0,0,0>(1, list, buffers, 0, off_end);
86 bin_newton<flt_t,acc_t,1,1,0,0,0>(0, list, buffers, host_start, nlocal,
87 off_end);
88 } else {
89 bin_newton<flt_t,acc_t,0,1,0,0,0>(1, list, buffers, 0, off_end);
90 bin_newton<flt_t,acc_t,0,1,0,0,0>(0, list, buffers, host_start, nlocal);
91 }
92 } else {
93 if (offload_noghost) {
94 bin_newton<flt_t,acc_t,1,0,0,0,0>(1, list, buffers, 0, off_end);
95 bin_newton<flt_t,acc_t,1,0,0,0,0>(0, list, buffers, host_start, nlocal,
96 off_end);
97 } else {
98 bin_newton<flt_t,acc_t,0,0,0,0,0>(1, list, buffers, 0, off_end);
99 bin_newton<flt_t,acc_t,0,0,0,0,0>(0, list, buffers, host_start, nlocal);
100 }
101 }
102 #else
103 if (need_ic)
104 bin_newton<flt_t,acc_t,0,1,0,0,0>(0, list, buffers, host_start, nlocal);
105 else
106 bin_newton<flt_t,acc_t,0,0,0,0,0>(0, list, buffers, host_start, nlocal);
107 #endif
108 }
109