2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014, 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.
36 * Note: this file was generated by the GROMACS avx_128_fma_double kernel generator.
42 #include "../nb_kernel.h"
43 #include "types/simple.h"
44 #include "gromacs/math/vec.h"
47 #include "gromacs/simd/math_x86_avx_128_fma_double.h"
48 #include "kernelutil_x86_avx_128_fma_double.h"
51 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwCSTab_GeomW4P1_VF_avx_128_fma_double
52 * Electrostatics interaction: ReactionField
53 * VdW interaction: CubicSplineTable
54 * Geometry: Water4-Particle
55 * Calculate force/pot: PotentialAndForce
58 nb_kernel_ElecRFCut_VdwCSTab_GeomW4P1_VF_avx_128_fma_double
59 (t_nblist * gmx_restrict nlist,
60 rvec * gmx_restrict xx,
61 rvec * gmx_restrict ff,
62 t_forcerec * gmx_restrict fr,
63 t_mdatoms * gmx_restrict mdatoms,
64 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
65 t_nrnb * gmx_restrict nrnb)
67 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
68 * just 0 for non-waters.
69 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
70 * jnr indices corresponding to data put in the four positions in the SIMD register.
72 int i_shift_offset,i_coord_offset,outeriter,inneriter;
73 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
75 int j_coord_offsetA,j_coord_offsetB;
76 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
78 real *shiftvec,*fshift,*x,*f;
79 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
81 __m128d ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
83 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
85 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
87 __m128d ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
88 int vdwjidx0A,vdwjidx0B;
89 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
90 __m128d dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
91 __m128d dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
92 __m128d dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
93 __m128d dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
94 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
97 __m128d rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
100 __m128d one_sixth = _mm_set1_pd(1.0/6.0);
101 __m128d one_twelfth = _mm_set1_pd(1.0/12.0);
103 __m128i ifour = _mm_set1_epi32(4);
104 __m128d rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF,twovfeps;
106 __m128d dummy_mask,cutoff_mask;
107 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
108 __m128d one = _mm_set1_pd(1.0);
109 __m128d two = _mm_set1_pd(2.0);
115 jindex = nlist->jindex;
117 shiftidx = nlist->shift;
119 shiftvec = fr->shift_vec[0];
120 fshift = fr->fshift[0];
121 facel = _mm_set1_pd(fr->epsfac);
122 charge = mdatoms->chargeA;
123 krf = _mm_set1_pd(fr->ic->k_rf);
124 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
125 crf = _mm_set1_pd(fr->ic->c_rf);
126 nvdwtype = fr->ntype;
128 vdwtype = mdatoms->typeA;
130 vftab = kernel_data->table_vdw->data;
131 vftabscale = _mm_set1_pd(kernel_data->table_vdw->scale);
133 /* Setup water-specific parameters */
134 inr = nlist->iinr[0];
135 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
136 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
137 iq3 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+3]));
138 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
140 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
141 rcutoff_scalar = fr->rcoulomb;
142 rcutoff = _mm_set1_pd(rcutoff_scalar);
143 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
145 /* Avoid stupid compiler warnings */
153 /* Start outer loop over neighborlists */
154 for(iidx=0; iidx<nri; iidx++)
156 /* Load shift vector for this list */
157 i_shift_offset = DIM*shiftidx[iidx];
159 /* Load limits for loop over neighbors */
160 j_index_start = jindex[iidx];
161 j_index_end = jindex[iidx+1];
163 /* Get outer coordinate index */
165 i_coord_offset = DIM*inr;
167 /* Load i particle coords and add shift vector */
168 gmx_mm_load_shift_and_4rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
169 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
171 fix0 = _mm_setzero_pd();
172 fiy0 = _mm_setzero_pd();
173 fiz0 = _mm_setzero_pd();
174 fix1 = _mm_setzero_pd();
175 fiy1 = _mm_setzero_pd();
176 fiz1 = _mm_setzero_pd();
177 fix2 = _mm_setzero_pd();
178 fiy2 = _mm_setzero_pd();
179 fiz2 = _mm_setzero_pd();
180 fix3 = _mm_setzero_pd();
181 fiy3 = _mm_setzero_pd();
182 fiz3 = _mm_setzero_pd();
184 /* Reset potential sums */
185 velecsum = _mm_setzero_pd();
186 vvdwsum = _mm_setzero_pd();
188 /* Start inner kernel loop */
189 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
192 /* Get j neighbor index, and coordinate index */
195 j_coord_offsetA = DIM*jnrA;
196 j_coord_offsetB = DIM*jnrB;
198 /* load j atom coordinates */
199 gmx_mm_load_1rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
202 /* Calculate displacement vector */
203 dx00 = _mm_sub_pd(ix0,jx0);
204 dy00 = _mm_sub_pd(iy0,jy0);
205 dz00 = _mm_sub_pd(iz0,jz0);
206 dx10 = _mm_sub_pd(ix1,jx0);
207 dy10 = _mm_sub_pd(iy1,jy0);
208 dz10 = _mm_sub_pd(iz1,jz0);
209 dx20 = _mm_sub_pd(ix2,jx0);
210 dy20 = _mm_sub_pd(iy2,jy0);
211 dz20 = _mm_sub_pd(iz2,jz0);
212 dx30 = _mm_sub_pd(ix3,jx0);
213 dy30 = _mm_sub_pd(iy3,jy0);
214 dz30 = _mm_sub_pd(iz3,jz0);
216 /* Calculate squared distance and things based on it */
217 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
218 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
219 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
220 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
222 rinv00 = gmx_mm_invsqrt_pd(rsq00);
223 rinv10 = gmx_mm_invsqrt_pd(rsq10);
224 rinv20 = gmx_mm_invsqrt_pd(rsq20);
225 rinv30 = gmx_mm_invsqrt_pd(rsq30);
227 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
228 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
229 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
231 /* Load parameters for j particles */
232 jq0 = gmx_mm_load_2real_swizzle_pd(charge+jnrA+0,charge+jnrB+0);
233 vdwjidx0A = 2*vdwtype[jnrA+0];
234 vdwjidx0B = 2*vdwtype[jnrB+0];
236 fjx0 = _mm_setzero_pd();
237 fjy0 = _mm_setzero_pd();
238 fjz0 = _mm_setzero_pd();
240 /**************************
241 * CALCULATE INTERACTIONS *
242 **************************/
244 r00 = _mm_mul_pd(rsq00,rinv00);
246 /* Compute parameters for interactions between i and j atoms */
247 gmx_mm_load_2pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,
248 vdwparam+vdwioffset0+vdwjidx0B,&c6_00,&c12_00);
250 /* Calculate table index by multiplying r with table scale and truncate to integer */
251 rt = _mm_mul_pd(r00,vftabscale);
252 vfitab = _mm_cvttpd_epi32(rt);
254 vfeps = _mm_frcz_pd(rt);
256 vfeps = _mm_sub_pd(rt,_mm_round_pd(rt, _MM_FROUND_FLOOR));
258 twovfeps = _mm_add_pd(vfeps,vfeps);
259 vfitab = _mm_slli_epi32(vfitab,3);
261 /* CUBIC SPLINE TABLE DISPERSION */
262 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
263 F = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
264 GMX_MM_TRANSPOSE2_PD(Y,F);
265 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
266 H = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) +2);
267 GMX_MM_TRANSPOSE2_PD(G,H);
268 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
269 VV = _mm_macc_pd(vfeps,Fp,Y);
270 vvdw6 = _mm_mul_pd(c6_00,VV);
271 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
272 fvdw6 = _mm_mul_pd(c6_00,FF);
274 /* CUBIC SPLINE TABLE REPULSION */
275 vfitab = _mm_add_epi32(vfitab,ifour);
276 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
277 F = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
278 GMX_MM_TRANSPOSE2_PD(Y,F);
279 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
280 H = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) +2);
281 GMX_MM_TRANSPOSE2_PD(G,H);
282 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
283 VV = _mm_macc_pd(vfeps,Fp,Y);
284 vvdw12 = _mm_mul_pd(c12_00,VV);
285 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
286 fvdw12 = _mm_mul_pd(c12_00,FF);
287 vvdw = _mm_add_pd(vvdw12,vvdw6);
288 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
290 /* Update potential sum for this i atom from the interaction with this j atom. */
291 vvdwsum = _mm_add_pd(vvdwsum,vvdw);
295 /* Update vectorial force */
296 fix0 = _mm_macc_pd(dx00,fscal,fix0);
297 fiy0 = _mm_macc_pd(dy00,fscal,fiy0);
298 fiz0 = _mm_macc_pd(dz00,fscal,fiz0);
300 fjx0 = _mm_macc_pd(dx00,fscal,fjx0);
301 fjy0 = _mm_macc_pd(dy00,fscal,fjy0);
302 fjz0 = _mm_macc_pd(dz00,fscal,fjz0);
304 /**************************
305 * CALCULATE INTERACTIONS *
306 **************************/
308 if (gmx_mm_any_lt(rsq10,rcutoff2))
311 /* Compute parameters for interactions between i and j atoms */
312 qq10 = _mm_mul_pd(iq1,jq0);
314 /* REACTION-FIELD ELECTROSTATICS */
315 velec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_macc_pd(krf,rsq10,rinv10),crf));
316 felec = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
318 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
320 /* Update potential sum for this i atom from the interaction with this j atom. */
321 velec = _mm_and_pd(velec,cutoff_mask);
322 velecsum = _mm_add_pd(velecsum,velec);
326 fscal = _mm_and_pd(fscal,cutoff_mask);
328 /* Update vectorial force */
329 fix1 = _mm_macc_pd(dx10,fscal,fix1);
330 fiy1 = _mm_macc_pd(dy10,fscal,fiy1);
331 fiz1 = _mm_macc_pd(dz10,fscal,fiz1);
333 fjx0 = _mm_macc_pd(dx10,fscal,fjx0);
334 fjy0 = _mm_macc_pd(dy10,fscal,fjy0);
335 fjz0 = _mm_macc_pd(dz10,fscal,fjz0);
339 /**************************
340 * CALCULATE INTERACTIONS *
341 **************************/
343 if (gmx_mm_any_lt(rsq20,rcutoff2))
346 /* Compute parameters for interactions between i and j atoms */
347 qq20 = _mm_mul_pd(iq2,jq0);
349 /* REACTION-FIELD ELECTROSTATICS */
350 velec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_macc_pd(krf,rsq20,rinv20),crf));
351 felec = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
353 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
355 /* Update potential sum for this i atom from the interaction with this j atom. */
356 velec = _mm_and_pd(velec,cutoff_mask);
357 velecsum = _mm_add_pd(velecsum,velec);
361 fscal = _mm_and_pd(fscal,cutoff_mask);
363 /* Update vectorial force */
364 fix2 = _mm_macc_pd(dx20,fscal,fix2);
365 fiy2 = _mm_macc_pd(dy20,fscal,fiy2);
366 fiz2 = _mm_macc_pd(dz20,fscal,fiz2);
368 fjx0 = _mm_macc_pd(dx20,fscal,fjx0);
369 fjy0 = _mm_macc_pd(dy20,fscal,fjy0);
370 fjz0 = _mm_macc_pd(dz20,fscal,fjz0);
374 /**************************
375 * CALCULATE INTERACTIONS *
376 **************************/
378 if (gmx_mm_any_lt(rsq30,rcutoff2))
381 /* Compute parameters for interactions between i and j atoms */
382 qq30 = _mm_mul_pd(iq3,jq0);
384 /* REACTION-FIELD ELECTROSTATICS */
385 velec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_macc_pd(krf,rsq30,rinv30),crf));
386 felec = _mm_mul_pd(qq30,_mm_msub_pd(rinv30,rinvsq30,krf2));
388 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
390 /* Update potential sum for this i atom from the interaction with this j atom. */
391 velec = _mm_and_pd(velec,cutoff_mask);
392 velecsum = _mm_add_pd(velecsum,velec);
396 fscal = _mm_and_pd(fscal,cutoff_mask);
398 /* Update vectorial force */
399 fix3 = _mm_macc_pd(dx30,fscal,fix3);
400 fiy3 = _mm_macc_pd(dy30,fscal,fiy3);
401 fiz3 = _mm_macc_pd(dz30,fscal,fiz3);
403 fjx0 = _mm_macc_pd(dx30,fscal,fjx0);
404 fjy0 = _mm_macc_pd(dy30,fscal,fjy0);
405 fjz0 = _mm_macc_pd(dz30,fscal,fjz0);
409 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0);
411 /* Inner loop uses 179 flops */
418 j_coord_offsetA = DIM*jnrA;
420 /* load j atom coordinates */
421 gmx_mm_load_1rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
424 /* Calculate displacement vector */
425 dx00 = _mm_sub_pd(ix0,jx0);
426 dy00 = _mm_sub_pd(iy0,jy0);
427 dz00 = _mm_sub_pd(iz0,jz0);
428 dx10 = _mm_sub_pd(ix1,jx0);
429 dy10 = _mm_sub_pd(iy1,jy0);
430 dz10 = _mm_sub_pd(iz1,jz0);
431 dx20 = _mm_sub_pd(ix2,jx0);
432 dy20 = _mm_sub_pd(iy2,jy0);
433 dz20 = _mm_sub_pd(iz2,jz0);
434 dx30 = _mm_sub_pd(ix3,jx0);
435 dy30 = _mm_sub_pd(iy3,jy0);
436 dz30 = _mm_sub_pd(iz3,jz0);
438 /* Calculate squared distance and things based on it */
439 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
440 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
441 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
442 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
444 rinv00 = gmx_mm_invsqrt_pd(rsq00);
445 rinv10 = gmx_mm_invsqrt_pd(rsq10);
446 rinv20 = gmx_mm_invsqrt_pd(rsq20);
447 rinv30 = gmx_mm_invsqrt_pd(rsq30);
449 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
450 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
451 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
453 /* Load parameters for j particles */
454 jq0 = _mm_load_sd(charge+jnrA+0);
455 vdwjidx0A = 2*vdwtype[jnrA+0];
457 fjx0 = _mm_setzero_pd();
458 fjy0 = _mm_setzero_pd();
459 fjz0 = _mm_setzero_pd();
461 /**************************
462 * CALCULATE INTERACTIONS *
463 **************************/
465 r00 = _mm_mul_pd(rsq00,rinv00);
467 /* Compute parameters for interactions between i and j atoms */
468 gmx_mm_load_1pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,&c6_00,&c12_00);
470 /* Calculate table index by multiplying r with table scale and truncate to integer */
471 rt = _mm_mul_pd(r00,vftabscale);
472 vfitab = _mm_cvttpd_epi32(rt);
474 vfeps = _mm_frcz_pd(rt);
476 vfeps = _mm_sub_pd(rt,_mm_round_pd(rt, _MM_FROUND_FLOOR));
478 twovfeps = _mm_add_pd(vfeps,vfeps);
479 vfitab = _mm_slli_epi32(vfitab,3);
481 /* CUBIC SPLINE TABLE DISPERSION */
482 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
483 F = _mm_setzero_pd();
484 GMX_MM_TRANSPOSE2_PD(Y,F);
485 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
486 H = _mm_setzero_pd();
487 GMX_MM_TRANSPOSE2_PD(G,H);
488 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
489 VV = _mm_macc_pd(vfeps,Fp,Y);
490 vvdw6 = _mm_mul_pd(c6_00,VV);
491 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
492 fvdw6 = _mm_mul_pd(c6_00,FF);
494 /* CUBIC SPLINE TABLE REPULSION */
495 vfitab = _mm_add_epi32(vfitab,ifour);
496 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
497 F = _mm_setzero_pd();
498 GMX_MM_TRANSPOSE2_PD(Y,F);
499 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
500 H = _mm_setzero_pd();
501 GMX_MM_TRANSPOSE2_PD(G,H);
502 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
503 VV = _mm_macc_pd(vfeps,Fp,Y);
504 vvdw12 = _mm_mul_pd(c12_00,VV);
505 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
506 fvdw12 = _mm_mul_pd(c12_00,FF);
507 vvdw = _mm_add_pd(vvdw12,vvdw6);
508 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
510 /* Update potential sum for this i atom from the interaction with this j atom. */
511 vvdw = _mm_unpacklo_pd(vvdw,_mm_setzero_pd());
512 vvdwsum = _mm_add_pd(vvdwsum,vvdw);
516 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
518 /* Update vectorial force */
519 fix0 = _mm_macc_pd(dx00,fscal,fix0);
520 fiy0 = _mm_macc_pd(dy00,fscal,fiy0);
521 fiz0 = _mm_macc_pd(dz00,fscal,fiz0);
523 fjx0 = _mm_macc_pd(dx00,fscal,fjx0);
524 fjy0 = _mm_macc_pd(dy00,fscal,fjy0);
525 fjz0 = _mm_macc_pd(dz00,fscal,fjz0);
527 /**************************
528 * CALCULATE INTERACTIONS *
529 **************************/
531 if (gmx_mm_any_lt(rsq10,rcutoff2))
534 /* Compute parameters for interactions between i and j atoms */
535 qq10 = _mm_mul_pd(iq1,jq0);
537 /* REACTION-FIELD ELECTROSTATICS */
538 velec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_macc_pd(krf,rsq10,rinv10),crf));
539 felec = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
541 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
543 /* Update potential sum for this i atom from the interaction with this j atom. */
544 velec = _mm_and_pd(velec,cutoff_mask);
545 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
546 velecsum = _mm_add_pd(velecsum,velec);
550 fscal = _mm_and_pd(fscal,cutoff_mask);
552 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
554 /* Update vectorial force */
555 fix1 = _mm_macc_pd(dx10,fscal,fix1);
556 fiy1 = _mm_macc_pd(dy10,fscal,fiy1);
557 fiz1 = _mm_macc_pd(dz10,fscal,fiz1);
559 fjx0 = _mm_macc_pd(dx10,fscal,fjx0);
560 fjy0 = _mm_macc_pd(dy10,fscal,fjy0);
561 fjz0 = _mm_macc_pd(dz10,fscal,fjz0);
565 /**************************
566 * CALCULATE INTERACTIONS *
567 **************************/
569 if (gmx_mm_any_lt(rsq20,rcutoff2))
572 /* Compute parameters for interactions between i and j atoms */
573 qq20 = _mm_mul_pd(iq2,jq0);
575 /* REACTION-FIELD ELECTROSTATICS */
576 velec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_macc_pd(krf,rsq20,rinv20),crf));
577 felec = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
579 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
581 /* Update potential sum for this i atom from the interaction with this j atom. */
582 velec = _mm_and_pd(velec,cutoff_mask);
583 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
584 velecsum = _mm_add_pd(velecsum,velec);
588 fscal = _mm_and_pd(fscal,cutoff_mask);
590 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
592 /* Update vectorial force */
593 fix2 = _mm_macc_pd(dx20,fscal,fix2);
594 fiy2 = _mm_macc_pd(dy20,fscal,fiy2);
595 fiz2 = _mm_macc_pd(dz20,fscal,fiz2);
597 fjx0 = _mm_macc_pd(dx20,fscal,fjx0);
598 fjy0 = _mm_macc_pd(dy20,fscal,fjy0);
599 fjz0 = _mm_macc_pd(dz20,fscal,fjz0);
603 /**************************
604 * CALCULATE INTERACTIONS *
605 **************************/
607 if (gmx_mm_any_lt(rsq30,rcutoff2))
610 /* Compute parameters for interactions between i and j atoms */
611 qq30 = _mm_mul_pd(iq3,jq0);
613 /* REACTION-FIELD ELECTROSTATICS */
614 velec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_macc_pd(krf,rsq30,rinv30),crf));
615 felec = _mm_mul_pd(qq30,_mm_msub_pd(rinv30,rinvsq30,krf2));
617 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
619 /* Update potential sum for this i atom from the interaction with this j atom. */
620 velec = _mm_and_pd(velec,cutoff_mask);
621 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
622 velecsum = _mm_add_pd(velecsum,velec);
626 fscal = _mm_and_pd(fscal,cutoff_mask);
628 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
630 /* Update vectorial force */
631 fix3 = _mm_macc_pd(dx30,fscal,fix3);
632 fiy3 = _mm_macc_pd(dy30,fscal,fiy3);
633 fiz3 = _mm_macc_pd(dz30,fscal,fiz3);
635 fjx0 = _mm_macc_pd(dx30,fscal,fjx0);
636 fjy0 = _mm_macc_pd(dy30,fscal,fjy0);
637 fjz0 = _mm_macc_pd(dz30,fscal,fjz0);
641 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0);
643 /* Inner loop uses 179 flops */
646 /* End of innermost loop */
648 gmx_mm_update_iforce_4atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
649 f+i_coord_offset,fshift+i_shift_offset);
652 /* Update potential energies */
653 gmx_mm_update_1pot_pd(velecsum,kernel_data->energygrp_elec+ggid);
654 gmx_mm_update_1pot_pd(vvdwsum,kernel_data->energygrp_vdw+ggid);
656 /* Increment number of inner iterations */
657 inneriter += j_index_end - j_index_start;
659 /* Outer loop uses 26 flops */
662 /* Increment number of outer iterations */
665 /* Update outer/inner flops */
667 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4_VF,outeriter*26 + inneriter*179);
670 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwCSTab_GeomW4P1_F_avx_128_fma_double
671 * Electrostatics interaction: ReactionField
672 * VdW interaction: CubicSplineTable
673 * Geometry: Water4-Particle
674 * Calculate force/pot: Force
677 nb_kernel_ElecRFCut_VdwCSTab_GeomW4P1_F_avx_128_fma_double
678 (t_nblist * gmx_restrict nlist,
679 rvec * gmx_restrict xx,
680 rvec * gmx_restrict ff,
681 t_forcerec * gmx_restrict fr,
682 t_mdatoms * gmx_restrict mdatoms,
683 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
684 t_nrnb * gmx_restrict nrnb)
686 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
687 * just 0 for non-waters.
688 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
689 * jnr indices corresponding to data put in the four positions in the SIMD register.
691 int i_shift_offset,i_coord_offset,outeriter,inneriter;
692 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
694 int j_coord_offsetA,j_coord_offsetB;
695 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
697 real *shiftvec,*fshift,*x,*f;
698 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
700 __m128d ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
702 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
704 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
706 __m128d ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
707 int vdwjidx0A,vdwjidx0B;
708 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
709 __m128d dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
710 __m128d dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
711 __m128d dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
712 __m128d dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
713 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
716 __m128d rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
719 __m128d one_sixth = _mm_set1_pd(1.0/6.0);
720 __m128d one_twelfth = _mm_set1_pd(1.0/12.0);
722 __m128i ifour = _mm_set1_epi32(4);
723 __m128d rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF,twovfeps;
725 __m128d dummy_mask,cutoff_mask;
726 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
727 __m128d one = _mm_set1_pd(1.0);
728 __m128d two = _mm_set1_pd(2.0);
734 jindex = nlist->jindex;
736 shiftidx = nlist->shift;
738 shiftvec = fr->shift_vec[0];
739 fshift = fr->fshift[0];
740 facel = _mm_set1_pd(fr->epsfac);
741 charge = mdatoms->chargeA;
742 krf = _mm_set1_pd(fr->ic->k_rf);
743 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
744 crf = _mm_set1_pd(fr->ic->c_rf);
745 nvdwtype = fr->ntype;
747 vdwtype = mdatoms->typeA;
749 vftab = kernel_data->table_vdw->data;
750 vftabscale = _mm_set1_pd(kernel_data->table_vdw->scale);
752 /* Setup water-specific parameters */
753 inr = nlist->iinr[0];
754 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
755 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
756 iq3 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+3]));
757 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
759 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
760 rcutoff_scalar = fr->rcoulomb;
761 rcutoff = _mm_set1_pd(rcutoff_scalar);
762 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
764 /* Avoid stupid compiler warnings */
772 /* Start outer loop over neighborlists */
773 for(iidx=0; iidx<nri; iidx++)
775 /* Load shift vector for this list */
776 i_shift_offset = DIM*shiftidx[iidx];
778 /* Load limits for loop over neighbors */
779 j_index_start = jindex[iidx];
780 j_index_end = jindex[iidx+1];
782 /* Get outer coordinate index */
784 i_coord_offset = DIM*inr;
786 /* Load i particle coords and add shift vector */
787 gmx_mm_load_shift_and_4rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
788 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
790 fix0 = _mm_setzero_pd();
791 fiy0 = _mm_setzero_pd();
792 fiz0 = _mm_setzero_pd();
793 fix1 = _mm_setzero_pd();
794 fiy1 = _mm_setzero_pd();
795 fiz1 = _mm_setzero_pd();
796 fix2 = _mm_setzero_pd();
797 fiy2 = _mm_setzero_pd();
798 fiz2 = _mm_setzero_pd();
799 fix3 = _mm_setzero_pd();
800 fiy3 = _mm_setzero_pd();
801 fiz3 = _mm_setzero_pd();
803 /* Start inner kernel loop */
804 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
807 /* Get j neighbor index, and coordinate index */
810 j_coord_offsetA = DIM*jnrA;
811 j_coord_offsetB = DIM*jnrB;
813 /* load j atom coordinates */
814 gmx_mm_load_1rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
817 /* Calculate displacement vector */
818 dx00 = _mm_sub_pd(ix0,jx0);
819 dy00 = _mm_sub_pd(iy0,jy0);
820 dz00 = _mm_sub_pd(iz0,jz0);
821 dx10 = _mm_sub_pd(ix1,jx0);
822 dy10 = _mm_sub_pd(iy1,jy0);
823 dz10 = _mm_sub_pd(iz1,jz0);
824 dx20 = _mm_sub_pd(ix2,jx0);
825 dy20 = _mm_sub_pd(iy2,jy0);
826 dz20 = _mm_sub_pd(iz2,jz0);
827 dx30 = _mm_sub_pd(ix3,jx0);
828 dy30 = _mm_sub_pd(iy3,jy0);
829 dz30 = _mm_sub_pd(iz3,jz0);
831 /* Calculate squared distance and things based on it */
832 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
833 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
834 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
835 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
837 rinv00 = gmx_mm_invsqrt_pd(rsq00);
838 rinv10 = gmx_mm_invsqrt_pd(rsq10);
839 rinv20 = gmx_mm_invsqrt_pd(rsq20);
840 rinv30 = gmx_mm_invsqrt_pd(rsq30);
842 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
843 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
844 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
846 /* Load parameters for j particles */
847 jq0 = gmx_mm_load_2real_swizzle_pd(charge+jnrA+0,charge+jnrB+0);
848 vdwjidx0A = 2*vdwtype[jnrA+0];
849 vdwjidx0B = 2*vdwtype[jnrB+0];
851 fjx0 = _mm_setzero_pd();
852 fjy0 = _mm_setzero_pd();
853 fjz0 = _mm_setzero_pd();
855 /**************************
856 * CALCULATE INTERACTIONS *
857 **************************/
859 r00 = _mm_mul_pd(rsq00,rinv00);
861 /* Compute parameters for interactions between i and j atoms */
862 gmx_mm_load_2pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,
863 vdwparam+vdwioffset0+vdwjidx0B,&c6_00,&c12_00);
865 /* Calculate table index by multiplying r with table scale and truncate to integer */
866 rt = _mm_mul_pd(r00,vftabscale);
867 vfitab = _mm_cvttpd_epi32(rt);
869 vfeps = _mm_frcz_pd(rt);
871 vfeps = _mm_sub_pd(rt,_mm_round_pd(rt, _MM_FROUND_FLOOR));
873 twovfeps = _mm_add_pd(vfeps,vfeps);
874 vfitab = _mm_slli_epi32(vfitab,3);
876 /* CUBIC SPLINE TABLE DISPERSION */
877 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
878 F = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
879 GMX_MM_TRANSPOSE2_PD(Y,F);
880 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
881 H = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) +2);
882 GMX_MM_TRANSPOSE2_PD(G,H);
883 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
884 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
885 fvdw6 = _mm_mul_pd(c6_00,FF);
887 /* CUBIC SPLINE TABLE REPULSION */
888 vfitab = _mm_add_epi32(vfitab,ifour);
889 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
890 F = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
891 GMX_MM_TRANSPOSE2_PD(Y,F);
892 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
893 H = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,1) +2);
894 GMX_MM_TRANSPOSE2_PD(G,H);
895 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
896 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
897 fvdw12 = _mm_mul_pd(c12_00,FF);
898 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
902 /* Update vectorial force */
903 fix0 = _mm_macc_pd(dx00,fscal,fix0);
904 fiy0 = _mm_macc_pd(dy00,fscal,fiy0);
905 fiz0 = _mm_macc_pd(dz00,fscal,fiz0);
907 fjx0 = _mm_macc_pd(dx00,fscal,fjx0);
908 fjy0 = _mm_macc_pd(dy00,fscal,fjy0);
909 fjz0 = _mm_macc_pd(dz00,fscal,fjz0);
911 /**************************
912 * CALCULATE INTERACTIONS *
913 **************************/
915 if (gmx_mm_any_lt(rsq10,rcutoff2))
918 /* Compute parameters for interactions between i and j atoms */
919 qq10 = _mm_mul_pd(iq1,jq0);
921 /* REACTION-FIELD ELECTROSTATICS */
922 felec = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
924 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
928 fscal = _mm_and_pd(fscal,cutoff_mask);
930 /* Update vectorial force */
931 fix1 = _mm_macc_pd(dx10,fscal,fix1);
932 fiy1 = _mm_macc_pd(dy10,fscal,fiy1);
933 fiz1 = _mm_macc_pd(dz10,fscal,fiz1);
935 fjx0 = _mm_macc_pd(dx10,fscal,fjx0);
936 fjy0 = _mm_macc_pd(dy10,fscal,fjy0);
937 fjz0 = _mm_macc_pd(dz10,fscal,fjz0);
941 /**************************
942 * CALCULATE INTERACTIONS *
943 **************************/
945 if (gmx_mm_any_lt(rsq20,rcutoff2))
948 /* Compute parameters for interactions between i and j atoms */
949 qq20 = _mm_mul_pd(iq2,jq0);
951 /* REACTION-FIELD ELECTROSTATICS */
952 felec = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
954 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
958 fscal = _mm_and_pd(fscal,cutoff_mask);
960 /* Update vectorial force */
961 fix2 = _mm_macc_pd(dx20,fscal,fix2);
962 fiy2 = _mm_macc_pd(dy20,fscal,fiy2);
963 fiz2 = _mm_macc_pd(dz20,fscal,fiz2);
965 fjx0 = _mm_macc_pd(dx20,fscal,fjx0);
966 fjy0 = _mm_macc_pd(dy20,fscal,fjy0);
967 fjz0 = _mm_macc_pd(dz20,fscal,fjz0);
971 /**************************
972 * CALCULATE INTERACTIONS *
973 **************************/
975 if (gmx_mm_any_lt(rsq30,rcutoff2))
978 /* Compute parameters for interactions between i and j atoms */
979 qq30 = _mm_mul_pd(iq3,jq0);
981 /* REACTION-FIELD ELECTROSTATICS */
982 felec = _mm_mul_pd(qq30,_mm_msub_pd(rinv30,rinvsq30,krf2));
984 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
988 fscal = _mm_and_pd(fscal,cutoff_mask);
990 /* Update vectorial force */
991 fix3 = _mm_macc_pd(dx30,fscal,fix3);
992 fiy3 = _mm_macc_pd(dy30,fscal,fiy3);
993 fiz3 = _mm_macc_pd(dz30,fscal,fiz3);
995 fjx0 = _mm_macc_pd(dx30,fscal,fjx0);
996 fjy0 = _mm_macc_pd(dy30,fscal,fjy0);
997 fjz0 = _mm_macc_pd(dz30,fscal,fjz0);
1001 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0);
1003 /* Inner loop uses 153 flops */
1006 if(jidx<j_index_end)
1010 j_coord_offsetA = DIM*jnrA;
1012 /* load j atom coordinates */
1013 gmx_mm_load_1rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
1016 /* Calculate displacement vector */
1017 dx00 = _mm_sub_pd(ix0,jx0);
1018 dy00 = _mm_sub_pd(iy0,jy0);
1019 dz00 = _mm_sub_pd(iz0,jz0);
1020 dx10 = _mm_sub_pd(ix1,jx0);
1021 dy10 = _mm_sub_pd(iy1,jy0);
1022 dz10 = _mm_sub_pd(iz1,jz0);
1023 dx20 = _mm_sub_pd(ix2,jx0);
1024 dy20 = _mm_sub_pd(iy2,jy0);
1025 dz20 = _mm_sub_pd(iz2,jz0);
1026 dx30 = _mm_sub_pd(ix3,jx0);
1027 dy30 = _mm_sub_pd(iy3,jy0);
1028 dz30 = _mm_sub_pd(iz3,jz0);
1030 /* Calculate squared distance and things based on it */
1031 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
1032 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
1033 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
1034 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
1036 rinv00 = gmx_mm_invsqrt_pd(rsq00);
1037 rinv10 = gmx_mm_invsqrt_pd(rsq10);
1038 rinv20 = gmx_mm_invsqrt_pd(rsq20);
1039 rinv30 = gmx_mm_invsqrt_pd(rsq30);
1041 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
1042 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
1043 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
1045 /* Load parameters for j particles */
1046 jq0 = _mm_load_sd(charge+jnrA+0);
1047 vdwjidx0A = 2*vdwtype[jnrA+0];
1049 fjx0 = _mm_setzero_pd();
1050 fjy0 = _mm_setzero_pd();
1051 fjz0 = _mm_setzero_pd();
1053 /**************************
1054 * CALCULATE INTERACTIONS *
1055 **************************/
1057 r00 = _mm_mul_pd(rsq00,rinv00);
1059 /* Compute parameters for interactions between i and j atoms */
1060 gmx_mm_load_1pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,&c6_00,&c12_00);
1062 /* Calculate table index by multiplying r with table scale and truncate to integer */
1063 rt = _mm_mul_pd(r00,vftabscale);
1064 vfitab = _mm_cvttpd_epi32(rt);
1066 vfeps = _mm_frcz_pd(rt);
1068 vfeps = _mm_sub_pd(rt,_mm_round_pd(rt, _MM_FROUND_FLOOR));
1070 twovfeps = _mm_add_pd(vfeps,vfeps);
1071 vfitab = _mm_slli_epi32(vfitab,3);
1073 /* CUBIC SPLINE TABLE DISPERSION */
1074 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1075 F = _mm_setzero_pd();
1076 GMX_MM_TRANSPOSE2_PD(Y,F);
1077 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
1078 H = _mm_setzero_pd();
1079 GMX_MM_TRANSPOSE2_PD(G,H);
1080 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
1081 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
1082 fvdw6 = _mm_mul_pd(c6_00,FF);
1084 /* CUBIC SPLINE TABLE REPULSION */
1085 vfitab = _mm_add_epi32(vfitab,ifour);
1086 Y = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1087 F = _mm_setzero_pd();
1088 GMX_MM_TRANSPOSE2_PD(Y,F);
1089 G = _mm_load_pd( vftab + _mm_extract_epi32(vfitab,0) +2);
1090 H = _mm_setzero_pd();
1091 GMX_MM_TRANSPOSE2_PD(G,H);
1092 Fp = _mm_macc_pd(vfeps,_mm_macc_pd(H,vfeps,G),F);
1093 FF = _mm_macc_pd(vfeps,_mm_macc_pd(twovfeps,H,G),Fp);
1094 fvdw12 = _mm_mul_pd(c12_00,FF);
1095 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
1099 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1101 /* Update vectorial force */
1102 fix0 = _mm_macc_pd(dx00,fscal,fix0);
1103 fiy0 = _mm_macc_pd(dy00,fscal,fiy0);
1104 fiz0 = _mm_macc_pd(dz00,fscal,fiz0);
1106 fjx0 = _mm_macc_pd(dx00,fscal,fjx0);
1107 fjy0 = _mm_macc_pd(dy00,fscal,fjy0);
1108 fjz0 = _mm_macc_pd(dz00,fscal,fjz0);
1110 /**************************
1111 * CALCULATE INTERACTIONS *
1112 **************************/
1114 if (gmx_mm_any_lt(rsq10,rcutoff2))
1117 /* Compute parameters for interactions between i and j atoms */
1118 qq10 = _mm_mul_pd(iq1,jq0);
1120 /* REACTION-FIELD ELECTROSTATICS */
1121 felec = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
1123 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
1127 fscal = _mm_and_pd(fscal,cutoff_mask);
1129 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1131 /* Update vectorial force */
1132 fix1 = _mm_macc_pd(dx10,fscal,fix1);
1133 fiy1 = _mm_macc_pd(dy10,fscal,fiy1);
1134 fiz1 = _mm_macc_pd(dz10,fscal,fiz1);
1136 fjx0 = _mm_macc_pd(dx10,fscal,fjx0);
1137 fjy0 = _mm_macc_pd(dy10,fscal,fjy0);
1138 fjz0 = _mm_macc_pd(dz10,fscal,fjz0);
1142 /**************************
1143 * CALCULATE INTERACTIONS *
1144 **************************/
1146 if (gmx_mm_any_lt(rsq20,rcutoff2))
1149 /* Compute parameters for interactions between i and j atoms */
1150 qq20 = _mm_mul_pd(iq2,jq0);
1152 /* REACTION-FIELD ELECTROSTATICS */
1153 felec = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
1155 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
1159 fscal = _mm_and_pd(fscal,cutoff_mask);
1161 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1163 /* Update vectorial force */
1164 fix2 = _mm_macc_pd(dx20,fscal,fix2);
1165 fiy2 = _mm_macc_pd(dy20,fscal,fiy2);
1166 fiz2 = _mm_macc_pd(dz20,fscal,fiz2);
1168 fjx0 = _mm_macc_pd(dx20,fscal,fjx0);
1169 fjy0 = _mm_macc_pd(dy20,fscal,fjy0);
1170 fjz0 = _mm_macc_pd(dz20,fscal,fjz0);
1174 /**************************
1175 * CALCULATE INTERACTIONS *
1176 **************************/
1178 if (gmx_mm_any_lt(rsq30,rcutoff2))
1181 /* Compute parameters for interactions between i and j atoms */
1182 qq30 = _mm_mul_pd(iq3,jq0);
1184 /* REACTION-FIELD ELECTROSTATICS */
1185 felec = _mm_mul_pd(qq30,_mm_msub_pd(rinv30,rinvsq30,krf2));
1187 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
1191 fscal = _mm_and_pd(fscal,cutoff_mask);
1193 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1195 /* Update vectorial force */
1196 fix3 = _mm_macc_pd(dx30,fscal,fix3);
1197 fiy3 = _mm_macc_pd(dy30,fscal,fiy3);
1198 fiz3 = _mm_macc_pd(dz30,fscal,fiz3);
1200 fjx0 = _mm_macc_pd(dx30,fscal,fjx0);
1201 fjy0 = _mm_macc_pd(dy30,fscal,fjy0);
1202 fjz0 = _mm_macc_pd(dz30,fscal,fjz0);
1206 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0);
1208 /* Inner loop uses 153 flops */
1211 /* End of innermost loop */
1213 gmx_mm_update_iforce_4atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
1214 f+i_coord_offset,fshift+i_shift_offset);
1216 /* Increment number of inner iterations */
1217 inneriter += j_index_end - j_index_start;
1219 /* Outer loop uses 24 flops */
1222 /* Increment number of outer iterations */
1225 /* Update outer/inner flops */
1227 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4_F,outeriter*24 + inneriter*153);