2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014,2015,2017, 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 sse2_double kernel generator.
44 #include "../nb_kernel.h"
45 #include "gromacs/gmxlib/nrnb.h"
47 #include "kernelutil_x86_sse2_double.h"
50 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwNone_GeomW3W3_VF_sse2_double
51 * Electrostatics interaction: ReactionField
52 * VdW interaction: None
53 * Geometry: Water3-Water3
54 * Calculate force/pot: PotentialAndForce
57 nb_kernel_ElecRFCut_VdwNone_GeomW3W3_VF_sse2_double
58 (t_nblist * gmx_restrict nlist,
59 rvec * gmx_restrict xx,
60 rvec * gmx_restrict ff,
61 struct t_forcerec * gmx_restrict fr,
62 t_mdatoms * gmx_restrict mdatoms,
63 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
64 t_nrnb * gmx_restrict nrnb)
66 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
67 * just 0 for non-waters.
68 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
69 * jnr indices corresponding to data put in the four positions in the SIMD register.
71 int i_shift_offset,i_coord_offset,outeriter,inneriter;
72 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
74 int j_coord_offsetA,j_coord_offsetB;
75 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
77 real *shiftvec,*fshift,*x,*f;
78 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
80 __m128d ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
82 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
84 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
85 int vdwjidx0A,vdwjidx0B;
86 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
87 int vdwjidx1A,vdwjidx1B;
88 __m128d jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
89 int vdwjidx2A,vdwjidx2B;
90 __m128d jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
91 __m128d dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
92 __m128d dx01,dy01,dz01,rsq01,rinv01,rinvsq01,r01,qq01,c6_01,c12_01;
93 __m128d dx02,dy02,dz02,rsq02,rinv02,rinvsq02,r02,qq02,c6_02,c12_02;
94 __m128d dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
95 __m128d dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
96 __m128d dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
97 __m128d dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
98 __m128d dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
99 __m128d dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
100 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
102 __m128d dummy_mask,cutoff_mask;
103 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
104 __m128d one = _mm_set1_pd(1.0);
105 __m128d two = _mm_set1_pd(2.0);
111 jindex = nlist->jindex;
113 shiftidx = nlist->shift;
115 shiftvec = fr->shift_vec[0];
116 fshift = fr->fshift[0];
117 facel = _mm_set1_pd(fr->ic->epsfac);
118 charge = mdatoms->chargeA;
119 krf = _mm_set1_pd(fr->ic->k_rf);
120 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
121 crf = _mm_set1_pd(fr->ic->c_rf);
123 /* Setup water-specific parameters */
124 inr = nlist->iinr[0];
125 iq0 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+0]));
126 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
127 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
129 jq0 = _mm_set1_pd(charge[inr+0]);
130 jq1 = _mm_set1_pd(charge[inr+1]);
131 jq2 = _mm_set1_pd(charge[inr+2]);
132 qq00 = _mm_mul_pd(iq0,jq0);
133 qq01 = _mm_mul_pd(iq0,jq1);
134 qq02 = _mm_mul_pd(iq0,jq2);
135 qq10 = _mm_mul_pd(iq1,jq0);
136 qq11 = _mm_mul_pd(iq1,jq1);
137 qq12 = _mm_mul_pd(iq1,jq2);
138 qq20 = _mm_mul_pd(iq2,jq0);
139 qq21 = _mm_mul_pd(iq2,jq1);
140 qq22 = _mm_mul_pd(iq2,jq2);
142 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
143 rcutoff_scalar = fr->ic->rcoulomb;
144 rcutoff = _mm_set1_pd(rcutoff_scalar);
145 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
147 /* Avoid stupid compiler warnings */
155 /* Start outer loop over neighborlists */
156 for(iidx=0; iidx<nri; iidx++)
158 /* Load shift vector for this list */
159 i_shift_offset = DIM*shiftidx[iidx];
161 /* Load limits for loop over neighbors */
162 j_index_start = jindex[iidx];
163 j_index_end = jindex[iidx+1];
165 /* Get outer coordinate index */
167 i_coord_offset = DIM*inr;
169 /* Load i particle coords and add shift vector */
170 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
171 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2);
173 fix0 = _mm_setzero_pd();
174 fiy0 = _mm_setzero_pd();
175 fiz0 = _mm_setzero_pd();
176 fix1 = _mm_setzero_pd();
177 fiy1 = _mm_setzero_pd();
178 fiz1 = _mm_setzero_pd();
179 fix2 = _mm_setzero_pd();
180 fiy2 = _mm_setzero_pd();
181 fiz2 = _mm_setzero_pd();
183 /* Reset potential sums */
184 velecsum = _mm_setzero_pd();
186 /* Start inner kernel loop */
187 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
190 /* Get j neighbor index, and coordinate index */
193 j_coord_offsetA = DIM*jnrA;
194 j_coord_offsetB = DIM*jnrB;
196 /* load j atom coordinates */
197 gmx_mm_load_3rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
198 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
200 /* Calculate displacement vector */
201 dx00 = _mm_sub_pd(ix0,jx0);
202 dy00 = _mm_sub_pd(iy0,jy0);
203 dz00 = _mm_sub_pd(iz0,jz0);
204 dx01 = _mm_sub_pd(ix0,jx1);
205 dy01 = _mm_sub_pd(iy0,jy1);
206 dz01 = _mm_sub_pd(iz0,jz1);
207 dx02 = _mm_sub_pd(ix0,jx2);
208 dy02 = _mm_sub_pd(iy0,jy2);
209 dz02 = _mm_sub_pd(iz0,jz2);
210 dx10 = _mm_sub_pd(ix1,jx0);
211 dy10 = _mm_sub_pd(iy1,jy0);
212 dz10 = _mm_sub_pd(iz1,jz0);
213 dx11 = _mm_sub_pd(ix1,jx1);
214 dy11 = _mm_sub_pd(iy1,jy1);
215 dz11 = _mm_sub_pd(iz1,jz1);
216 dx12 = _mm_sub_pd(ix1,jx2);
217 dy12 = _mm_sub_pd(iy1,jy2);
218 dz12 = _mm_sub_pd(iz1,jz2);
219 dx20 = _mm_sub_pd(ix2,jx0);
220 dy20 = _mm_sub_pd(iy2,jy0);
221 dz20 = _mm_sub_pd(iz2,jz0);
222 dx21 = _mm_sub_pd(ix2,jx1);
223 dy21 = _mm_sub_pd(iy2,jy1);
224 dz21 = _mm_sub_pd(iz2,jz1);
225 dx22 = _mm_sub_pd(ix2,jx2);
226 dy22 = _mm_sub_pd(iy2,jy2);
227 dz22 = _mm_sub_pd(iz2,jz2);
229 /* Calculate squared distance and things based on it */
230 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
231 rsq01 = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
232 rsq02 = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
233 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
234 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
235 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
236 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
237 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
238 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
240 rinv00 = sse2_invsqrt_d(rsq00);
241 rinv01 = sse2_invsqrt_d(rsq01);
242 rinv02 = sse2_invsqrt_d(rsq02);
243 rinv10 = sse2_invsqrt_d(rsq10);
244 rinv11 = sse2_invsqrt_d(rsq11);
245 rinv12 = sse2_invsqrt_d(rsq12);
246 rinv20 = sse2_invsqrt_d(rsq20);
247 rinv21 = sse2_invsqrt_d(rsq21);
248 rinv22 = sse2_invsqrt_d(rsq22);
250 rinvsq00 = _mm_mul_pd(rinv00,rinv00);
251 rinvsq01 = _mm_mul_pd(rinv01,rinv01);
252 rinvsq02 = _mm_mul_pd(rinv02,rinv02);
253 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
254 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
255 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
256 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
257 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
258 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
260 fjx0 = _mm_setzero_pd();
261 fjy0 = _mm_setzero_pd();
262 fjz0 = _mm_setzero_pd();
263 fjx1 = _mm_setzero_pd();
264 fjy1 = _mm_setzero_pd();
265 fjz1 = _mm_setzero_pd();
266 fjx2 = _mm_setzero_pd();
267 fjy2 = _mm_setzero_pd();
268 fjz2 = _mm_setzero_pd();
270 /**************************
271 * CALCULATE INTERACTIONS *
272 **************************/
274 if (gmx_mm_any_lt(rsq00,rcutoff2))
277 /* REACTION-FIELD ELECTROSTATICS */
278 velec = _mm_mul_pd(qq00,_mm_sub_pd(_mm_add_pd(rinv00,_mm_mul_pd(krf,rsq00)),crf));
279 felec = _mm_mul_pd(qq00,_mm_sub_pd(_mm_mul_pd(rinv00,rinvsq00),krf2));
281 cutoff_mask = _mm_cmplt_pd(rsq00,rcutoff2);
283 /* Update potential sum for this i atom from the interaction with this j atom. */
284 velec = _mm_and_pd(velec,cutoff_mask);
285 velecsum = _mm_add_pd(velecsum,velec);
289 fscal = _mm_and_pd(fscal,cutoff_mask);
291 /* Calculate temporary vectorial force */
292 tx = _mm_mul_pd(fscal,dx00);
293 ty = _mm_mul_pd(fscal,dy00);
294 tz = _mm_mul_pd(fscal,dz00);
296 /* Update vectorial force */
297 fix0 = _mm_add_pd(fix0,tx);
298 fiy0 = _mm_add_pd(fiy0,ty);
299 fiz0 = _mm_add_pd(fiz0,tz);
301 fjx0 = _mm_add_pd(fjx0,tx);
302 fjy0 = _mm_add_pd(fjy0,ty);
303 fjz0 = _mm_add_pd(fjz0,tz);
307 /**************************
308 * CALCULATE INTERACTIONS *
309 **************************/
311 if (gmx_mm_any_lt(rsq01,rcutoff2))
314 /* REACTION-FIELD ELECTROSTATICS */
315 velec = _mm_mul_pd(qq01,_mm_sub_pd(_mm_add_pd(rinv01,_mm_mul_pd(krf,rsq01)),crf));
316 felec = _mm_mul_pd(qq01,_mm_sub_pd(_mm_mul_pd(rinv01,rinvsq01),krf2));
318 cutoff_mask = _mm_cmplt_pd(rsq01,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 /* Calculate temporary vectorial force */
329 tx = _mm_mul_pd(fscal,dx01);
330 ty = _mm_mul_pd(fscal,dy01);
331 tz = _mm_mul_pd(fscal,dz01);
333 /* Update vectorial force */
334 fix0 = _mm_add_pd(fix0,tx);
335 fiy0 = _mm_add_pd(fiy0,ty);
336 fiz0 = _mm_add_pd(fiz0,tz);
338 fjx1 = _mm_add_pd(fjx1,tx);
339 fjy1 = _mm_add_pd(fjy1,ty);
340 fjz1 = _mm_add_pd(fjz1,tz);
344 /**************************
345 * CALCULATE INTERACTIONS *
346 **************************/
348 if (gmx_mm_any_lt(rsq02,rcutoff2))
351 /* REACTION-FIELD ELECTROSTATICS */
352 velec = _mm_mul_pd(qq02,_mm_sub_pd(_mm_add_pd(rinv02,_mm_mul_pd(krf,rsq02)),crf));
353 felec = _mm_mul_pd(qq02,_mm_sub_pd(_mm_mul_pd(rinv02,rinvsq02),krf2));
355 cutoff_mask = _mm_cmplt_pd(rsq02,rcutoff2);
357 /* Update potential sum for this i atom from the interaction with this j atom. */
358 velec = _mm_and_pd(velec,cutoff_mask);
359 velecsum = _mm_add_pd(velecsum,velec);
363 fscal = _mm_and_pd(fscal,cutoff_mask);
365 /* Calculate temporary vectorial force */
366 tx = _mm_mul_pd(fscal,dx02);
367 ty = _mm_mul_pd(fscal,dy02);
368 tz = _mm_mul_pd(fscal,dz02);
370 /* Update vectorial force */
371 fix0 = _mm_add_pd(fix0,tx);
372 fiy0 = _mm_add_pd(fiy0,ty);
373 fiz0 = _mm_add_pd(fiz0,tz);
375 fjx2 = _mm_add_pd(fjx2,tx);
376 fjy2 = _mm_add_pd(fjy2,ty);
377 fjz2 = _mm_add_pd(fjz2,tz);
381 /**************************
382 * CALCULATE INTERACTIONS *
383 **************************/
385 if (gmx_mm_any_lt(rsq10,rcutoff2))
388 /* REACTION-FIELD ELECTROSTATICS */
389 velec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_add_pd(rinv10,_mm_mul_pd(krf,rsq10)),crf));
390 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
392 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
394 /* Update potential sum for this i atom from the interaction with this j atom. */
395 velec = _mm_and_pd(velec,cutoff_mask);
396 velecsum = _mm_add_pd(velecsum,velec);
400 fscal = _mm_and_pd(fscal,cutoff_mask);
402 /* Calculate temporary vectorial force */
403 tx = _mm_mul_pd(fscal,dx10);
404 ty = _mm_mul_pd(fscal,dy10);
405 tz = _mm_mul_pd(fscal,dz10);
407 /* Update vectorial force */
408 fix1 = _mm_add_pd(fix1,tx);
409 fiy1 = _mm_add_pd(fiy1,ty);
410 fiz1 = _mm_add_pd(fiz1,tz);
412 fjx0 = _mm_add_pd(fjx0,tx);
413 fjy0 = _mm_add_pd(fjy0,ty);
414 fjz0 = _mm_add_pd(fjz0,tz);
418 /**************************
419 * CALCULATE INTERACTIONS *
420 **************************/
422 if (gmx_mm_any_lt(rsq11,rcutoff2))
425 /* REACTION-FIELD ELECTROSTATICS */
426 velec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_add_pd(rinv11,_mm_mul_pd(krf,rsq11)),crf));
427 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
429 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
431 /* Update potential sum for this i atom from the interaction with this j atom. */
432 velec = _mm_and_pd(velec,cutoff_mask);
433 velecsum = _mm_add_pd(velecsum,velec);
437 fscal = _mm_and_pd(fscal,cutoff_mask);
439 /* Calculate temporary vectorial force */
440 tx = _mm_mul_pd(fscal,dx11);
441 ty = _mm_mul_pd(fscal,dy11);
442 tz = _mm_mul_pd(fscal,dz11);
444 /* Update vectorial force */
445 fix1 = _mm_add_pd(fix1,tx);
446 fiy1 = _mm_add_pd(fiy1,ty);
447 fiz1 = _mm_add_pd(fiz1,tz);
449 fjx1 = _mm_add_pd(fjx1,tx);
450 fjy1 = _mm_add_pd(fjy1,ty);
451 fjz1 = _mm_add_pd(fjz1,tz);
455 /**************************
456 * CALCULATE INTERACTIONS *
457 **************************/
459 if (gmx_mm_any_lt(rsq12,rcutoff2))
462 /* REACTION-FIELD ELECTROSTATICS */
463 velec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_add_pd(rinv12,_mm_mul_pd(krf,rsq12)),crf));
464 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
466 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
468 /* Update potential sum for this i atom from the interaction with this j atom. */
469 velec = _mm_and_pd(velec,cutoff_mask);
470 velecsum = _mm_add_pd(velecsum,velec);
474 fscal = _mm_and_pd(fscal,cutoff_mask);
476 /* Calculate temporary vectorial force */
477 tx = _mm_mul_pd(fscal,dx12);
478 ty = _mm_mul_pd(fscal,dy12);
479 tz = _mm_mul_pd(fscal,dz12);
481 /* Update vectorial force */
482 fix1 = _mm_add_pd(fix1,tx);
483 fiy1 = _mm_add_pd(fiy1,ty);
484 fiz1 = _mm_add_pd(fiz1,tz);
486 fjx2 = _mm_add_pd(fjx2,tx);
487 fjy2 = _mm_add_pd(fjy2,ty);
488 fjz2 = _mm_add_pd(fjz2,tz);
492 /**************************
493 * CALCULATE INTERACTIONS *
494 **************************/
496 if (gmx_mm_any_lt(rsq20,rcutoff2))
499 /* REACTION-FIELD ELECTROSTATICS */
500 velec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_add_pd(rinv20,_mm_mul_pd(krf,rsq20)),crf));
501 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
503 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
505 /* Update potential sum for this i atom from the interaction with this j atom. */
506 velec = _mm_and_pd(velec,cutoff_mask);
507 velecsum = _mm_add_pd(velecsum,velec);
511 fscal = _mm_and_pd(fscal,cutoff_mask);
513 /* Calculate temporary vectorial force */
514 tx = _mm_mul_pd(fscal,dx20);
515 ty = _mm_mul_pd(fscal,dy20);
516 tz = _mm_mul_pd(fscal,dz20);
518 /* Update vectorial force */
519 fix2 = _mm_add_pd(fix2,tx);
520 fiy2 = _mm_add_pd(fiy2,ty);
521 fiz2 = _mm_add_pd(fiz2,tz);
523 fjx0 = _mm_add_pd(fjx0,tx);
524 fjy0 = _mm_add_pd(fjy0,ty);
525 fjz0 = _mm_add_pd(fjz0,tz);
529 /**************************
530 * CALCULATE INTERACTIONS *
531 **************************/
533 if (gmx_mm_any_lt(rsq21,rcutoff2))
536 /* REACTION-FIELD ELECTROSTATICS */
537 velec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_add_pd(rinv21,_mm_mul_pd(krf,rsq21)),crf));
538 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
540 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
542 /* Update potential sum for this i atom from the interaction with this j atom. */
543 velec = _mm_and_pd(velec,cutoff_mask);
544 velecsum = _mm_add_pd(velecsum,velec);
548 fscal = _mm_and_pd(fscal,cutoff_mask);
550 /* Calculate temporary vectorial force */
551 tx = _mm_mul_pd(fscal,dx21);
552 ty = _mm_mul_pd(fscal,dy21);
553 tz = _mm_mul_pd(fscal,dz21);
555 /* Update vectorial force */
556 fix2 = _mm_add_pd(fix2,tx);
557 fiy2 = _mm_add_pd(fiy2,ty);
558 fiz2 = _mm_add_pd(fiz2,tz);
560 fjx1 = _mm_add_pd(fjx1,tx);
561 fjy1 = _mm_add_pd(fjy1,ty);
562 fjz1 = _mm_add_pd(fjz1,tz);
566 /**************************
567 * CALCULATE INTERACTIONS *
568 **************************/
570 if (gmx_mm_any_lt(rsq22,rcutoff2))
573 /* REACTION-FIELD ELECTROSTATICS */
574 velec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_add_pd(rinv22,_mm_mul_pd(krf,rsq22)),crf));
575 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
577 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
579 /* Update potential sum for this i atom from the interaction with this j atom. */
580 velec = _mm_and_pd(velec,cutoff_mask);
581 velecsum = _mm_add_pd(velecsum,velec);
585 fscal = _mm_and_pd(fscal,cutoff_mask);
587 /* Calculate temporary vectorial force */
588 tx = _mm_mul_pd(fscal,dx22);
589 ty = _mm_mul_pd(fscal,dy22);
590 tz = _mm_mul_pd(fscal,dz22);
592 /* Update vectorial force */
593 fix2 = _mm_add_pd(fix2,tx);
594 fiy2 = _mm_add_pd(fiy2,ty);
595 fiz2 = _mm_add_pd(fiz2,tz);
597 fjx2 = _mm_add_pd(fjx2,tx);
598 fjy2 = _mm_add_pd(fjy2,ty);
599 fjz2 = _mm_add_pd(fjz2,tz);
603 gmx_mm_decrement_3rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
605 /* Inner loop uses 324 flops */
612 j_coord_offsetA = DIM*jnrA;
614 /* load j atom coordinates */
615 gmx_mm_load_3rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
616 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
618 /* Calculate displacement vector */
619 dx00 = _mm_sub_pd(ix0,jx0);
620 dy00 = _mm_sub_pd(iy0,jy0);
621 dz00 = _mm_sub_pd(iz0,jz0);
622 dx01 = _mm_sub_pd(ix0,jx1);
623 dy01 = _mm_sub_pd(iy0,jy1);
624 dz01 = _mm_sub_pd(iz0,jz1);
625 dx02 = _mm_sub_pd(ix0,jx2);
626 dy02 = _mm_sub_pd(iy0,jy2);
627 dz02 = _mm_sub_pd(iz0,jz2);
628 dx10 = _mm_sub_pd(ix1,jx0);
629 dy10 = _mm_sub_pd(iy1,jy0);
630 dz10 = _mm_sub_pd(iz1,jz0);
631 dx11 = _mm_sub_pd(ix1,jx1);
632 dy11 = _mm_sub_pd(iy1,jy1);
633 dz11 = _mm_sub_pd(iz1,jz1);
634 dx12 = _mm_sub_pd(ix1,jx2);
635 dy12 = _mm_sub_pd(iy1,jy2);
636 dz12 = _mm_sub_pd(iz1,jz2);
637 dx20 = _mm_sub_pd(ix2,jx0);
638 dy20 = _mm_sub_pd(iy2,jy0);
639 dz20 = _mm_sub_pd(iz2,jz0);
640 dx21 = _mm_sub_pd(ix2,jx1);
641 dy21 = _mm_sub_pd(iy2,jy1);
642 dz21 = _mm_sub_pd(iz2,jz1);
643 dx22 = _mm_sub_pd(ix2,jx2);
644 dy22 = _mm_sub_pd(iy2,jy2);
645 dz22 = _mm_sub_pd(iz2,jz2);
647 /* Calculate squared distance and things based on it */
648 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
649 rsq01 = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
650 rsq02 = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
651 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
652 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
653 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
654 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
655 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
656 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
658 rinv00 = sse2_invsqrt_d(rsq00);
659 rinv01 = sse2_invsqrt_d(rsq01);
660 rinv02 = sse2_invsqrt_d(rsq02);
661 rinv10 = sse2_invsqrt_d(rsq10);
662 rinv11 = sse2_invsqrt_d(rsq11);
663 rinv12 = sse2_invsqrt_d(rsq12);
664 rinv20 = sse2_invsqrt_d(rsq20);
665 rinv21 = sse2_invsqrt_d(rsq21);
666 rinv22 = sse2_invsqrt_d(rsq22);
668 rinvsq00 = _mm_mul_pd(rinv00,rinv00);
669 rinvsq01 = _mm_mul_pd(rinv01,rinv01);
670 rinvsq02 = _mm_mul_pd(rinv02,rinv02);
671 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
672 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
673 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
674 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
675 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
676 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
678 fjx0 = _mm_setzero_pd();
679 fjy0 = _mm_setzero_pd();
680 fjz0 = _mm_setzero_pd();
681 fjx1 = _mm_setzero_pd();
682 fjy1 = _mm_setzero_pd();
683 fjz1 = _mm_setzero_pd();
684 fjx2 = _mm_setzero_pd();
685 fjy2 = _mm_setzero_pd();
686 fjz2 = _mm_setzero_pd();
688 /**************************
689 * CALCULATE INTERACTIONS *
690 **************************/
692 if (gmx_mm_any_lt(rsq00,rcutoff2))
695 /* REACTION-FIELD ELECTROSTATICS */
696 velec = _mm_mul_pd(qq00,_mm_sub_pd(_mm_add_pd(rinv00,_mm_mul_pd(krf,rsq00)),crf));
697 felec = _mm_mul_pd(qq00,_mm_sub_pd(_mm_mul_pd(rinv00,rinvsq00),krf2));
699 cutoff_mask = _mm_cmplt_pd(rsq00,rcutoff2);
701 /* Update potential sum for this i atom from the interaction with this j atom. */
702 velec = _mm_and_pd(velec,cutoff_mask);
703 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
704 velecsum = _mm_add_pd(velecsum,velec);
708 fscal = _mm_and_pd(fscal,cutoff_mask);
710 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
712 /* Calculate temporary vectorial force */
713 tx = _mm_mul_pd(fscal,dx00);
714 ty = _mm_mul_pd(fscal,dy00);
715 tz = _mm_mul_pd(fscal,dz00);
717 /* Update vectorial force */
718 fix0 = _mm_add_pd(fix0,tx);
719 fiy0 = _mm_add_pd(fiy0,ty);
720 fiz0 = _mm_add_pd(fiz0,tz);
722 fjx0 = _mm_add_pd(fjx0,tx);
723 fjy0 = _mm_add_pd(fjy0,ty);
724 fjz0 = _mm_add_pd(fjz0,tz);
728 /**************************
729 * CALCULATE INTERACTIONS *
730 **************************/
732 if (gmx_mm_any_lt(rsq01,rcutoff2))
735 /* REACTION-FIELD ELECTROSTATICS */
736 velec = _mm_mul_pd(qq01,_mm_sub_pd(_mm_add_pd(rinv01,_mm_mul_pd(krf,rsq01)),crf));
737 felec = _mm_mul_pd(qq01,_mm_sub_pd(_mm_mul_pd(rinv01,rinvsq01),krf2));
739 cutoff_mask = _mm_cmplt_pd(rsq01,rcutoff2);
741 /* Update potential sum for this i atom from the interaction with this j atom. */
742 velec = _mm_and_pd(velec,cutoff_mask);
743 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
744 velecsum = _mm_add_pd(velecsum,velec);
748 fscal = _mm_and_pd(fscal,cutoff_mask);
750 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
752 /* Calculate temporary vectorial force */
753 tx = _mm_mul_pd(fscal,dx01);
754 ty = _mm_mul_pd(fscal,dy01);
755 tz = _mm_mul_pd(fscal,dz01);
757 /* Update vectorial force */
758 fix0 = _mm_add_pd(fix0,tx);
759 fiy0 = _mm_add_pd(fiy0,ty);
760 fiz0 = _mm_add_pd(fiz0,tz);
762 fjx1 = _mm_add_pd(fjx1,tx);
763 fjy1 = _mm_add_pd(fjy1,ty);
764 fjz1 = _mm_add_pd(fjz1,tz);
768 /**************************
769 * CALCULATE INTERACTIONS *
770 **************************/
772 if (gmx_mm_any_lt(rsq02,rcutoff2))
775 /* REACTION-FIELD ELECTROSTATICS */
776 velec = _mm_mul_pd(qq02,_mm_sub_pd(_mm_add_pd(rinv02,_mm_mul_pd(krf,rsq02)),crf));
777 felec = _mm_mul_pd(qq02,_mm_sub_pd(_mm_mul_pd(rinv02,rinvsq02),krf2));
779 cutoff_mask = _mm_cmplt_pd(rsq02,rcutoff2);
781 /* Update potential sum for this i atom from the interaction with this j atom. */
782 velec = _mm_and_pd(velec,cutoff_mask);
783 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
784 velecsum = _mm_add_pd(velecsum,velec);
788 fscal = _mm_and_pd(fscal,cutoff_mask);
790 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
792 /* Calculate temporary vectorial force */
793 tx = _mm_mul_pd(fscal,dx02);
794 ty = _mm_mul_pd(fscal,dy02);
795 tz = _mm_mul_pd(fscal,dz02);
797 /* Update vectorial force */
798 fix0 = _mm_add_pd(fix0,tx);
799 fiy0 = _mm_add_pd(fiy0,ty);
800 fiz0 = _mm_add_pd(fiz0,tz);
802 fjx2 = _mm_add_pd(fjx2,tx);
803 fjy2 = _mm_add_pd(fjy2,ty);
804 fjz2 = _mm_add_pd(fjz2,tz);
808 /**************************
809 * CALCULATE INTERACTIONS *
810 **************************/
812 if (gmx_mm_any_lt(rsq10,rcutoff2))
815 /* REACTION-FIELD ELECTROSTATICS */
816 velec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_add_pd(rinv10,_mm_mul_pd(krf,rsq10)),crf));
817 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
819 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
821 /* Update potential sum for this i atom from the interaction with this j atom. */
822 velec = _mm_and_pd(velec,cutoff_mask);
823 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
824 velecsum = _mm_add_pd(velecsum,velec);
828 fscal = _mm_and_pd(fscal,cutoff_mask);
830 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
832 /* Calculate temporary vectorial force */
833 tx = _mm_mul_pd(fscal,dx10);
834 ty = _mm_mul_pd(fscal,dy10);
835 tz = _mm_mul_pd(fscal,dz10);
837 /* Update vectorial force */
838 fix1 = _mm_add_pd(fix1,tx);
839 fiy1 = _mm_add_pd(fiy1,ty);
840 fiz1 = _mm_add_pd(fiz1,tz);
842 fjx0 = _mm_add_pd(fjx0,tx);
843 fjy0 = _mm_add_pd(fjy0,ty);
844 fjz0 = _mm_add_pd(fjz0,tz);
848 /**************************
849 * CALCULATE INTERACTIONS *
850 **************************/
852 if (gmx_mm_any_lt(rsq11,rcutoff2))
855 /* REACTION-FIELD ELECTROSTATICS */
856 velec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_add_pd(rinv11,_mm_mul_pd(krf,rsq11)),crf));
857 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
859 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
861 /* Update potential sum for this i atom from the interaction with this j atom. */
862 velec = _mm_and_pd(velec,cutoff_mask);
863 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
864 velecsum = _mm_add_pd(velecsum,velec);
868 fscal = _mm_and_pd(fscal,cutoff_mask);
870 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
872 /* Calculate temporary vectorial force */
873 tx = _mm_mul_pd(fscal,dx11);
874 ty = _mm_mul_pd(fscal,dy11);
875 tz = _mm_mul_pd(fscal,dz11);
877 /* Update vectorial force */
878 fix1 = _mm_add_pd(fix1,tx);
879 fiy1 = _mm_add_pd(fiy1,ty);
880 fiz1 = _mm_add_pd(fiz1,tz);
882 fjx1 = _mm_add_pd(fjx1,tx);
883 fjy1 = _mm_add_pd(fjy1,ty);
884 fjz1 = _mm_add_pd(fjz1,tz);
888 /**************************
889 * CALCULATE INTERACTIONS *
890 **************************/
892 if (gmx_mm_any_lt(rsq12,rcutoff2))
895 /* REACTION-FIELD ELECTROSTATICS */
896 velec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_add_pd(rinv12,_mm_mul_pd(krf,rsq12)),crf));
897 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
899 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
901 /* Update potential sum for this i atom from the interaction with this j atom. */
902 velec = _mm_and_pd(velec,cutoff_mask);
903 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
904 velecsum = _mm_add_pd(velecsum,velec);
908 fscal = _mm_and_pd(fscal,cutoff_mask);
910 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
912 /* Calculate temporary vectorial force */
913 tx = _mm_mul_pd(fscal,dx12);
914 ty = _mm_mul_pd(fscal,dy12);
915 tz = _mm_mul_pd(fscal,dz12);
917 /* Update vectorial force */
918 fix1 = _mm_add_pd(fix1,tx);
919 fiy1 = _mm_add_pd(fiy1,ty);
920 fiz1 = _mm_add_pd(fiz1,tz);
922 fjx2 = _mm_add_pd(fjx2,tx);
923 fjy2 = _mm_add_pd(fjy2,ty);
924 fjz2 = _mm_add_pd(fjz2,tz);
928 /**************************
929 * CALCULATE INTERACTIONS *
930 **************************/
932 if (gmx_mm_any_lt(rsq20,rcutoff2))
935 /* REACTION-FIELD ELECTROSTATICS */
936 velec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_add_pd(rinv20,_mm_mul_pd(krf,rsq20)),crf));
937 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
939 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
941 /* Update potential sum for this i atom from the interaction with this j atom. */
942 velec = _mm_and_pd(velec,cutoff_mask);
943 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
944 velecsum = _mm_add_pd(velecsum,velec);
948 fscal = _mm_and_pd(fscal,cutoff_mask);
950 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
952 /* Calculate temporary vectorial force */
953 tx = _mm_mul_pd(fscal,dx20);
954 ty = _mm_mul_pd(fscal,dy20);
955 tz = _mm_mul_pd(fscal,dz20);
957 /* Update vectorial force */
958 fix2 = _mm_add_pd(fix2,tx);
959 fiy2 = _mm_add_pd(fiy2,ty);
960 fiz2 = _mm_add_pd(fiz2,tz);
962 fjx0 = _mm_add_pd(fjx0,tx);
963 fjy0 = _mm_add_pd(fjy0,ty);
964 fjz0 = _mm_add_pd(fjz0,tz);
968 /**************************
969 * CALCULATE INTERACTIONS *
970 **************************/
972 if (gmx_mm_any_lt(rsq21,rcutoff2))
975 /* REACTION-FIELD ELECTROSTATICS */
976 velec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_add_pd(rinv21,_mm_mul_pd(krf,rsq21)),crf));
977 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
979 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
981 /* Update potential sum for this i atom from the interaction with this j atom. */
982 velec = _mm_and_pd(velec,cutoff_mask);
983 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
984 velecsum = _mm_add_pd(velecsum,velec);
988 fscal = _mm_and_pd(fscal,cutoff_mask);
990 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
992 /* Calculate temporary vectorial force */
993 tx = _mm_mul_pd(fscal,dx21);
994 ty = _mm_mul_pd(fscal,dy21);
995 tz = _mm_mul_pd(fscal,dz21);
997 /* Update vectorial force */
998 fix2 = _mm_add_pd(fix2,tx);
999 fiy2 = _mm_add_pd(fiy2,ty);
1000 fiz2 = _mm_add_pd(fiz2,tz);
1002 fjx1 = _mm_add_pd(fjx1,tx);
1003 fjy1 = _mm_add_pd(fjy1,ty);
1004 fjz1 = _mm_add_pd(fjz1,tz);
1008 /**************************
1009 * CALCULATE INTERACTIONS *
1010 **************************/
1012 if (gmx_mm_any_lt(rsq22,rcutoff2))
1015 /* REACTION-FIELD ELECTROSTATICS */
1016 velec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_add_pd(rinv22,_mm_mul_pd(krf,rsq22)),crf));
1017 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
1019 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
1021 /* Update potential sum for this i atom from the interaction with this j atom. */
1022 velec = _mm_and_pd(velec,cutoff_mask);
1023 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
1024 velecsum = _mm_add_pd(velecsum,velec);
1028 fscal = _mm_and_pd(fscal,cutoff_mask);
1030 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1032 /* Calculate temporary vectorial force */
1033 tx = _mm_mul_pd(fscal,dx22);
1034 ty = _mm_mul_pd(fscal,dy22);
1035 tz = _mm_mul_pd(fscal,dz22);
1037 /* Update vectorial force */
1038 fix2 = _mm_add_pd(fix2,tx);
1039 fiy2 = _mm_add_pd(fiy2,ty);
1040 fiz2 = _mm_add_pd(fiz2,tz);
1042 fjx2 = _mm_add_pd(fjx2,tx);
1043 fjy2 = _mm_add_pd(fjy2,ty);
1044 fjz2 = _mm_add_pd(fjz2,tz);
1048 gmx_mm_decrement_3rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
1050 /* Inner loop uses 324 flops */
1053 /* End of innermost loop */
1055 gmx_mm_update_iforce_3atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,
1056 f+i_coord_offset,fshift+i_shift_offset);
1059 /* Update potential energies */
1060 gmx_mm_update_1pot_pd(velecsum,kernel_data->energygrp_elec+ggid);
1062 /* Increment number of inner iterations */
1063 inneriter += j_index_end - j_index_start;
1065 /* Outer loop uses 19 flops */
1068 /* Increment number of outer iterations */
1071 /* Update outer/inner flops */
1073 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W3W3_VF,outeriter*19 + inneriter*324);
1076 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwNone_GeomW3W3_F_sse2_double
1077 * Electrostatics interaction: ReactionField
1078 * VdW interaction: None
1079 * Geometry: Water3-Water3
1080 * Calculate force/pot: Force
1083 nb_kernel_ElecRFCut_VdwNone_GeomW3W3_F_sse2_double
1084 (t_nblist * gmx_restrict nlist,
1085 rvec * gmx_restrict xx,
1086 rvec * gmx_restrict ff,
1087 struct t_forcerec * gmx_restrict fr,
1088 t_mdatoms * gmx_restrict mdatoms,
1089 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
1090 t_nrnb * gmx_restrict nrnb)
1092 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
1093 * just 0 for non-waters.
1094 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
1095 * jnr indices corresponding to data put in the four positions in the SIMD register.
1097 int i_shift_offset,i_coord_offset,outeriter,inneriter;
1098 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
1100 int j_coord_offsetA,j_coord_offsetB;
1101 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
1102 real rcutoff_scalar;
1103 real *shiftvec,*fshift,*x,*f;
1104 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
1106 __m128d ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
1108 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
1110 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
1111 int vdwjidx0A,vdwjidx0B;
1112 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
1113 int vdwjidx1A,vdwjidx1B;
1114 __m128d jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
1115 int vdwjidx2A,vdwjidx2B;
1116 __m128d jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
1117 __m128d dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
1118 __m128d dx01,dy01,dz01,rsq01,rinv01,rinvsq01,r01,qq01,c6_01,c12_01;
1119 __m128d dx02,dy02,dz02,rsq02,rinv02,rinvsq02,r02,qq02,c6_02,c12_02;
1120 __m128d dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
1121 __m128d dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
1122 __m128d dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
1123 __m128d dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
1124 __m128d dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
1125 __m128d dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
1126 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
1128 __m128d dummy_mask,cutoff_mask;
1129 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
1130 __m128d one = _mm_set1_pd(1.0);
1131 __m128d two = _mm_set1_pd(2.0);
1137 jindex = nlist->jindex;
1139 shiftidx = nlist->shift;
1141 shiftvec = fr->shift_vec[0];
1142 fshift = fr->fshift[0];
1143 facel = _mm_set1_pd(fr->ic->epsfac);
1144 charge = mdatoms->chargeA;
1145 krf = _mm_set1_pd(fr->ic->k_rf);
1146 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
1147 crf = _mm_set1_pd(fr->ic->c_rf);
1149 /* Setup water-specific parameters */
1150 inr = nlist->iinr[0];
1151 iq0 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+0]));
1152 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
1153 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
1155 jq0 = _mm_set1_pd(charge[inr+0]);
1156 jq1 = _mm_set1_pd(charge[inr+1]);
1157 jq2 = _mm_set1_pd(charge[inr+2]);
1158 qq00 = _mm_mul_pd(iq0,jq0);
1159 qq01 = _mm_mul_pd(iq0,jq1);
1160 qq02 = _mm_mul_pd(iq0,jq2);
1161 qq10 = _mm_mul_pd(iq1,jq0);
1162 qq11 = _mm_mul_pd(iq1,jq1);
1163 qq12 = _mm_mul_pd(iq1,jq2);
1164 qq20 = _mm_mul_pd(iq2,jq0);
1165 qq21 = _mm_mul_pd(iq2,jq1);
1166 qq22 = _mm_mul_pd(iq2,jq2);
1168 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
1169 rcutoff_scalar = fr->ic->rcoulomb;
1170 rcutoff = _mm_set1_pd(rcutoff_scalar);
1171 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
1173 /* Avoid stupid compiler warnings */
1175 j_coord_offsetA = 0;
1176 j_coord_offsetB = 0;
1181 /* Start outer loop over neighborlists */
1182 for(iidx=0; iidx<nri; iidx++)
1184 /* Load shift vector for this list */
1185 i_shift_offset = DIM*shiftidx[iidx];
1187 /* Load limits for loop over neighbors */
1188 j_index_start = jindex[iidx];
1189 j_index_end = jindex[iidx+1];
1191 /* Get outer coordinate index */
1193 i_coord_offset = DIM*inr;
1195 /* Load i particle coords and add shift vector */
1196 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
1197 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2);
1199 fix0 = _mm_setzero_pd();
1200 fiy0 = _mm_setzero_pd();
1201 fiz0 = _mm_setzero_pd();
1202 fix1 = _mm_setzero_pd();
1203 fiy1 = _mm_setzero_pd();
1204 fiz1 = _mm_setzero_pd();
1205 fix2 = _mm_setzero_pd();
1206 fiy2 = _mm_setzero_pd();
1207 fiz2 = _mm_setzero_pd();
1209 /* Start inner kernel loop */
1210 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
1213 /* Get j neighbor index, and coordinate index */
1215 jnrB = jjnr[jidx+1];
1216 j_coord_offsetA = DIM*jnrA;
1217 j_coord_offsetB = DIM*jnrB;
1219 /* load j atom coordinates */
1220 gmx_mm_load_3rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
1221 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
1223 /* Calculate displacement vector */
1224 dx00 = _mm_sub_pd(ix0,jx0);
1225 dy00 = _mm_sub_pd(iy0,jy0);
1226 dz00 = _mm_sub_pd(iz0,jz0);
1227 dx01 = _mm_sub_pd(ix0,jx1);
1228 dy01 = _mm_sub_pd(iy0,jy1);
1229 dz01 = _mm_sub_pd(iz0,jz1);
1230 dx02 = _mm_sub_pd(ix0,jx2);
1231 dy02 = _mm_sub_pd(iy0,jy2);
1232 dz02 = _mm_sub_pd(iz0,jz2);
1233 dx10 = _mm_sub_pd(ix1,jx0);
1234 dy10 = _mm_sub_pd(iy1,jy0);
1235 dz10 = _mm_sub_pd(iz1,jz0);
1236 dx11 = _mm_sub_pd(ix1,jx1);
1237 dy11 = _mm_sub_pd(iy1,jy1);
1238 dz11 = _mm_sub_pd(iz1,jz1);
1239 dx12 = _mm_sub_pd(ix1,jx2);
1240 dy12 = _mm_sub_pd(iy1,jy2);
1241 dz12 = _mm_sub_pd(iz1,jz2);
1242 dx20 = _mm_sub_pd(ix2,jx0);
1243 dy20 = _mm_sub_pd(iy2,jy0);
1244 dz20 = _mm_sub_pd(iz2,jz0);
1245 dx21 = _mm_sub_pd(ix2,jx1);
1246 dy21 = _mm_sub_pd(iy2,jy1);
1247 dz21 = _mm_sub_pd(iz2,jz1);
1248 dx22 = _mm_sub_pd(ix2,jx2);
1249 dy22 = _mm_sub_pd(iy2,jy2);
1250 dz22 = _mm_sub_pd(iz2,jz2);
1252 /* Calculate squared distance and things based on it */
1253 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
1254 rsq01 = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
1255 rsq02 = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
1256 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
1257 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
1258 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
1259 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
1260 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
1261 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
1263 rinv00 = sse2_invsqrt_d(rsq00);
1264 rinv01 = sse2_invsqrt_d(rsq01);
1265 rinv02 = sse2_invsqrt_d(rsq02);
1266 rinv10 = sse2_invsqrt_d(rsq10);
1267 rinv11 = sse2_invsqrt_d(rsq11);
1268 rinv12 = sse2_invsqrt_d(rsq12);
1269 rinv20 = sse2_invsqrt_d(rsq20);
1270 rinv21 = sse2_invsqrt_d(rsq21);
1271 rinv22 = sse2_invsqrt_d(rsq22);
1273 rinvsq00 = _mm_mul_pd(rinv00,rinv00);
1274 rinvsq01 = _mm_mul_pd(rinv01,rinv01);
1275 rinvsq02 = _mm_mul_pd(rinv02,rinv02);
1276 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
1277 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
1278 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
1279 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
1280 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
1281 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
1283 fjx0 = _mm_setzero_pd();
1284 fjy0 = _mm_setzero_pd();
1285 fjz0 = _mm_setzero_pd();
1286 fjx1 = _mm_setzero_pd();
1287 fjy1 = _mm_setzero_pd();
1288 fjz1 = _mm_setzero_pd();
1289 fjx2 = _mm_setzero_pd();
1290 fjy2 = _mm_setzero_pd();
1291 fjz2 = _mm_setzero_pd();
1293 /**************************
1294 * CALCULATE INTERACTIONS *
1295 **************************/
1297 if (gmx_mm_any_lt(rsq00,rcutoff2))
1300 /* REACTION-FIELD ELECTROSTATICS */
1301 felec = _mm_mul_pd(qq00,_mm_sub_pd(_mm_mul_pd(rinv00,rinvsq00),krf2));
1303 cutoff_mask = _mm_cmplt_pd(rsq00,rcutoff2);
1307 fscal = _mm_and_pd(fscal,cutoff_mask);
1309 /* Calculate temporary vectorial force */
1310 tx = _mm_mul_pd(fscal,dx00);
1311 ty = _mm_mul_pd(fscal,dy00);
1312 tz = _mm_mul_pd(fscal,dz00);
1314 /* Update vectorial force */
1315 fix0 = _mm_add_pd(fix0,tx);
1316 fiy0 = _mm_add_pd(fiy0,ty);
1317 fiz0 = _mm_add_pd(fiz0,tz);
1319 fjx0 = _mm_add_pd(fjx0,tx);
1320 fjy0 = _mm_add_pd(fjy0,ty);
1321 fjz0 = _mm_add_pd(fjz0,tz);
1325 /**************************
1326 * CALCULATE INTERACTIONS *
1327 **************************/
1329 if (gmx_mm_any_lt(rsq01,rcutoff2))
1332 /* REACTION-FIELD ELECTROSTATICS */
1333 felec = _mm_mul_pd(qq01,_mm_sub_pd(_mm_mul_pd(rinv01,rinvsq01),krf2));
1335 cutoff_mask = _mm_cmplt_pd(rsq01,rcutoff2);
1339 fscal = _mm_and_pd(fscal,cutoff_mask);
1341 /* Calculate temporary vectorial force */
1342 tx = _mm_mul_pd(fscal,dx01);
1343 ty = _mm_mul_pd(fscal,dy01);
1344 tz = _mm_mul_pd(fscal,dz01);
1346 /* Update vectorial force */
1347 fix0 = _mm_add_pd(fix0,tx);
1348 fiy0 = _mm_add_pd(fiy0,ty);
1349 fiz0 = _mm_add_pd(fiz0,tz);
1351 fjx1 = _mm_add_pd(fjx1,tx);
1352 fjy1 = _mm_add_pd(fjy1,ty);
1353 fjz1 = _mm_add_pd(fjz1,tz);
1357 /**************************
1358 * CALCULATE INTERACTIONS *
1359 **************************/
1361 if (gmx_mm_any_lt(rsq02,rcutoff2))
1364 /* REACTION-FIELD ELECTROSTATICS */
1365 felec = _mm_mul_pd(qq02,_mm_sub_pd(_mm_mul_pd(rinv02,rinvsq02),krf2));
1367 cutoff_mask = _mm_cmplt_pd(rsq02,rcutoff2);
1371 fscal = _mm_and_pd(fscal,cutoff_mask);
1373 /* Calculate temporary vectorial force */
1374 tx = _mm_mul_pd(fscal,dx02);
1375 ty = _mm_mul_pd(fscal,dy02);
1376 tz = _mm_mul_pd(fscal,dz02);
1378 /* Update vectorial force */
1379 fix0 = _mm_add_pd(fix0,tx);
1380 fiy0 = _mm_add_pd(fiy0,ty);
1381 fiz0 = _mm_add_pd(fiz0,tz);
1383 fjx2 = _mm_add_pd(fjx2,tx);
1384 fjy2 = _mm_add_pd(fjy2,ty);
1385 fjz2 = _mm_add_pd(fjz2,tz);
1389 /**************************
1390 * CALCULATE INTERACTIONS *
1391 **************************/
1393 if (gmx_mm_any_lt(rsq10,rcutoff2))
1396 /* REACTION-FIELD ELECTROSTATICS */
1397 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
1399 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
1403 fscal = _mm_and_pd(fscal,cutoff_mask);
1405 /* Calculate temporary vectorial force */
1406 tx = _mm_mul_pd(fscal,dx10);
1407 ty = _mm_mul_pd(fscal,dy10);
1408 tz = _mm_mul_pd(fscal,dz10);
1410 /* Update vectorial force */
1411 fix1 = _mm_add_pd(fix1,tx);
1412 fiy1 = _mm_add_pd(fiy1,ty);
1413 fiz1 = _mm_add_pd(fiz1,tz);
1415 fjx0 = _mm_add_pd(fjx0,tx);
1416 fjy0 = _mm_add_pd(fjy0,ty);
1417 fjz0 = _mm_add_pd(fjz0,tz);
1421 /**************************
1422 * CALCULATE INTERACTIONS *
1423 **************************/
1425 if (gmx_mm_any_lt(rsq11,rcutoff2))
1428 /* REACTION-FIELD ELECTROSTATICS */
1429 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
1431 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
1435 fscal = _mm_and_pd(fscal,cutoff_mask);
1437 /* Calculate temporary vectorial force */
1438 tx = _mm_mul_pd(fscal,dx11);
1439 ty = _mm_mul_pd(fscal,dy11);
1440 tz = _mm_mul_pd(fscal,dz11);
1442 /* Update vectorial force */
1443 fix1 = _mm_add_pd(fix1,tx);
1444 fiy1 = _mm_add_pd(fiy1,ty);
1445 fiz1 = _mm_add_pd(fiz1,tz);
1447 fjx1 = _mm_add_pd(fjx1,tx);
1448 fjy1 = _mm_add_pd(fjy1,ty);
1449 fjz1 = _mm_add_pd(fjz1,tz);
1453 /**************************
1454 * CALCULATE INTERACTIONS *
1455 **************************/
1457 if (gmx_mm_any_lt(rsq12,rcutoff2))
1460 /* REACTION-FIELD ELECTROSTATICS */
1461 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
1463 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
1467 fscal = _mm_and_pd(fscal,cutoff_mask);
1469 /* Calculate temporary vectorial force */
1470 tx = _mm_mul_pd(fscal,dx12);
1471 ty = _mm_mul_pd(fscal,dy12);
1472 tz = _mm_mul_pd(fscal,dz12);
1474 /* Update vectorial force */
1475 fix1 = _mm_add_pd(fix1,tx);
1476 fiy1 = _mm_add_pd(fiy1,ty);
1477 fiz1 = _mm_add_pd(fiz1,tz);
1479 fjx2 = _mm_add_pd(fjx2,tx);
1480 fjy2 = _mm_add_pd(fjy2,ty);
1481 fjz2 = _mm_add_pd(fjz2,tz);
1485 /**************************
1486 * CALCULATE INTERACTIONS *
1487 **************************/
1489 if (gmx_mm_any_lt(rsq20,rcutoff2))
1492 /* REACTION-FIELD ELECTROSTATICS */
1493 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
1495 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
1499 fscal = _mm_and_pd(fscal,cutoff_mask);
1501 /* Calculate temporary vectorial force */
1502 tx = _mm_mul_pd(fscal,dx20);
1503 ty = _mm_mul_pd(fscal,dy20);
1504 tz = _mm_mul_pd(fscal,dz20);
1506 /* Update vectorial force */
1507 fix2 = _mm_add_pd(fix2,tx);
1508 fiy2 = _mm_add_pd(fiy2,ty);
1509 fiz2 = _mm_add_pd(fiz2,tz);
1511 fjx0 = _mm_add_pd(fjx0,tx);
1512 fjy0 = _mm_add_pd(fjy0,ty);
1513 fjz0 = _mm_add_pd(fjz0,tz);
1517 /**************************
1518 * CALCULATE INTERACTIONS *
1519 **************************/
1521 if (gmx_mm_any_lt(rsq21,rcutoff2))
1524 /* REACTION-FIELD ELECTROSTATICS */
1525 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
1527 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
1531 fscal = _mm_and_pd(fscal,cutoff_mask);
1533 /* Calculate temporary vectorial force */
1534 tx = _mm_mul_pd(fscal,dx21);
1535 ty = _mm_mul_pd(fscal,dy21);
1536 tz = _mm_mul_pd(fscal,dz21);
1538 /* Update vectorial force */
1539 fix2 = _mm_add_pd(fix2,tx);
1540 fiy2 = _mm_add_pd(fiy2,ty);
1541 fiz2 = _mm_add_pd(fiz2,tz);
1543 fjx1 = _mm_add_pd(fjx1,tx);
1544 fjy1 = _mm_add_pd(fjy1,ty);
1545 fjz1 = _mm_add_pd(fjz1,tz);
1549 /**************************
1550 * CALCULATE INTERACTIONS *
1551 **************************/
1553 if (gmx_mm_any_lt(rsq22,rcutoff2))
1556 /* REACTION-FIELD ELECTROSTATICS */
1557 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
1559 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
1563 fscal = _mm_and_pd(fscal,cutoff_mask);
1565 /* Calculate temporary vectorial force */
1566 tx = _mm_mul_pd(fscal,dx22);
1567 ty = _mm_mul_pd(fscal,dy22);
1568 tz = _mm_mul_pd(fscal,dz22);
1570 /* Update vectorial force */
1571 fix2 = _mm_add_pd(fix2,tx);
1572 fiy2 = _mm_add_pd(fiy2,ty);
1573 fiz2 = _mm_add_pd(fiz2,tz);
1575 fjx2 = _mm_add_pd(fjx2,tx);
1576 fjy2 = _mm_add_pd(fjy2,ty);
1577 fjz2 = _mm_add_pd(fjz2,tz);
1581 gmx_mm_decrement_3rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
1583 /* Inner loop uses 270 flops */
1586 if(jidx<j_index_end)
1590 j_coord_offsetA = DIM*jnrA;
1592 /* load j atom coordinates */
1593 gmx_mm_load_3rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
1594 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
1596 /* Calculate displacement vector */
1597 dx00 = _mm_sub_pd(ix0,jx0);
1598 dy00 = _mm_sub_pd(iy0,jy0);
1599 dz00 = _mm_sub_pd(iz0,jz0);
1600 dx01 = _mm_sub_pd(ix0,jx1);
1601 dy01 = _mm_sub_pd(iy0,jy1);
1602 dz01 = _mm_sub_pd(iz0,jz1);
1603 dx02 = _mm_sub_pd(ix0,jx2);
1604 dy02 = _mm_sub_pd(iy0,jy2);
1605 dz02 = _mm_sub_pd(iz0,jz2);
1606 dx10 = _mm_sub_pd(ix1,jx0);
1607 dy10 = _mm_sub_pd(iy1,jy0);
1608 dz10 = _mm_sub_pd(iz1,jz0);
1609 dx11 = _mm_sub_pd(ix1,jx1);
1610 dy11 = _mm_sub_pd(iy1,jy1);
1611 dz11 = _mm_sub_pd(iz1,jz1);
1612 dx12 = _mm_sub_pd(ix1,jx2);
1613 dy12 = _mm_sub_pd(iy1,jy2);
1614 dz12 = _mm_sub_pd(iz1,jz2);
1615 dx20 = _mm_sub_pd(ix2,jx0);
1616 dy20 = _mm_sub_pd(iy2,jy0);
1617 dz20 = _mm_sub_pd(iz2,jz0);
1618 dx21 = _mm_sub_pd(ix2,jx1);
1619 dy21 = _mm_sub_pd(iy2,jy1);
1620 dz21 = _mm_sub_pd(iz2,jz1);
1621 dx22 = _mm_sub_pd(ix2,jx2);
1622 dy22 = _mm_sub_pd(iy2,jy2);
1623 dz22 = _mm_sub_pd(iz2,jz2);
1625 /* Calculate squared distance and things based on it */
1626 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
1627 rsq01 = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
1628 rsq02 = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
1629 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
1630 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
1631 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
1632 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
1633 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
1634 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
1636 rinv00 = sse2_invsqrt_d(rsq00);
1637 rinv01 = sse2_invsqrt_d(rsq01);
1638 rinv02 = sse2_invsqrt_d(rsq02);
1639 rinv10 = sse2_invsqrt_d(rsq10);
1640 rinv11 = sse2_invsqrt_d(rsq11);
1641 rinv12 = sse2_invsqrt_d(rsq12);
1642 rinv20 = sse2_invsqrt_d(rsq20);
1643 rinv21 = sse2_invsqrt_d(rsq21);
1644 rinv22 = sse2_invsqrt_d(rsq22);
1646 rinvsq00 = _mm_mul_pd(rinv00,rinv00);
1647 rinvsq01 = _mm_mul_pd(rinv01,rinv01);
1648 rinvsq02 = _mm_mul_pd(rinv02,rinv02);
1649 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
1650 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
1651 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
1652 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
1653 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
1654 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
1656 fjx0 = _mm_setzero_pd();
1657 fjy0 = _mm_setzero_pd();
1658 fjz0 = _mm_setzero_pd();
1659 fjx1 = _mm_setzero_pd();
1660 fjy1 = _mm_setzero_pd();
1661 fjz1 = _mm_setzero_pd();
1662 fjx2 = _mm_setzero_pd();
1663 fjy2 = _mm_setzero_pd();
1664 fjz2 = _mm_setzero_pd();
1666 /**************************
1667 * CALCULATE INTERACTIONS *
1668 **************************/
1670 if (gmx_mm_any_lt(rsq00,rcutoff2))
1673 /* REACTION-FIELD ELECTROSTATICS */
1674 felec = _mm_mul_pd(qq00,_mm_sub_pd(_mm_mul_pd(rinv00,rinvsq00),krf2));
1676 cutoff_mask = _mm_cmplt_pd(rsq00,rcutoff2);
1680 fscal = _mm_and_pd(fscal,cutoff_mask);
1682 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1684 /* Calculate temporary vectorial force */
1685 tx = _mm_mul_pd(fscal,dx00);
1686 ty = _mm_mul_pd(fscal,dy00);
1687 tz = _mm_mul_pd(fscal,dz00);
1689 /* Update vectorial force */
1690 fix0 = _mm_add_pd(fix0,tx);
1691 fiy0 = _mm_add_pd(fiy0,ty);
1692 fiz0 = _mm_add_pd(fiz0,tz);
1694 fjx0 = _mm_add_pd(fjx0,tx);
1695 fjy0 = _mm_add_pd(fjy0,ty);
1696 fjz0 = _mm_add_pd(fjz0,tz);
1700 /**************************
1701 * CALCULATE INTERACTIONS *
1702 **************************/
1704 if (gmx_mm_any_lt(rsq01,rcutoff2))
1707 /* REACTION-FIELD ELECTROSTATICS */
1708 felec = _mm_mul_pd(qq01,_mm_sub_pd(_mm_mul_pd(rinv01,rinvsq01),krf2));
1710 cutoff_mask = _mm_cmplt_pd(rsq01,rcutoff2);
1714 fscal = _mm_and_pd(fscal,cutoff_mask);
1716 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1718 /* Calculate temporary vectorial force */
1719 tx = _mm_mul_pd(fscal,dx01);
1720 ty = _mm_mul_pd(fscal,dy01);
1721 tz = _mm_mul_pd(fscal,dz01);
1723 /* Update vectorial force */
1724 fix0 = _mm_add_pd(fix0,tx);
1725 fiy0 = _mm_add_pd(fiy0,ty);
1726 fiz0 = _mm_add_pd(fiz0,tz);
1728 fjx1 = _mm_add_pd(fjx1,tx);
1729 fjy1 = _mm_add_pd(fjy1,ty);
1730 fjz1 = _mm_add_pd(fjz1,tz);
1734 /**************************
1735 * CALCULATE INTERACTIONS *
1736 **************************/
1738 if (gmx_mm_any_lt(rsq02,rcutoff2))
1741 /* REACTION-FIELD ELECTROSTATICS */
1742 felec = _mm_mul_pd(qq02,_mm_sub_pd(_mm_mul_pd(rinv02,rinvsq02),krf2));
1744 cutoff_mask = _mm_cmplt_pd(rsq02,rcutoff2);
1748 fscal = _mm_and_pd(fscal,cutoff_mask);
1750 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1752 /* Calculate temporary vectorial force */
1753 tx = _mm_mul_pd(fscal,dx02);
1754 ty = _mm_mul_pd(fscal,dy02);
1755 tz = _mm_mul_pd(fscal,dz02);
1757 /* Update vectorial force */
1758 fix0 = _mm_add_pd(fix0,tx);
1759 fiy0 = _mm_add_pd(fiy0,ty);
1760 fiz0 = _mm_add_pd(fiz0,tz);
1762 fjx2 = _mm_add_pd(fjx2,tx);
1763 fjy2 = _mm_add_pd(fjy2,ty);
1764 fjz2 = _mm_add_pd(fjz2,tz);
1768 /**************************
1769 * CALCULATE INTERACTIONS *
1770 **************************/
1772 if (gmx_mm_any_lt(rsq10,rcutoff2))
1775 /* REACTION-FIELD ELECTROSTATICS */
1776 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
1778 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
1782 fscal = _mm_and_pd(fscal,cutoff_mask);
1784 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1786 /* Calculate temporary vectorial force */
1787 tx = _mm_mul_pd(fscal,dx10);
1788 ty = _mm_mul_pd(fscal,dy10);
1789 tz = _mm_mul_pd(fscal,dz10);
1791 /* Update vectorial force */
1792 fix1 = _mm_add_pd(fix1,tx);
1793 fiy1 = _mm_add_pd(fiy1,ty);
1794 fiz1 = _mm_add_pd(fiz1,tz);
1796 fjx0 = _mm_add_pd(fjx0,tx);
1797 fjy0 = _mm_add_pd(fjy0,ty);
1798 fjz0 = _mm_add_pd(fjz0,tz);
1802 /**************************
1803 * CALCULATE INTERACTIONS *
1804 **************************/
1806 if (gmx_mm_any_lt(rsq11,rcutoff2))
1809 /* REACTION-FIELD ELECTROSTATICS */
1810 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
1812 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
1816 fscal = _mm_and_pd(fscal,cutoff_mask);
1818 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1820 /* Calculate temporary vectorial force */
1821 tx = _mm_mul_pd(fscal,dx11);
1822 ty = _mm_mul_pd(fscal,dy11);
1823 tz = _mm_mul_pd(fscal,dz11);
1825 /* Update vectorial force */
1826 fix1 = _mm_add_pd(fix1,tx);
1827 fiy1 = _mm_add_pd(fiy1,ty);
1828 fiz1 = _mm_add_pd(fiz1,tz);
1830 fjx1 = _mm_add_pd(fjx1,tx);
1831 fjy1 = _mm_add_pd(fjy1,ty);
1832 fjz1 = _mm_add_pd(fjz1,tz);
1836 /**************************
1837 * CALCULATE INTERACTIONS *
1838 **************************/
1840 if (gmx_mm_any_lt(rsq12,rcutoff2))
1843 /* REACTION-FIELD ELECTROSTATICS */
1844 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
1846 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
1850 fscal = _mm_and_pd(fscal,cutoff_mask);
1852 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1854 /* Calculate temporary vectorial force */
1855 tx = _mm_mul_pd(fscal,dx12);
1856 ty = _mm_mul_pd(fscal,dy12);
1857 tz = _mm_mul_pd(fscal,dz12);
1859 /* Update vectorial force */
1860 fix1 = _mm_add_pd(fix1,tx);
1861 fiy1 = _mm_add_pd(fiy1,ty);
1862 fiz1 = _mm_add_pd(fiz1,tz);
1864 fjx2 = _mm_add_pd(fjx2,tx);
1865 fjy2 = _mm_add_pd(fjy2,ty);
1866 fjz2 = _mm_add_pd(fjz2,tz);
1870 /**************************
1871 * CALCULATE INTERACTIONS *
1872 **************************/
1874 if (gmx_mm_any_lt(rsq20,rcutoff2))
1877 /* REACTION-FIELD ELECTROSTATICS */
1878 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
1880 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
1884 fscal = _mm_and_pd(fscal,cutoff_mask);
1886 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1888 /* Calculate temporary vectorial force */
1889 tx = _mm_mul_pd(fscal,dx20);
1890 ty = _mm_mul_pd(fscal,dy20);
1891 tz = _mm_mul_pd(fscal,dz20);
1893 /* Update vectorial force */
1894 fix2 = _mm_add_pd(fix2,tx);
1895 fiy2 = _mm_add_pd(fiy2,ty);
1896 fiz2 = _mm_add_pd(fiz2,tz);
1898 fjx0 = _mm_add_pd(fjx0,tx);
1899 fjy0 = _mm_add_pd(fjy0,ty);
1900 fjz0 = _mm_add_pd(fjz0,tz);
1904 /**************************
1905 * CALCULATE INTERACTIONS *
1906 **************************/
1908 if (gmx_mm_any_lt(rsq21,rcutoff2))
1911 /* REACTION-FIELD ELECTROSTATICS */
1912 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
1914 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
1918 fscal = _mm_and_pd(fscal,cutoff_mask);
1920 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1922 /* Calculate temporary vectorial force */
1923 tx = _mm_mul_pd(fscal,dx21);
1924 ty = _mm_mul_pd(fscal,dy21);
1925 tz = _mm_mul_pd(fscal,dz21);
1927 /* Update vectorial force */
1928 fix2 = _mm_add_pd(fix2,tx);
1929 fiy2 = _mm_add_pd(fiy2,ty);
1930 fiz2 = _mm_add_pd(fiz2,tz);
1932 fjx1 = _mm_add_pd(fjx1,tx);
1933 fjy1 = _mm_add_pd(fjy1,ty);
1934 fjz1 = _mm_add_pd(fjz1,tz);
1938 /**************************
1939 * CALCULATE INTERACTIONS *
1940 **************************/
1942 if (gmx_mm_any_lt(rsq22,rcutoff2))
1945 /* REACTION-FIELD ELECTROSTATICS */
1946 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
1948 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
1952 fscal = _mm_and_pd(fscal,cutoff_mask);
1954 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1956 /* Calculate temporary vectorial force */
1957 tx = _mm_mul_pd(fscal,dx22);
1958 ty = _mm_mul_pd(fscal,dy22);
1959 tz = _mm_mul_pd(fscal,dz22);
1961 /* Update vectorial force */
1962 fix2 = _mm_add_pd(fix2,tx);
1963 fiy2 = _mm_add_pd(fiy2,ty);
1964 fiz2 = _mm_add_pd(fiz2,tz);
1966 fjx2 = _mm_add_pd(fjx2,tx);
1967 fjy2 = _mm_add_pd(fjy2,ty);
1968 fjz2 = _mm_add_pd(fjz2,tz);
1972 gmx_mm_decrement_3rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
1974 /* Inner loop uses 270 flops */
1977 /* End of innermost loop */
1979 gmx_mm_update_iforce_3atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,
1980 f+i_coord_offset,fshift+i_shift_offset);
1982 /* Increment number of inner iterations */
1983 inneriter += j_index_end - j_index_start;
1985 /* Outer loop uses 18 flops */
1988 /* Increment number of outer iterations */
1991 /* Update outer/inner flops */
1993 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W3W3_F,outeriter*18 + inneriter*270);