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 "npair.h"
16 #include <cmath>
17 #include "neighbor.h"
18 #include "neigh_request.h"
19 #include "nbin.h"
20 #include "nstencil.h"
21 #include "atom.h"
22 #include "update.h"
23 #include "memory.h"
24 #include "error.h"
25
26 using namespace LAMMPS_NS;
27
28 /* ---------------------------------------------------------------------- */
29
NPair(LAMMPS * lmp)30 NPair::NPair(LAMMPS *lmp)
31 : Pointers(lmp), nb(nullptr), ns(nullptr), bins(nullptr), stencil(nullptr)
32 {
33 last_build = -1;
34 mycutneighsq = nullptr;
35 molecular = atom->molecular;
36 copymode = 0;
37 execution_space = Host;
38 }
39
40 /* ---------------------------------------------------------------------- */
41
~NPair()42 NPair::~NPair()
43 {
44 if (copymode) return;
45
46 memory->destroy(mycutneighsq);
47 }
48
49 /* ---------------------------------------------------------------------- */
50
post_constructor(NeighRequest * nrq)51 void NPair::post_constructor(NeighRequest *nrq)
52 {
53 cutoff_custom = 0.0;
54 if (nrq->cut) cutoff_custom = nrq->cutoff;
55 }
56
57 /* ----------------------------------------------------------------------
58 copy needed info from Neighbor class to this build class
59 done once per run
60 ------------------------------------------------------------------------- */
61
copy_neighbor_info()62 void NPair::copy_neighbor_info()
63 {
64 // general params
65
66 includegroup = neighbor->includegroup;
67 exclude = neighbor->exclude;
68 skin = neighbor->skin;
69 cutneighsq = neighbor->cutneighsq;
70 cutneighghostsq = neighbor->cutneighghostsq;
71 cut_inner_sq = neighbor->cut_inner_sq;
72 cut_middle_sq = neighbor->cut_middle_sq;
73 cut_middle_inside_sq = neighbor->cut_middle_inside_sq;
74 bboxlo = neighbor->bboxlo;
75 bboxhi = neighbor->bboxhi;
76
77 // exclusion info
78
79 nex_type = neighbor->nex_type;
80 ex1_type = neighbor->ex1_type;
81 ex2_type = neighbor->ex2_type;
82 ex_type = neighbor->ex_type;
83
84 nex_group = neighbor->nex_group;
85 ex1_group = neighbor->ex1_group;
86 ex2_group = neighbor->ex2_group;
87 ex1_bit = neighbor->ex1_bit;
88 ex2_bit = neighbor->ex2_bit;
89
90 nex_mol = neighbor->nex_mol;
91 ex_mol_group = neighbor->ex_mol_group;
92 ex_mol_bit = neighbor->ex_mol_bit;
93 ex_mol_intra = neighbor->ex_mol_intra;
94
95 // special info
96
97 special_flag = neighbor->special_flag;
98
99 // multi info
100
101 ncollections = neighbor->ncollections;
102 cutcollectionsq = neighbor->cutcollectionsq;
103
104 // overwrite per-type Neighbor cutoffs with custom value set by requestor
105 // only works for style = BIN (checked by Neighbor class)
106
107 if (cutoff_custom > 0.0) {
108 memory->destroy(mycutneighsq);
109 int n = atom->ntypes;
110 memory->create(mycutneighsq,n+1,n+1,"npair:cutneighsq");
111 int i,j;
112 for (i = 1; i <= n; i++)
113 for (j = 1; j <= n; j++)
114 mycutneighsq[i][j] = cutoff_custom * cutoff_custom;
115 cutneighsq = mycutneighsq;
116 }
117 }
118
119 /* ----------------------------------------------------------------------
120 copy info from NBin class to this build class
121 ------------------------------------------------------------------------- */
122
copy_bin_info()123 void NPair::copy_bin_info()
124 {
125 nbinx = nb->nbinx;
126 nbiny = nb->nbiny;
127 nbinz = nb->nbinz;
128 mbins = nb->mbins;
129 mbinx = nb->mbinx;
130 mbiny = nb->mbiny;
131 mbinz = nb->mbinz;
132 mbinxlo = nb->mbinxlo;
133 mbinylo = nb->mbinylo;
134 mbinzlo = nb->mbinzlo;
135
136 bininvx = nb->bininvx;
137 bininvy = nb->bininvy;
138 bininvz = nb->bininvz;
139
140 atom2bin = nb->atom2bin;
141 bins = nb->bins;
142 binhead = nb->binhead;
143
144 nbinx_multi = nb->nbinx_multi;
145 nbiny_multi = nb->nbiny_multi;
146 nbinz_multi = nb->nbinz_multi;
147 mbins_multi = nb->mbins_multi;
148 mbinx_multi = nb->mbinx_multi;
149 mbiny_multi = nb->mbiny_multi;
150 mbinz_multi = nb->mbinz_multi;
151 mbinxlo_multi = nb->mbinxlo_multi;
152 mbinylo_multi = nb->mbinylo_multi;
153 mbinzlo_multi = nb->mbinzlo_multi;
154
155 bininvx_multi = nb->bininvx_multi;
156 bininvy_multi = nb->bininvy_multi;
157 bininvz_multi = nb->bininvz_multi;
158
159 binhead_multi = nb->binhead_multi;
160 }
161
162 /* ----------------------------------------------------------------------
163 copy info from NStencil class to this build class
164 ------------------------------------------------------------------------- */
165
copy_stencil_info()166 void NPair::copy_stencil_info()
167 {
168 nstencil = ns->nstencil;
169 stencil = ns->stencil;
170 stencilxyz = ns->stencilxyz;
171 nstencil_multi_old = ns->nstencil_multi_old;
172 stencil_multi_old = ns->stencil_multi_old;
173 distsq_multi_old = ns->distsq_multi_old;
174
175 nstencil_multi = ns->nstencil_multi;
176 stencil_multi = ns->stencil_multi;
177 }
178
179 /* ----------------------------------------------------------------------
180 copy info from NBin and NStencil classes to this build class
181 ------------------------------------------------------------------------- */
182
build_setup()183 void NPair::build_setup()
184 {
185 if (nb) copy_bin_info();
186 if (ns) copy_stencil_info();
187
188 // set here, since build_setup() always called before build()
189 last_build = update->ntimestep;
190 }
191
192 /* ----------------------------------------------------------------------
193 test if atom pair i,j is excluded from neighbor list
194 due to type, group, molecule settings from neigh_modify command
195 return 1 if should be excluded, 0 if included
196 ------------------------------------------------------------------------- */
197
exclusion(int i,int j,int itype,int jtype,int * mask,tagint * molecule) const198 int NPair::exclusion(int i, int j, int itype, int jtype,
199 int *mask, tagint *molecule) const {
200 int m;
201
202 if (nex_type && ex_type[itype][jtype]) return 1;
203
204 if (nex_group) {
205 for (m = 0; m < nex_group; m++) {
206 if (mask[i] & ex1_bit[m] && mask[j] & ex2_bit[m]) return 1;
207 if (mask[i] & ex2_bit[m] && mask[j] & ex1_bit[m]) return 1;
208 }
209 }
210
211 if (nex_mol) {
212 for (m = 0; m < nex_mol; m++)
213
214 // intra-chain: exclude i-j pair if in same molecule
215 // inter-chain: exclude i-j pair if in different molecules
216
217 if (ex_mol_intra[m]) {
218 if (mask[i] & ex_mol_bit[m] && mask[j] & ex_mol_bit[m] &&
219 molecule[i] == molecule[j]) return 1;
220 } else {
221 if (mask[i] & ex_mol_bit[m] && mask[j] & ex_mol_bit[m] &&
222 molecule[i] != molecule[j]) return 1;
223 }
224 }
225
226 return 0;
227 }
228
229 /* ----------------------------------------------------------------------
230 same as coord2bin in Nbin, but also return ix,iy,iz offsets in each dim
231 used by some of the ghost neighbor lists
232 ------------------------------------------------------------------------- */
233
coord2bin(double * x,int & ix,int & iy,int & iz)234 int NPair::coord2bin(double *x, int &ix, int &iy, int &iz)
235 {
236 if (!std::isfinite(x[0]) || !std::isfinite(x[1]) || !std::isfinite(x[2]))
237 error->one(FLERR,"Non-numeric positions - simulation unstable");
238
239 if (x[0] >= bboxhi[0])
240 ix = static_cast<int> ((x[0]-bboxhi[0])*bininvx) + nbinx;
241 else if (x[0] >= bboxlo[0]) {
242 ix = static_cast<int> ((x[0]-bboxlo[0])*bininvx);
243 ix = MIN(ix,nbinx-1);
244 } else
245 ix = static_cast<int> ((x[0]-bboxlo[0])*bininvx) - 1;
246
247 if (x[1] >= bboxhi[1])
248 iy = static_cast<int> ((x[1]-bboxhi[1])*bininvy) + nbiny;
249 else if (x[1] >= bboxlo[1]) {
250 iy = static_cast<int> ((x[1]-bboxlo[1])*bininvy);
251 iy = MIN(iy,nbiny-1);
252 } else
253 iy = static_cast<int> ((x[1]-bboxlo[1])*bininvy) - 1;
254
255 if (x[2] >= bboxhi[2])
256 iz = static_cast<int> ((x[2]-bboxhi[2])*bininvz) + nbinz;
257 else if (x[2] >= bboxlo[2]) {
258 iz = static_cast<int> ((x[2]-bboxlo[2])*bininvz);
259 iz = MIN(iz,nbinz-1);
260 } else
261 iz = static_cast<int> ((x[2]-bboxlo[2])*bininvz) - 1;
262
263 ix -= mbinxlo;
264 iy -= mbinylo;
265 iz -= mbinzlo;
266 return iz*mbiny*mbinx + iy*mbinx + ix;
267 }
268
269
270 /* ----------------------------------------------------------------------
271 multi version of coord2bin for a given collection
272 ------------------------------------------------------------------------- */
273
coord2bin(double * x,int ic)274 int NPair::coord2bin(double *x, int ic)
275 {
276 int ix,iy,iz;
277 int ibin;
278
279 if (!std::isfinite(x[0]) || !std::isfinite(x[1]) || !std::isfinite(x[2]))
280 error->one(FLERR,"Non-numeric positions - simulation unstable");
281
282 if (x[0] >= bboxhi[0])
283 ix = static_cast<int> ((x[0]-bboxhi[0])*bininvx_multi[ic]) + nbinx_multi[ic];
284 else if (x[0] >= bboxlo[0]) {
285 ix = static_cast<int> ((x[0]-bboxlo[0])*bininvx_multi[ic]);
286 ix = MIN(ix,nbinx_multi[ic]-1);
287 } else
288 ix = static_cast<int> ((x[0]-bboxlo[0])*bininvx_multi[ic]) - 1;
289
290 if (x[1] >= bboxhi[1])
291 iy = static_cast<int> ((x[1]-bboxhi[1])*bininvy_multi[ic]) + nbiny_multi[ic];
292 else if (x[1] >= bboxlo[1]) {
293 iy = static_cast<int> ((x[1]-bboxlo[1])*bininvy_multi[ic]);
294 iy = MIN(iy,nbiny_multi[ic]-1);
295 } else
296 iy = static_cast<int> ((x[1]-bboxlo[1])*bininvy_multi[ic]) - 1;
297
298 if (x[2] >= bboxhi[2])
299 iz = static_cast<int> ((x[2]-bboxhi[2])*bininvz_multi[ic]) + nbinz_multi[ic];
300 else if (x[2] >= bboxlo[2]) {
301 iz = static_cast<int> ((x[2]-bboxlo[2])*bininvz_multi[ic]);
302 iz = MIN(iz,nbinz_multi[ic]-1);
303 } else
304 iz = static_cast<int> ((x[2]-bboxlo[2])*bininvz_multi[ic]) - 1;
305
306 ix -= mbinxlo_multi[ic];
307 iy -= mbinylo_multi[ic];
308 iz -= mbinzlo_multi[ic];
309 ibin = iz*mbiny_multi[ic]*mbinx_multi[ic] + iy*mbinx_multi[ic] + ix;
310 return ibin;
311 }
312