Introduce SimulatorBuilder
[gromacs.git] / src / gromacs / nbnxm / kernels_simd_4xm / kernel_outer.h
blob0e7f68b53dedfb1f783860504f3188931ea99087
1 /*
2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
38 using namespace gmx;
40 /* Unpack pointers for output */
41 real *f = out->f.data();
42 real *fshift = out->fshift.data();
43 #ifdef CALC_ENERGIES
44 #ifdef ENERGY_GROUPS
45 real *Vvdw = out->VSvdw.data();
46 real *Vc = out->VSc.data();
47 #else
48 real *Vvdw = out->Vvdw.data();
49 real *Vc = out->Vc.data();
50 #endif
51 #endif
53 const nbnxn_cj_t *l_cj;
54 int ci, ci_sh;
55 int ish, ish3;
56 gmx_bool do_LJ, half_LJ, do_coul;
57 int cjind0, cjind1, cjind;
59 #ifdef ENERGY_GROUPS
60 int Vstride_i;
61 int egps_ishift, egps_imask;
62 int egps_jshift, egps_jmask, egps_jstride;
63 int egps_i;
64 real *vvdwtp[UNROLLI];
65 real *vctp[UNROLLI];
66 #endif
68 SimdReal shX_S;
69 SimdReal shY_S;
70 SimdReal shZ_S;
71 SimdReal ix_S0, iy_S0, iz_S0;
72 SimdReal ix_S1, iy_S1, iz_S1;
73 SimdReal ix_S2, iy_S2, iz_S2;
74 SimdReal ix_S3, iy_S3, iz_S3;
75 SimdReal fix_S0, fiy_S0, fiz_S0;
76 SimdReal fix_S1, fiy_S1, fiz_S1;
77 SimdReal fix_S2, fiy_S2, fiz_S2;
78 SimdReal fix_S3, fiy_S3, fiz_S3;
80 SimdReal diagonal_jmi_S;
81 #if UNROLLI == UNROLLJ
82 SimdBool diagonal_mask_S0, diagonal_mask_S1, diagonal_mask_S2, diagonal_mask_S3;
83 #else
84 SimdBool diagonal_mask0_S0, diagonal_mask0_S1, diagonal_mask0_S2, diagonal_mask0_S3;
85 SimdBool diagonal_mask1_S0, diagonal_mask1_S1, diagonal_mask1_S2, diagonal_mask1_S3;
86 #endif
88 SimdBitMask filter_S0, filter_S1, filter_S2, filter_S3;
90 SimdReal zero_S(0.0);
92 SimdReal one_S(1.0);
93 SimdReal iq_S0 = setZero();
94 SimdReal iq_S1 = setZero();
95 SimdReal iq_S2 = setZero();
96 SimdReal iq_S3 = setZero();
98 #ifdef CALC_COUL_RF
99 SimdReal mrc_3_S;
100 #ifdef CALC_ENERGIES
101 SimdReal hrc_3_S, moh_rc_S;
102 #endif
103 #endif
105 #ifdef CALC_COUL_TAB
106 /* Coulomb table variables */
107 SimdReal invtsp_S;
108 const real * tab_coul_F;
109 #if defined CALC_ENERGIES && !defined TAB_FDV0
110 const real gmx_unused *tab_coul_V;
111 #endif
112 #ifdef CALC_ENERGIES
113 SimdReal mhalfsp_S;
114 #endif
115 #endif
117 #ifdef CALC_COUL_EWALD
118 SimdReal beta2_S, beta_S;
119 #endif
121 #if defined CALC_ENERGIES && (defined CALC_COUL_EWALD || defined CALC_COUL_TAB)
122 SimdReal sh_ewald_S;
123 #endif
125 #if defined LJ_CUT && defined CALC_ENERGIES
126 SimdReal p6_cpot_S, p12_cpot_S;
127 #endif
128 #ifdef LJ_POT_SWITCH
129 SimdReal rswitch_S;
130 SimdReal swV3_S, swV4_S, swV5_S;
131 SimdReal swF2_S, swF3_S, swF4_S;
132 #endif
133 #ifdef LJ_FORCE_SWITCH
134 SimdReal rswitch_S;
135 SimdReal p6_fc2_S, p6_fc3_S;
136 SimdReal p12_fc2_S, p12_fc3_S;
137 #ifdef CALC_ENERGIES
138 SimdReal p6_vc3_S, p6_vc4_S;
139 SimdReal p12_vc3_S, p12_vc4_S;
140 SimdReal p6_6cpot_S, p12_12cpot_S;
141 #endif
142 #endif
143 #ifdef LJ_EWALD_GEOM
144 real lj_ewaldcoeff2, lj_ewaldcoeff6_6;
145 SimdReal mone_S, half_S, lje_c2_S, lje_c6_6_S;
146 #endif
148 #ifdef LJ_COMB_LB
149 SimdReal hsig_i_S0, seps_i_S0;
150 SimdReal hsig_i_S1, seps_i_S1;
151 SimdReal hsig_i_S2, seps_i_S2;
152 SimdReal hsig_i_S3, seps_i_S3;
153 #endif /* LJ_COMB_LB */
155 SimdReal minRsq_S;
156 SimdReal rc2_S;
157 #ifdef VDW_CUTOFF_CHECK
158 SimdReal rcvdw2_S;
159 #endif
161 int ninner;
163 #ifdef COUNT_PAIRS
164 int npair = 0;
165 #endif
167 const nbnxn_atomdata_t::Params &nbatParams = nbat->params();
169 #if defined LJ_COMB_GEOM || defined LJ_COMB_LB || defined LJ_EWALD_GEOM
170 const real * gmx_restrict ljc = nbatParams.lj_comb.data();
171 #endif
172 #if !(defined LJ_COMB_GEOM || defined LJ_COMB_LB || defined FIX_LJ_C)
173 /* No combination rule used */
174 const real * gmx_restrict nbfp_ptr = nbatParams.nbfp_aligned.data();
175 const int * gmx_restrict type = nbatParams.type.data();
176 #endif
178 /* Load j-i for the first i */
179 diagonal_jmi_S = load<SimdReal>(nbat->simdMasks.diagonal_4xn_j_minus_i.data());
180 /* Generate all the diagonal masks as comparison results */
181 #if UNROLLI == UNROLLJ
182 diagonal_mask_S0 = (zero_S < diagonal_jmi_S);
183 diagonal_jmi_S = diagonal_jmi_S - one_S;
184 diagonal_mask_S1 = (zero_S < diagonal_jmi_S);
185 diagonal_jmi_S = diagonal_jmi_S - one_S;
186 diagonal_mask_S2 = (zero_S < diagonal_jmi_S);
187 diagonal_jmi_S = diagonal_jmi_S - one_S;
188 diagonal_mask_S3 = (zero_S < diagonal_jmi_S);
189 #else
190 #if UNROLLI == 2*UNROLLJ || 2*UNROLLI == UNROLLJ
191 diagonal_mask0_S0 = (zero_S < diagonal_jmi_S);
192 diagonal_jmi_S = diagonal_jmi_S - one_S;
193 diagonal_mask0_S1 = (zero_S < diagonal_jmi_S);
194 diagonal_jmi_S = diagonal_jmi_S - one_S;
195 diagonal_mask0_S2 = (zero_S < diagonal_jmi_S);
196 diagonal_jmi_S = diagonal_jmi_S - one_S;
197 diagonal_mask0_S3 = (zero_S < diagonal_jmi_S);
198 diagonal_jmi_S = diagonal_jmi_S - one_S;
200 #if UNROLLI == 2*UNROLLJ
201 /* Load j-i for the second half of the j-cluster */
202 diagonal_jmi_S = load<SimdReal>(nbat->simdMasks.diagonal_4xn_j_minus_i.data() + UNROLLJ);
203 #endif
205 diagonal_mask1_S0 = (zero_S < diagonal_jmi_S);
206 diagonal_jmi_S = diagonal_jmi_S - one_S;
207 diagonal_mask1_S1 = (zero_S < diagonal_jmi_S);
208 diagonal_jmi_S = diagonal_jmi_S - one_S;
209 diagonal_mask1_S2 = (zero_S < diagonal_jmi_S);
210 diagonal_jmi_S = diagonal_jmi_S - one_S;
211 diagonal_mask1_S3 = (zero_S < diagonal_jmi_S);
212 #endif
213 #endif
215 #if GMX_DOUBLE && !GMX_SIMD_HAVE_INT32_LOGICAL
216 const std::uint64_t * gmx_restrict exclusion_filter = nbat->simdMasks.exclusion_filter64.data();
217 #else
218 const std::uint32_t * gmx_restrict exclusion_filter = nbat->simdMasks.exclusion_filter.data();
219 #endif
221 /* Here we cast the exclusion filters from unsigned * to int * or real *.
222 * Since we only check bits, the actual value they represent does not
223 * matter, as long as both filter and mask data are treated the same way.
225 #if GMX_SIMD_HAVE_INT32_LOGICAL
226 filter_S0 = load<SimdBitMask>(reinterpret_cast<const int *>(exclusion_filter + 0*UNROLLJ));
227 filter_S1 = load<SimdBitMask>(reinterpret_cast<const int *>(exclusion_filter + 1*UNROLLJ));
228 filter_S2 = load<SimdBitMask>(reinterpret_cast<const int *>(exclusion_filter + 2*UNROLLJ));
229 filter_S3 = load<SimdBitMask>(reinterpret_cast<const int *>(exclusion_filter + 3*UNROLLJ));
230 #else
231 filter_S0 = load<SimdBitMask>(reinterpret_cast<const real *>(exclusion_filter + 0*UNROLLJ));
232 filter_S1 = load<SimdBitMask>(reinterpret_cast<const real *>(exclusion_filter + 1*UNROLLJ));
233 filter_S2 = load<SimdBitMask>(reinterpret_cast<const real *>(exclusion_filter + 2*UNROLLJ));
234 filter_S3 = load<SimdBitMask>(reinterpret_cast<const real *>(exclusion_filter + 3*UNROLLJ));
235 #endif
237 #ifdef CALC_COUL_RF
238 /* Reaction-field constants */
239 mrc_3_S = SimdReal(-2*ic->k_rf);
240 #ifdef CALC_ENERGIES
241 hrc_3_S = SimdReal(ic->k_rf);
242 moh_rc_S = SimdReal(-ic->c_rf);
243 #endif
244 #endif
246 #ifdef CALC_COUL_TAB
248 invtsp_S = SimdReal(ic->coulombEwaldTables->scale);
249 #ifdef CALC_ENERGIES
250 mhalfsp_S = SimdReal(-0.5/ic->coulombEwaldTables->scale);
251 #endif
253 #ifdef TAB_FDV0
254 tab_coul_F = ic->coulombEwaldTables->tableFDV0.data();
255 #else
256 tab_coul_F = ic->coulombEwaldTables->tableF.data();
257 #ifdef CALC_ENERGIES
258 tab_coul_V = ic->coulombEwaldTables->tableV.data();
259 #endif
260 #endif
261 #endif /* CALC_COUL_TAB */
263 #ifdef CALC_COUL_EWALD
264 beta2_S = SimdReal(ic->ewaldcoeff_q*ic->ewaldcoeff_q);
265 beta_S = SimdReal(ic->ewaldcoeff_q);
266 #endif
268 #if (defined CALC_COUL_TAB || defined CALC_COUL_EWALD) && defined CALC_ENERGIES
269 sh_ewald_S = SimdReal(ic->sh_ewald);
270 #endif
272 /* LJ function constants */
273 #if defined CALC_ENERGIES || defined LJ_POT_SWITCH
274 SimdReal sixth_S(1.0/6.0);
275 SimdReal twelveth_S(1.0/12.0);
276 #endif
278 #if defined LJ_CUT && defined CALC_ENERGIES
279 /* We shift the potential by cpot, which can be zero */
280 p6_cpot_S = SimdReal(ic->dispersion_shift.cpot);
281 p12_cpot_S = SimdReal(ic->repulsion_shift.cpot);
282 #endif
283 #ifdef LJ_POT_SWITCH
284 rswitch_S = SimdReal(ic->rvdw_switch);
285 swV3_S = SimdReal(ic->vdw_switch.c3);
286 swV4_S = SimdReal(ic->vdw_switch.c4);
287 swV5_S = SimdReal(ic->vdw_switch.c5);
288 swF2_S = SimdReal(3*ic->vdw_switch.c3);
289 swF3_S = SimdReal(4*ic->vdw_switch.c4);
290 swF4_S = SimdReal(5*ic->vdw_switch.c5);
291 #endif
292 #ifdef LJ_FORCE_SWITCH
293 rswitch_S = SimdReal(ic->rvdw_switch);
294 p6_fc2_S = SimdReal(ic->dispersion_shift.c2);
295 p6_fc3_S = SimdReal(ic->dispersion_shift.c3);
296 p12_fc2_S = SimdReal(ic->repulsion_shift.c2);
297 p12_fc3_S = SimdReal(ic->repulsion_shift.c3);
298 #ifdef CALC_ENERGIES
300 SimdReal mthird_S(-1.0/3.0);
301 SimdReal mfourth_S(-1.0/4.0);
303 p6_vc3_S = mthird_S * p6_fc2_S;
304 p6_vc4_S = mfourth_S * p6_fc3_S;
305 p6_6cpot_S = SimdReal(ic->dispersion_shift.cpot/6);
306 p12_vc3_S = mthird_S * p12_fc2_S;
307 p12_vc4_S = mfourth_S * p12_fc3_S;
308 p12_12cpot_S = SimdReal(ic->repulsion_shift.cpot/12);
310 #endif
311 #endif
312 #ifdef LJ_EWALD_GEOM
313 mone_S = SimdReal(-1.0);
314 half_S = SimdReal(0.5);
315 lj_ewaldcoeff2 = ic->ewaldcoeff_lj*ic->ewaldcoeff_lj;
316 lj_ewaldcoeff6_6 = lj_ewaldcoeff2*lj_ewaldcoeff2*lj_ewaldcoeff2/6;
317 lje_c2_S = SimdReal(lj_ewaldcoeff2);
318 lje_c6_6_S = SimdReal(lj_ewaldcoeff6_6);
319 #ifdef CALC_ENERGIES
320 /* Determine the grid potential at the cut-off */
321 SimdReal lje_vc_S(ic->sh_lj_ewald);
322 #endif
323 #endif
325 /* The kernel either supports rcoulomb = rvdw or rcoulomb >= rvdw */
326 rc2_S = SimdReal(ic->rcoulomb*ic->rcoulomb);
327 #ifdef VDW_CUTOFF_CHECK
328 rcvdw2_S = SimdReal(ic->rvdw*ic->rvdw);
329 #endif
331 minRsq_S = SimdReal(NBNXN_MIN_RSQ);
333 const real * gmx_restrict q = nbatParams.q.data();
334 const real facel = ic->epsfac;
335 const real * gmx_restrict shiftvec = shift_vec[0];
336 const real * gmx_restrict x = nbat->x().data();
338 #ifdef FIX_LJ_C
339 alignas(GMX_SIMD_ALIGNMENT) real pvdw_c6[2*UNROLLI*UNROLLJ];
340 real *pvdw_c12 = pvdw_c6 + UNROLLI*UNROLLJ;
342 for (int jp = 0; jp < UNROLLJ; jp++)
344 pvdw_c6 [0*UNROLLJ+jp] = nbat->nbfp[0*2];
345 pvdw_c6 [1*UNROLLJ+jp] = nbat->nbfp[0*2];
346 pvdw_c6 [2*UNROLLJ+jp] = nbat->nbfp[0*2];
347 pvdw_c6 [3*UNROLLJ+jp] = nbat->nbfp[0*2];
349 pvdw_c12[0*UNROLLJ+jp] = nbat->nbfp[0*2+1];
350 pvdw_c12[1*UNROLLJ+jp] = nbat->nbfp[0*2+1];
351 pvdw_c12[2*UNROLLJ+jp] = nbat->nbfp[0*2+1];
352 pvdw_c12[3*UNROLLJ+jp] = nbat->nbfp[0*2+1];
354 SimdReal c6_S0 = load<SimdReal>(pvdw_c6 +0*UNROLLJ);
355 SimdReal c6_S1 = load<SimdReal>(pvdw_c6 +1*UNROLLJ);
356 SimdReal c6_S2 = load<SimdReal>(pvdw_c6 +2*UNROLLJ);
357 SimdReal c6_S3 = load<SimdReal>(pvdw_c6 +3*UNROLLJ);
359 SimdReal c12_S0 = load<SimdReal>(pvdw_c12+0*UNROLLJ);
360 SimdReal c12_S1 = load<SimdReal>(pvdw_c12+1*UNROLLJ);
361 SimdReal c12_S2 = load<SimdReal>(pvdw_c12+2*UNROLLJ);
362 SimdReal c12_S3 = load<SimdReal>(pvdw_c12+3*UNROLLJ);
363 #endif /* FIX_LJ_C */
365 #ifdef ENERGY_GROUPS
366 egps_ishift = nbatParams.neg_2log;
367 egps_imask = (1<<egps_ishift) - 1;
368 egps_jshift = 2*nbatParams.neg_2log;
369 egps_jmask = (1<<egps_jshift) - 1;
370 egps_jstride = (UNROLLJ>>1)*UNROLLJ;
371 /* Major division is over i-particle energy groups, determine the stride */
372 Vstride_i = nbatParams.nenergrp*(1 << nbatParams.neg_2log)*egps_jstride;
373 #endif
375 l_cj = nbl->cj.data();
377 ninner = 0;
379 for (const nbnxn_ci_t &ciEntry : nbl->ci)
381 ish = (ciEntry.shift & NBNXN_CI_SHIFT);
382 ish3 = ish*3;
383 cjind0 = ciEntry.cj_ind_start;
384 cjind1 = ciEntry.cj_ind_end;
385 ci = ciEntry.ci;
386 ci_sh = (ish == CENTRAL ? ci : -1);
388 shX_S = SimdReal(shiftvec[ish3]);
389 shY_S = SimdReal(shiftvec[ish3+1]);
390 shZ_S = SimdReal(shiftvec[ish3+2]);
392 #if UNROLLJ <= 4
393 int sci = ci*STRIDE;
394 int scix = sci*DIM;
395 #if defined LJ_COMB_LB || defined LJ_COMB_GEOM || defined LJ_EWALD_GEOM
396 int sci2 = sci*2;
397 #endif
398 #else
399 int sci = (ci>>1)*STRIDE;
400 int scix = sci*DIM + (ci & 1)*(STRIDE>>1);
401 #if defined LJ_COMB_LB || defined LJ_COMB_GEOM || defined LJ_EWALD_GEOM
402 int sci2 = sci*2 + (ci & 1)*(STRIDE>>1);
403 #endif
404 sci += (ci & 1)*(STRIDE>>1);
405 #endif
407 /* We have 5 LJ/C combinations, but use only three inner loops,
408 * as the other combinations are unlikely and/or not much faster:
409 * inner half-LJ + C for half-LJ + C / no-LJ + C
410 * inner LJ + C for full-LJ + C
411 * inner LJ for full-LJ + no-C / half-LJ + no-C
413 do_LJ = ((ciEntry.shift & NBNXN_CI_DO_LJ(0)) != 0);
414 do_coul = ((ciEntry.shift & NBNXN_CI_DO_COUL(0)) != 0);
415 half_LJ = (((ciEntry.shift & NBNXN_CI_HALF_LJ(0)) != 0) || !do_LJ) && do_coul;
417 #ifdef ENERGY_GROUPS
418 egps_i = nbatParams.energrp[ci];
420 int ia, egp_ia;
422 for (ia = 0; ia < UNROLLI; ia++)
424 egp_ia = (egps_i >> (ia*egps_ishift)) & egps_imask;
425 vvdwtp[ia] = Vvdw + egp_ia*Vstride_i;
426 vctp[ia] = Vc + egp_ia*Vstride_i;
429 #endif
431 #ifdef CALC_ENERGIES
432 #ifdef LJ_EWALD_GEOM
433 gmx_bool do_self = TRUE;
434 #else
435 gmx_bool do_self = do_coul;
436 #endif
437 #if UNROLLJ == 4
438 if (do_self && l_cj[ciEntry.cj_ind_start].cj == ci_sh)
439 #endif
440 #if UNROLLJ == 2
441 if (do_self && l_cj[ciEntry.cj_ind_start].cj == (ci_sh<<1))
442 #endif
443 #if UNROLLJ == 8
444 if (do_self && l_cj[ciEntry.cj_ind_start].cj == (ci_sh>>1))
445 #endif
447 if (do_coul)
449 real Vc_sub_self;
450 int ia;
452 #ifdef CALC_COUL_RF
453 Vc_sub_self = 0.5*ic->c_rf;
454 #endif
455 #ifdef CALC_COUL_TAB
456 #ifdef TAB_FDV0
457 Vc_sub_self = 0.5*tab_coul_F[2];
458 #else
459 Vc_sub_self = 0.5*tab_coul_V[0];
460 #endif
461 #endif
462 #ifdef CALC_COUL_EWALD
463 /* beta/sqrt(pi) */
464 Vc_sub_self = 0.5*ic->ewaldcoeff_q*M_2_SQRTPI;
465 #endif
467 for (ia = 0; ia < UNROLLI; ia++)
469 real qi;
471 qi = q[sci+ia];
472 #ifdef ENERGY_GROUPS
473 vctp[ia][((egps_i>>(ia*egps_ishift)) & egps_imask)*egps_jstride]
474 #else
475 Vc[0]
476 #endif
477 -= facel*qi*qi*Vc_sub_self;
481 #ifdef LJ_EWALD_GEOM
483 int ia;
485 for (ia = 0; ia < UNROLLI; ia++)
487 real c6_i;
489 c6_i = nbatParams.nbfp[nbatParams.type[sci+ia]*(nbatParams.numTypes + 1)*2]/6;
490 #ifdef ENERGY_GROUPS
491 vvdwtp[ia][((egps_i>>(ia*egps_ishift)) & egps_imask)*egps_jstride]
492 #else
493 Vvdw[0]
494 #endif
495 += 0.5*c6_i*lj_ewaldcoeff6_6;
498 #endif /* LJ_EWALD_GEOM */
500 #endif
502 /* Load i atom data */
503 int sciy = scix + STRIDE;
504 int sciz = sciy + STRIDE;
505 ix_S0 = SimdReal(x[scix]) + shX_S;
506 ix_S1 = SimdReal(x[scix+1]) + shX_S;
507 ix_S2 = SimdReal(x[scix+2]) + shX_S;
508 ix_S3 = SimdReal(x[scix+3]) + shX_S;
509 iy_S0 = SimdReal(x[sciy]) + shY_S;
510 iy_S1 = SimdReal(x[sciy+1]) + shY_S;
511 iy_S2 = SimdReal(x[sciy+2]) + shY_S;
512 iy_S3 = SimdReal(x[sciy+3]) + shY_S;
513 iz_S0 = SimdReal(x[sciz]) + shZ_S;
514 iz_S1 = SimdReal(x[sciz+1]) + shZ_S;
515 iz_S2 = SimdReal(x[sciz+2]) + shZ_S;
516 iz_S3 = SimdReal(x[sciz+3]) + shZ_S;
518 if (do_coul)
520 iq_S0 = SimdReal(facel*q[sci]);
521 iq_S1 = SimdReal(facel*q[sci+1]);
522 iq_S2 = SimdReal(facel*q[sci+2]);
523 iq_S3 = SimdReal(facel*q[sci+3]);
526 #ifdef LJ_COMB_LB
527 hsig_i_S0 = SimdReal(ljc[sci2+0]);
528 hsig_i_S1 = SimdReal(ljc[sci2+1]);
529 hsig_i_S2 = SimdReal(ljc[sci2+2]);
530 hsig_i_S3 = SimdReal(ljc[sci2+3]);
531 seps_i_S0 = SimdReal(ljc[sci2+STRIDE+0]);
532 seps_i_S1 = SimdReal(ljc[sci2+STRIDE+1]);
533 seps_i_S2 = SimdReal(ljc[sci2+STRIDE+2]);
534 seps_i_S3 = SimdReal(ljc[sci2+STRIDE+3]);
535 #else
536 #ifdef LJ_COMB_GEOM
537 SimdReal c6s_S0 = SimdReal(ljc[sci2+0]);
538 SimdReal c6s_S1 = SimdReal(ljc[sci2+1]);
539 SimdReal c6s_S2, c6s_S3;
540 if (!half_LJ)
542 c6s_S2 = SimdReal(ljc[sci2+2]);
543 c6s_S3 = SimdReal(ljc[sci2+3]);
545 else
547 c6s_S2 = setZero();
548 c6s_S3 = setZero();
550 SimdReal c12s_S0 = SimdReal(ljc[sci2+STRIDE+0]);
551 SimdReal c12s_S1 = SimdReal(ljc[sci2+STRIDE+1]);
552 SimdReal c12s_S2, c12s_S3;
553 if (!half_LJ)
555 c12s_S2 = SimdReal(ljc[sci2+STRIDE+2]);
556 c12s_S3 = SimdReal(ljc[sci2+STRIDE+3]);
558 else
560 c12s_S2 = setZero();
561 c12s_S3 = setZero();
563 #else
564 const int numTypes = nbatParams.numTypes;
565 const real *nbfp0 = nbfp_ptr + type[sci ]*numTypes*c_simdBestPairAlignment;
566 const real *nbfp1 = nbfp_ptr + type[sci+1]*numTypes*c_simdBestPairAlignment;
567 const real *nbfp2 = nullptr, *nbfp3 = nullptr;
568 if (!half_LJ)
570 nbfp2 = nbfp_ptr + type[sci+2]*numTypes*c_simdBestPairAlignment;
571 nbfp3 = nbfp_ptr + type[sci+3]*numTypes*c_simdBestPairAlignment;
573 #endif
574 #endif
575 #ifdef LJ_EWALD_GEOM
576 /* We need the geometrically combined C6 for the PME grid correction */
577 SimdReal c6s_S0 = SimdReal(ljc[sci2+0]);
578 SimdReal c6s_S1 = SimdReal(ljc[sci2+1]);
579 SimdReal c6s_S2, c6s_S3;
580 if (!half_LJ)
582 c6s_S2 = SimdReal(ljc[sci2+2]);
583 c6s_S3 = SimdReal(ljc[sci2+3]);
585 #endif
587 /* Zero the potential energy for this list */
588 #ifdef CALC_ENERGIES
589 SimdReal Vvdwtot_S = setZero();
590 SimdReal vctot_S = setZero();
591 #endif
593 /* Clear i atom forces */
594 fix_S0 = setZero();
595 fix_S1 = setZero();
596 fix_S2 = setZero();
597 fix_S3 = setZero();
598 fiy_S0 = setZero();
599 fiy_S1 = setZero();
600 fiy_S2 = setZero();
601 fiy_S3 = setZero();
602 fiz_S0 = setZero();
603 fiz_S1 = setZero();
604 fiz_S2 = setZero();
605 fiz_S3 = setZero();
607 cjind = cjind0;
609 /* Currently all kernels use (at least half) LJ */
610 #define CALC_LJ
611 if (half_LJ)
613 /* Coulomb: all i-atoms, LJ: first half i-atoms */
614 #define CALC_COULOMB
615 #define HALF_LJ
616 #define CHECK_EXCLS
617 while (cjind < cjind1 && nbl->cj[cjind].excl != NBNXN_INTERACTION_MASK_ALL)
619 #include "gromacs/nbnxm/kernels_simd_4xm/kernel_inner.h"
620 cjind++;
622 #undef CHECK_EXCLS
623 for (; (cjind < cjind1); cjind++)
625 #include "gromacs/nbnxm/kernels_simd_4xm/kernel_inner.h"
627 #undef HALF_LJ
628 #undef CALC_COULOMB
630 else if (do_coul)
632 /* Coulomb: all i-atoms, LJ: all i-atoms */
633 #define CALC_COULOMB
634 #define CHECK_EXCLS
635 while (cjind < cjind1 && nbl->cj[cjind].excl != NBNXN_INTERACTION_MASK_ALL)
637 #include "gromacs/nbnxm/kernels_simd_4xm/kernel_inner.h"
638 cjind++;
640 #undef CHECK_EXCLS
641 for (; (cjind < cjind1); cjind++)
643 #include "gromacs/nbnxm/kernels_simd_4xm/kernel_inner.h"
645 #undef CALC_COULOMB
647 else
649 /* Coulomb: none, LJ: all i-atoms */
650 #define CHECK_EXCLS
651 while (cjind < cjind1 && nbl->cj[cjind].excl != NBNXN_INTERACTION_MASK_ALL)
653 #include "gromacs/nbnxm/kernels_simd_4xm/kernel_inner.h"
654 cjind++;
656 #undef CHECK_EXCLS
657 for (; (cjind < cjind1); cjind++)
659 #include "gromacs/nbnxm/kernels_simd_4xm/kernel_inner.h"
662 #undef CALC_LJ
663 ninner += cjind1 - cjind0;
665 /* Add accumulated i-forces to the force array */
666 real fShiftX = reduceIncr4ReturnSum(f+scix, fix_S0, fix_S1, fix_S2, fix_S3);
667 real fShiftY = reduceIncr4ReturnSum(f+sciy, fiy_S0, fiy_S1, fiy_S2, fiy_S3);
668 real fShiftZ = reduceIncr4ReturnSum(f+sciz, fiz_S0, fiz_S1, fiz_S2, fiz_S3);
671 #ifdef CALC_SHIFTFORCES
672 fshift[ish3+0] += fShiftX;
673 fshift[ish3+1] += fShiftY;
674 fshift[ish3+2] += fShiftZ;
675 #endif
677 #ifdef CALC_ENERGIES
678 if (do_coul)
680 *Vc += reduce(vctot_S);
683 *Vvdw += reduce(Vvdwtot_S);
684 #endif
686 /* Outer loop uses 6 flops/iteration */
689 #ifdef COUNT_PAIRS
690 printf("atom pairs %d\n", npair);
691 #endif