1 | |
2 | |
3 | |
4 | |
5 | |
6 | |
7 | |
8 | |
9 | |
10 | |
11 | |
12 | |
13 | |
14 | |
15 | |
16 | |
17 | |
18 | |
19 | |
20 | |
21 | |
22 | |
23 | |
24 | |
25 | |
26 | |
27 | |
28 | |
29 | |
30 | |
31 | |
32 | |
33 | |
34 | |
35 | |
36 | #define UNROLLI NBNXN_CPU_CLUSTER_I_SIZE4 |
37 | #define UNROLLJ NBNXN_CPU_CLUSTER_I_SIZE4 |
38 | |
39 | |
40 | #define X_STRIDE 3 |
41 | #define F_STRIDE 3 |
42 | |
43 | #define XI_STRIDE 3 |
44 | #define FI_STRIDE 3 |
45 | |
46 | |
47 | |
48 | |
49 | |
50 | |
51 | |
52 | |
53 | #define CALC_SHIFTFORCES |
54 | |
55 | #ifdef CALC_COUL_RF |
56 | #define NBK_FUNC_NAME2(ljt, feg) nbnxn_kernel ## _ElecRF ## ljt ## feg ## _ref |
57 | #endif |
58 | #ifdef CALC_COUL_TAB |
59 | #ifndef VDW_CUTOFF_CHECK |
60 | #define NBK_FUNC_NAME2(ljt, feg) nbnxn_kernel ## _ElecQSTab ## ljt ## feg ## _ref |
61 | #else |
62 | #define NBK_FUNC_NAME2(ljt, feg) nbnxn_kernel ## _ElecQSTabTwinCut ## ljt ## feg ## _ref |
63 | #endif |
64 | #endif |
65 | |
66 | #if defined LJ_CUT && !defined LJ_EWALD |
67 | #define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJ, feg) |
68 | #elif defined LJ_FORCE_SWITCH |
69 | #define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJFsw, feg) |
70 | #elif defined LJ_POT_SWITCH |
71 | #define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJPsw, feg) |
72 | #elif defined LJ_EWALD |
73 | #ifdef LJ_EWALD_COMB_GEOM |
74 | #define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJEwCombGeom, feg) |
75 | #else |
76 | #define NBK_FUNC_NAME(feg) NBK_FUNC_NAME2(_VdwLJEwCombLB, feg) |
77 | #endif |
78 | #else |
79 | #error "No VdW type defined" |
80 | #endif |
81 | |
82 | static void |
83 | #ifndef CALC_ENERGIES |
84 | NBK_FUNC_NAME(_F) |
85 | #else |
86 | #ifndef ENERGY_GROUPS |
87 | NBK_FUNC_NAME(_VF) |
88 | #else |
89 | NBK_FUNC_NAME(_VgrpF) |
90 | #endif |
91 | #endif |
92 | #undef NBK_FUNC_NAME |
93 | #undef NBK_FUNC_NAME2 |
94 | (const nbnxn_pairlist_t *nbl, |
95 | const nbnxn_atomdata_t *nbat, |
96 | const interaction_const_t *ic, |
97 | rvec *shift_vec, |
98 | real *f |
99 | #ifdef CALC_SHIFTFORCES |
100 | , |
101 | real *fshift |
102 | #endif |
103 | #ifdef CALC_ENERGIES |
104 | , |
105 | real *Vvdw, |
106 | real *Vc |
107 | #endif |
108 | ) |
109 | { |
110 | const nbnxn_ci_t *nbln; |
111 | const nbnxn_cj_t *l_cj; |
112 | const int *type; |
113 | const real *q; |
114 | const real *shiftvec; |
115 | const real *x; |
116 | const real *nbfp; |
117 | real rcut2; |
118 | #ifdef VDW_CUTOFF_CHECK |
119 | real rvdw2; |
120 | #endif |
121 | int ntype2; |
122 | real facel; |
123 | real *nbfp_i; |
124 | int n, ci, ci_sh; |
125 | int ish, ishf; |
126 | gmx_bool do_LJ, half_LJ, do_coul, do_self; |
127 | int cjind0, cjind1, cjind; |
128 | int ip, jp; |
129 | |
130 | real xi[UNROLLI*XI_STRIDE]; |
131 | real fi[UNROLLI*FI_STRIDE]; |
132 | real qi[UNROLLI]; |
133 | |
134 | #ifdef CALC_ENERGIES |
135 | #ifndef ENERGY_GROUPS |
136 | |
137 | real Vvdw_ci, Vc_ci; |
138 | #else |
139 | int egp_mask; |
140 | int egp_sh_i[UNROLLI]; |
141 | #endif |
142 | #endif |
143 | #ifdef LJ_POT_SWITCH |
144 | real swV3, swV4, swV5; |
145 | real swF2, swF3, swF4; |
146 | #endif |
147 | #ifdef LJ_EWALD |
148 | real lje_coeff2, lje_coeff6_6, lje_vc; |
149 | const real *ljc; |
150 | #endif |
151 | |
152 | #ifdef CALC_COUL_RF |
153 | real k_rf2; |
154 | #ifdef CALC_ENERGIES |
155 | real k_rf, c_rf; |
156 | #endif |
157 | #endif |
158 | #ifdef CALC_COUL_TAB |
159 | real tabscale; |
160 | #ifdef CALC_ENERGIES |
161 | real halfsp; |
162 | #endif |
163 | #ifndef GMX_DOUBLE |
164 | const real *tab_coul_FDV0; |
165 | #else |
166 | const real *tab_coul_F; |
167 | const real *tab_coul_V; |
168 | #endif |
169 | #endif |
170 | |
171 | int ninner; |
172 | |
173 | #ifdef COUNT_PAIRS |
174 | int npair = 0; |
175 | #endif |
176 | |
177 | #ifdef LJ_POT_SWITCH |
178 | swV3 = ic->vdw_switch.c3; |
179 | swV4 = ic->vdw_switch.c4; |
180 | swV5 = ic->vdw_switch.c5; |
181 | swF2 = 3*ic->vdw_switch.c3; |
182 | swF3 = 4*ic->vdw_switch.c4; |
183 | swF4 = 5*ic->vdw_switch.c5; |
184 | #endif |
185 | |
186 | #ifdef LJ_EWALD |
187 | lje_coeff2 = ic->ewaldcoeff_lj*ic->ewaldcoeff_lj; |
188 | lje_coeff6_6 = lje_coeff2*lje_coeff2*lje_coeff2/6.0; |
189 | lje_vc = ic->sh_lj_ewald; |
190 | |
191 | ljc = nbat->nbfp_comb; |
192 | #endif |
193 | |
194 | #ifdef CALC_COUL_RF |
195 | k_rf2 = 2*ic->k_rf; |
196 | #ifdef CALC_ENERGIES |
197 | k_rf = ic->k_rf; |
198 | c_rf = ic->c_rf; |
199 | #endif |
200 | #endif |
201 | #ifdef CALC_COUL_TAB |
202 | tabscale = ic->tabq_scale; |
203 | #ifdef CALC_ENERGIES |
204 | halfsp = 0.5/ic->tabq_scale; |
205 | #endif |
206 | |
207 | #ifndef GMX_DOUBLE |
208 | tab_coul_FDV0 = ic->tabq_coul_FDV0; |
209 | #else |
210 | tab_coul_F = ic->tabq_coul_F; |
211 | tab_coul_V = ic->tabq_coul_V; |
212 | #endif |
213 | #endif |
214 | |
215 | #ifdef ENERGY_GROUPS |
216 | egp_mask = (1<<nbat->neg_2log) - 1; |
217 | #endif |
218 | |
219 | |
220 | rcut2 = ic->rcoulomb*ic->rcoulomb; |
221 | #ifdef VDW_CUTOFF_CHECK |
222 | rvdw2 = ic->rvdw*ic->rvdw; |
223 | #endif |
224 | |
225 | ntype2 = nbat->ntype*2; |
226 | nbfp = nbat->nbfp; |
227 | q = nbat->q; |
228 | type = nbat->type; |
229 | facel = ic->epsfac; |
230 | shiftvec = shift_vec[0]; |
231 | x = nbat->x; |
232 | |
233 | l_cj = nbl->cj; |
234 | |
235 | ninner = 0; |
236 | for (n = 0; n < nbl->nci; n++) |
237 | { |
238 | int i, d; |
239 | |
240 | nbln = &nbl->ci[n]; |
241 | |
242 | ish = (nbln->shift & NBNXN_CI_SHIFT127); |
243 | |
244 | ishf = ish*DIM3; |
245 | cjind0 = nbln->cj_ind_start; |
246 | cjind1 = nbln->cj_ind_end; |
247 | |
248 | ci = nbln->ci; |
249 | ci_sh = (ish == CENTRAL(((2*1 +1)*(2*1 +1)*(2*2 +1))/2) ? ci : -1); |
250 | |
251 | |
252 | |
253 | |
254 | |
255 | |
256 | |
257 | do_LJ = (nbln->shift & NBNXN_CI_DO_LJ(0)(1<<(7+3*(0)))); |
258 | do_coul = (nbln->shift & NBNXN_CI_DO_COUL(0)(1<<(9+3*(0)))); |
259 | half_LJ = ((nbln->shift & NBNXN_CI_HALF_LJ(0)(1<<(8+3*(0)))) || !do_LJ) && do_coul; |
260 | #ifdef LJ_EWALD |
261 | do_self = TRUE1; |
262 | #else |
263 | do_self = do_coul; |
| Value stored to 'do_self' is never read |
264 | #endif |
265 | |
266 | #ifdef CALC_ENERGIES |
267 | #ifndef ENERGY_GROUPS |
268 | Vvdw_ci = 0; |
269 | Vc_ci = 0; |
270 | #else |
271 | for (i = 0; i < UNROLLI; i++) |
272 | { |
273 | egp_sh_i[i] = ((nbat->energrp[ci]>>(i*nbat->neg_2log)) & egp_mask)*nbat->nenergrp; |
274 | } |
275 | #endif |
276 | #endif |
277 | |
278 | for (i = 0; i < UNROLLI; i++) |
279 | { |
280 | for (d = 0; d < DIM3; d++) |
281 | { |
282 | xi[i*XI_STRIDE+d] = x[(ci*UNROLLI+i)*X_STRIDE+d] + shiftvec[ishf+d]; |
283 | fi[i*FI_STRIDE+d] = 0; |
284 | } |
285 | |
286 | qi[i] = facel*q[ci*UNROLLI+i]; |
287 | } |
288 | |
289 | #ifdef CALC_ENERGIES |
290 | if (do_self) |
291 | { |
292 | real Vc_sub_self; |
293 | |
294 | #ifdef CALC_COUL_RF |
295 | Vc_sub_self = 0.5*c_rf; |
296 | #endif |
297 | #ifdef CALC_COUL_TAB |
298 | #ifdef GMX_DOUBLE |
299 | Vc_sub_self = 0.5*tab_coul_V[0]; |
300 | #else |
301 | Vc_sub_self = 0.5*tab_coul_FDV0[2]; |
302 | #endif |
303 | #endif |
304 | |
305 | if (l_cj[nbln->cj_ind_start].cj == ci_sh) |
306 | { |
307 | for (i = 0; i < UNROLLI; i++) |
308 | { |
309 | int egp_ind; |
310 | #ifdef ENERGY_GROUPS |
311 | egp_ind = egp_sh_i[i] + ((nbat->energrp[ci]>>(i*nbat->neg_2log)) & egp_mask); |
312 | #else |
313 | egp_ind = 0; |
314 | #endif |
315 | |
316 | Vc[egp_ind] -= qi[i]*q[ci*UNROLLI+i]*Vc_sub_self; |
317 | |
318 | #ifdef LJ_EWALD |
319 | |
320 | Vvdw[egp_ind] += 0.5*nbat->nbfp[nbat->type[ci*UNROLLI+i]*(nbat->ntype + 1)*2]/6*lje_coeff6_6; |
321 | #endif |
322 | } |
323 | } |
324 | } |
325 | #endif /* CALC_ENERGIES */ |
326 | |
327 | cjind = cjind0; |
328 | while (cjind < cjind1 && nbl->cj[cjind].excl != 0xffff) |
329 | { |
330 | #define CHECK_EXCLS |
331 | if (half_LJ) |
332 | { |
333 | #define CALC_COULOMB |
334 | #define HALF_LJ |
335 | #include "nbnxn_kernel_ref_inner.h" |
336 | #undef HALF_LJ |
337 | #undef CALC_COULOMB |
338 | } |
339 | |
340 | else if (do_coul) |
341 | { |
342 | #define CALC_COULOMB |
343 | #include "nbnxn_kernel_ref_inner.h" |
344 | #undef CALC_COULOMB |
345 | } |
346 | else |
347 | { |
348 | #include "nbnxn_kernel_ref_inner.h" |
349 | } |
350 | #undef CHECK_EXCLS |
351 | cjind++; |
352 | } |
353 | |
354 | for (; (cjind < cjind1); cjind++) |
355 | { |
356 | if (half_LJ) |
357 | { |
358 | #define CALC_COULOMB |
359 | #define HALF_LJ |
360 | #include "nbnxn_kernel_ref_inner.h" |
361 | #undef HALF_LJ |
362 | #undef CALC_COULOMB |
363 | } |
364 | |
365 | else if (do_coul) |
366 | { |
367 | #define CALC_COULOMB |
368 | #include "nbnxn_kernel_ref_inner.h" |
369 | #undef CALC_COULOMB |
370 | } |
371 | else |
372 | { |
373 | #include "nbnxn_kernel_ref_inner.h" |
374 | } |
375 | } |
376 | ninner += cjind1 - cjind0; |
377 | |
378 | |
379 | for (i = 0; i < UNROLLI; i++) |
380 | { |
381 | for (d = 0; d < DIM3; d++) |
382 | { |
383 | f[(ci*UNROLLI+i)*F_STRIDE+d] += fi[i*FI_STRIDE+d]; |
384 | } |
385 | } |
386 | #ifdef CALC_SHIFTFORCES |
387 | if (fshift != NULL((void*)0)) |
388 | { |
389 | |
390 | for (i = 0; i < UNROLLI; i++) |
391 | { |
392 | for (d = 0; d < DIM3; d++) |
393 | { |
394 | fshift[ishf+d] += fi[i*FI_STRIDE+d]; |
395 | } |
396 | } |
397 | } |
398 | #endif |
399 | |
400 | #ifdef CALC_ENERGIES |
401 | #ifndef ENERGY_GROUPS |
402 | *Vvdw += Vvdw_ci; |
403 | *Vc += Vc_ci; |
404 | #endif |
405 | #endif |
406 | } |
407 | |
408 | #ifdef COUNT_PAIRS |
409 | printf("atom pairs %d\n", npair); |
410 | #endif |
411 | } |
412 | |
413 | #undef CALC_SHIFTFORCES |
414 | |
415 | #undef X_STRIDE |
416 | #undef F_STRIDE |
417 | #undef XI_STRIDE |
418 | #undef FI_STRIDE |
419 | |
420 | #undef UNROLLI |
421 | #undef UNROLLJ |