?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Remove nb-parameters from t_forcerec
[gromacs.git] / src / gromacs / gmxlib / nonbonded / nb_kernel_sse2_double / nb_kernel_ElecRFCut_VdwCSTab_GeomW4W4_sse2_double.c
blob05cbb7b5eb5c7a7ecac678ee422c2f4f4b06c3bc
1 /*
2 * This file is part of the GROMACS molecular simulation package.
3 *
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.
8 *
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.
13 *
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.
18 *
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.
23 *
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.
31 *
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.
34 */
35 /*
36 * Note: this file was generated by the GROMACS sse2_double kernel generator.
37 */
38 #include "gmxpre.h"
39
40 #include "config.h"
41
42 #include <math.h>
43
44 #include "../nb_kernel.h"
45 #include "gromacs/gmxlib/nrnb.h"
46
47 #include "kernelutil_x86_sse2_double.h"
48
49 /*
50 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwCSTab_GeomW4W4_VF_sse2_double
51 * Electrostatics interaction: ReactionField
52 * VdW interaction: CubicSplineTable
53 * Geometry: Water4-Water4
54 * Calculate force/pot: PotentialAndForce
55 */
56 void
57 nb_kernel_ElecRFCut_VdwCSTab_GeomW4W4_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)
65 {
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.
70 */
71 int i_shift_offset,i_coord_offset,outeriter,inneriter;
72 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
73 int jnrA,jnrB;
74 int j_coord_offsetA,j_coord_offsetB;
75 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
76 real rcutoff_scalar;
77 real *shiftvec,*fshift,*x,*f;
78 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
79 int vdwioffset0;
80 __m128d ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
81 int vdwioffset1;
82 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
83 int vdwioffset2;
84 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
85 int vdwioffset3;
86 __m128d ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
87 int vdwjidx0A,vdwjidx0B;
88 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
89 int vdwjidx1A,vdwjidx1B;
90 __m128d jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
91 int vdwjidx2A,vdwjidx2B;
92 __m128d jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
93 int vdwjidx3A,vdwjidx3B;
94 __m128d jx3,jy3,jz3,fjx3,fjy3,fjz3,jq3,isaj3;
95 __m128d dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
96 __m128d dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
97 __m128d dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
98 __m128d dx13,dy13,dz13,rsq13,rinv13,rinvsq13,r13,qq13,c6_13,c12_13;
99 __m128d dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
100 __m128d dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
101 __m128d dx23,dy23,dz23,rsq23,rinv23,rinvsq23,r23,qq23,c6_23,c12_23;
102 __m128d dx31,dy31,dz31,rsq31,rinv31,rinvsq31,r31,qq31,c6_31,c12_31;
103 __m128d dx32,dy32,dz32,rsq32,rinv32,rinvsq32,r32,qq32,c6_32,c12_32;
104 __m128d dx33,dy33,dz33,rsq33,rinv33,rinvsq33,r33,qq33,c6_33,c12_33;
105 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
106 real *charge;
107 int nvdwtype;
108 __m128d rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
109 int *vdwtype;
110 real *vdwparam;
111 __m128d one_sixth = _mm_set1_pd(1.0/6.0);
112 __m128d one_twelfth = _mm_set1_pd(1.0/12.0);
113 __m128i vfitab;
114 __m128i ifour = _mm_set1_epi32(4);
115 __m128d rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF;
116 real *vftab;
117 __m128d dummy_mask,cutoff_mask;
118 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
119 __m128d one = _mm_set1_pd(1.0);
120 __m128d two = _mm_set1_pd(2.0);
121 x = xx[0];
122 f = ff[0];
123
124 nri = nlist->nri;
125 iinr = nlist->iinr;
126 jindex = nlist->jindex;
127 jjnr = nlist->jjnr;
128 shiftidx = nlist->shift;
129 gid = nlist->gid;
130 shiftvec = fr->shift_vec[0];
131 fshift = fr->fshift[0];
132 facel = _mm_set1_pd(fr->ic->epsfac);
133 charge = mdatoms->chargeA;
134 krf = _mm_set1_pd(fr->ic->k_rf);
135 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
136 crf = _mm_set1_pd(fr->ic->c_rf);
137 nvdwtype = fr->ntype;
138 vdwparam = fr->nbfp;
139 vdwtype = mdatoms->typeA;
140
141 vftab = kernel_data->table_vdw->data;
142 vftabscale = _mm_set1_pd(kernel_data->table_vdw->scale);
143
144 /* Setup water-specific parameters */
145 inr = nlist->iinr[0];
146 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
147 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
148 iq3 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+3]));
149 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
150
151 jq1 = _mm_set1_pd(charge[inr+1]);
152 jq2 = _mm_set1_pd(charge[inr+2]);
153 jq3 = _mm_set1_pd(charge[inr+3]);
154 vdwjidx0A = 2*vdwtype[inr+0];
155 c6_00 = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A]);
156 c12_00 = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A+1]);
157 qq11 = _mm_mul_pd(iq1,jq1);
158 qq12 = _mm_mul_pd(iq1,jq2);
159 qq13 = _mm_mul_pd(iq1,jq3);
160 qq21 = _mm_mul_pd(iq2,jq1);
161 qq22 = _mm_mul_pd(iq2,jq2);
162 qq23 = _mm_mul_pd(iq2,jq3);
163 qq31 = _mm_mul_pd(iq3,jq1);
164 qq32 = _mm_mul_pd(iq3,jq2);
165 qq33 = _mm_mul_pd(iq3,jq3);
166
167 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
168 rcutoff_scalar = fr->ic->rcoulomb;
169 rcutoff = _mm_set1_pd(rcutoff_scalar);
170 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
171
172 /* Avoid stupid compiler warnings */
173 jnrA = jnrB = 0;
174 j_coord_offsetA = 0;
175 j_coord_offsetB = 0;
176
177 outeriter = 0;
178 inneriter = 0;
179
180 /* Start outer loop over neighborlists */
181 for(iidx=0; iidx<nri; iidx++)
182 {
183 /* Load shift vector for this list */
184 i_shift_offset = DIM*shiftidx[iidx];
185
186 /* Load limits for loop over neighbors */
187 j_index_start = jindex[iidx];
188 j_index_end = jindex[iidx+1];
189
190 /* Get outer coordinate index */
191 inr = iinr[iidx];
192 i_coord_offset = DIM*inr;
193
194 /* Load i particle coords and add shift vector */
195 gmx_mm_load_shift_and_4rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
196 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
197
198 fix0 = _mm_setzero_pd();
199 fiy0 = _mm_setzero_pd();
200 fiz0 = _mm_setzero_pd();
201 fix1 = _mm_setzero_pd();
202 fiy1 = _mm_setzero_pd();
203 fiz1 = _mm_setzero_pd();
204 fix2 = _mm_setzero_pd();
205 fiy2 = _mm_setzero_pd();
206 fiz2 = _mm_setzero_pd();
207 fix3 = _mm_setzero_pd();
208 fiy3 = _mm_setzero_pd();
209 fiz3 = _mm_setzero_pd();
210
211 /* Reset potential sums */
212 velecsum = _mm_setzero_pd();
213 vvdwsum = _mm_setzero_pd();
214
215 /* Start inner kernel loop */
216 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
217 {
218
219 /* Get j neighbor index, and coordinate index */
220 jnrA = jjnr[jidx];
221 jnrB = jjnr[jidx+1];
222 j_coord_offsetA = DIM*jnrA;
223 j_coord_offsetB = DIM*jnrB;
224
225 /* load j atom coordinates */
226 gmx_mm_load_4rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
227 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
228 &jy2,&jz2,&jx3,&jy3,&jz3);
229
230 /* Calculate displacement vector */
231 dx00 = _mm_sub_pd(ix0,jx0);
232 dy00 = _mm_sub_pd(iy0,jy0);
233 dz00 = _mm_sub_pd(iz0,jz0);
234 dx11 = _mm_sub_pd(ix1,jx1);
235 dy11 = _mm_sub_pd(iy1,jy1);
236 dz11 = _mm_sub_pd(iz1,jz1);
237 dx12 = _mm_sub_pd(ix1,jx2);
238 dy12 = _mm_sub_pd(iy1,jy2);
239 dz12 = _mm_sub_pd(iz1,jz2);
240 dx13 = _mm_sub_pd(ix1,jx3);
241 dy13 = _mm_sub_pd(iy1,jy3);
242 dz13 = _mm_sub_pd(iz1,jz3);
243 dx21 = _mm_sub_pd(ix2,jx1);
244 dy21 = _mm_sub_pd(iy2,jy1);
245 dz21 = _mm_sub_pd(iz2,jz1);
246 dx22 = _mm_sub_pd(ix2,jx2);
247 dy22 = _mm_sub_pd(iy2,jy2);
248 dz22 = _mm_sub_pd(iz2,jz2);
249 dx23 = _mm_sub_pd(ix2,jx3);
250 dy23 = _mm_sub_pd(iy2,jy3);
251 dz23 = _mm_sub_pd(iz2,jz3);
252 dx31 = _mm_sub_pd(ix3,jx1);
253 dy31 = _mm_sub_pd(iy3,jy1);
254 dz31 = _mm_sub_pd(iz3,jz1);
255 dx32 = _mm_sub_pd(ix3,jx2);
256 dy32 = _mm_sub_pd(iy3,jy2);
257 dz32 = _mm_sub_pd(iz3,jz2);
258 dx33 = _mm_sub_pd(ix3,jx3);
259 dy33 = _mm_sub_pd(iy3,jy3);
260 dz33 = _mm_sub_pd(iz3,jz3);
261
262 /* Calculate squared distance and things based on it */
263 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
264 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
265 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
266 rsq13 = gmx_mm_calc_rsq_pd(dx13,dy13,dz13);
267 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
268 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
269 rsq23 = gmx_mm_calc_rsq_pd(dx23,dy23,dz23);
270 rsq31 = gmx_mm_calc_rsq_pd(dx31,dy31,dz31);
271 rsq32 = gmx_mm_calc_rsq_pd(dx32,dy32,dz32);
272 rsq33 = gmx_mm_calc_rsq_pd(dx33,dy33,dz33);
273
274 rinv00 = sse2_invsqrt_d(rsq00);
275 rinv11 = sse2_invsqrt_d(rsq11);
276 rinv12 = sse2_invsqrt_d(rsq12);
277 rinv13 = sse2_invsqrt_d(rsq13);
278 rinv21 = sse2_invsqrt_d(rsq21);
279 rinv22 = sse2_invsqrt_d(rsq22);
280 rinv23 = sse2_invsqrt_d(rsq23);
281 rinv31 = sse2_invsqrt_d(rsq31);
282 rinv32 = sse2_invsqrt_d(rsq32);
283 rinv33 = sse2_invsqrt_d(rsq33);
284
285 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
286 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
287 rinvsq13 = _mm_mul_pd(rinv13,rinv13);
288 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
289 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
290 rinvsq23 = _mm_mul_pd(rinv23,rinv23);
291 rinvsq31 = _mm_mul_pd(rinv31,rinv31);
292 rinvsq32 = _mm_mul_pd(rinv32,rinv32);
293 rinvsq33 = _mm_mul_pd(rinv33,rinv33);
294
295 fjx0 = _mm_setzero_pd();
296 fjy0 = _mm_setzero_pd();
297 fjz0 = _mm_setzero_pd();
298 fjx1 = _mm_setzero_pd();
299 fjy1 = _mm_setzero_pd();
300 fjz1 = _mm_setzero_pd();
301 fjx2 = _mm_setzero_pd();
302 fjy2 = _mm_setzero_pd();
303 fjz2 = _mm_setzero_pd();
304 fjx3 = _mm_setzero_pd();
305 fjy3 = _mm_setzero_pd();
306 fjz3 = _mm_setzero_pd();
307
308 /**************************
309 * CALCULATE INTERACTIONS *
310 **************************/
311
312 r00 = _mm_mul_pd(rsq00,rinv00);
313
314 /* Calculate table index by multiplying r with table scale and truncate to integer */
315 rt = _mm_mul_pd(r00,vftabscale);
316 vfitab = _mm_cvttpd_epi32(rt);
317 vfeps = _mm_sub_pd(rt,_mm_cvtepi32_pd(vfitab));
318 vfitab = _mm_slli_epi32(vfitab,3);
319
320 /* CUBIC SPLINE TABLE DISPERSION */
321 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
322 F = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) );
323 GMX_MM_TRANSPOSE2_PD(Y,F);
324 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
325 H = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) +2);
326 GMX_MM_TRANSPOSE2_PD(G,H);
327 Heps = _mm_mul_pd(vfeps,H);
328 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
329 VV = _mm_add_pd(Y,_mm_mul_pd(vfeps,Fp));
330 vvdw6 = _mm_mul_pd(c6_00,VV);
331 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
332 fvdw6 = _mm_mul_pd(c6_00,FF);
333
334 /* CUBIC SPLINE TABLE REPULSION */
335 vfitab = _mm_add_epi32(vfitab,ifour);
336 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
337 F = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) );
338 GMX_MM_TRANSPOSE2_PD(Y,F);
339 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
340 H = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) +2);
341 GMX_MM_TRANSPOSE2_PD(G,H);
342 Heps = _mm_mul_pd(vfeps,H);
343 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
344 VV = _mm_add_pd(Y,_mm_mul_pd(vfeps,Fp));
345 vvdw12 = _mm_mul_pd(c12_00,VV);
346 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
347 fvdw12 = _mm_mul_pd(c12_00,FF);
348 vvdw = _mm_add_pd(vvdw12,vvdw6);
349 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
350
351 /* Update potential sum for this i atom from the interaction with this j atom. */
352 vvdwsum = _mm_add_pd(vvdwsum,vvdw);
353
354 fscal = fvdw;
355
356 /* Calculate temporary vectorial force */
357 tx = _mm_mul_pd(fscal,dx00);
358 ty = _mm_mul_pd(fscal,dy00);
359 tz = _mm_mul_pd(fscal,dz00);
360
361 /* Update vectorial force */
362 fix0 = _mm_add_pd(fix0,tx);
363 fiy0 = _mm_add_pd(fiy0,ty);
364 fiz0 = _mm_add_pd(fiz0,tz);
365
366 fjx0 = _mm_add_pd(fjx0,tx);
367 fjy0 = _mm_add_pd(fjy0,ty);
368 fjz0 = _mm_add_pd(fjz0,tz);
369
370 /**************************
371 * CALCULATE INTERACTIONS *
372 **************************/
373
374 if (gmx_mm_any_lt(rsq11,rcutoff2))
375 {
376
377 /* REACTION-FIELD ELECTROSTATICS */
378 velec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_add_pd(rinv11,_mm_mul_pd(krf,rsq11)),crf));
379 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
380
381 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
382
383 /* Update potential sum for this i atom from the interaction with this j atom. */
384 velec = _mm_and_pd(velec,cutoff_mask);
385 velecsum = _mm_add_pd(velecsum,velec);
386
387 fscal = felec;
388
389 fscal = _mm_and_pd(fscal,cutoff_mask);
390
391 /* Calculate temporary vectorial force */
392 tx = _mm_mul_pd(fscal,dx11);
393 ty = _mm_mul_pd(fscal,dy11);
394 tz = _mm_mul_pd(fscal,dz11);
395
396 /* Update vectorial force */
397 fix1 = _mm_add_pd(fix1,tx);
398 fiy1 = _mm_add_pd(fiy1,ty);
399 fiz1 = _mm_add_pd(fiz1,tz);
400
401 fjx1 = _mm_add_pd(fjx1,tx);
402 fjy1 = _mm_add_pd(fjy1,ty);
403 fjz1 = _mm_add_pd(fjz1,tz);
404
405 }
406
407 /**************************
408 * CALCULATE INTERACTIONS *
409 **************************/
410
411 if (gmx_mm_any_lt(rsq12,rcutoff2))
412 {
413
414 /* REACTION-FIELD ELECTROSTATICS */
415 velec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_add_pd(rinv12,_mm_mul_pd(krf,rsq12)),crf));
416 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
417
418 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
419
420 /* Update potential sum for this i atom from the interaction with this j atom. */
421 velec = _mm_and_pd(velec,cutoff_mask);
422 velecsum = _mm_add_pd(velecsum,velec);
423
424 fscal = felec;
425
426 fscal = _mm_and_pd(fscal,cutoff_mask);
427
428 /* Calculate temporary vectorial force */
429 tx = _mm_mul_pd(fscal,dx12);
430 ty = _mm_mul_pd(fscal,dy12);
431 tz = _mm_mul_pd(fscal,dz12);
432
433 /* Update vectorial force */
434 fix1 = _mm_add_pd(fix1,tx);
435 fiy1 = _mm_add_pd(fiy1,ty);
436 fiz1 = _mm_add_pd(fiz1,tz);
437
438 fjx2 = _mm_add_pd(fjx2,tx);
439 fjy2 = _mm_add_pd(fjy2,ty);
440 fjz2 = _mm_add_pd(fjz2,tz);
441
442 }
443
444 /**************************
445 * CALCULATE INTERACTIONS *
446 **************************/
447
448 if (gmx_mm_any_lt(rsq13,rcutoff2))
449 {
450
451 /* REACTION-FIELD ELECTROSTATICS */
452 velec = _mm_mul_pd(qq13,_mm_sub_pd(_mm_add_pd(rinv13,_mm_mul_pd(krf,rsq13)),crf));
453 felec = _mm_mul_pd(qq13,_mm_sub_pd(_mm_mul_pd(rinv13,rinvsq13),krf2));
454
455 cutoff_mask = _mm_cmplt_pd(rsq13,rcutoff2);
456
457 /* Update potential sum for this i atom from the interaction with this j atom. */
458 velec = _mm_and_pd(velec,cutoff_mask);
459 velecsum = _mm_add_pd(velecsum,velec);
460
461 fscal = felec;
462
463 fscal = _mm_and_pd(fscal,cutoff_mask);
464
465 /* Calculate temporary vectorial force */
466 tx = _mm_mul_pd(fscal,dx13);
467 ty = _mm_mul_pd(fscal,dy13);
468 tz = _mm_mul_pd(fscal,dz13);
469
470 /* Update vectorial force */
471 fix1 = _mm_add_pd(fix1,tx);
472 fiy1 = _mm_add_pd(fiy1,ty);
473 fiz1 = _mm_add_pd(fiz1,tz);
474
475 fjx3 = _mm_add_pd(fjx3,tx);
476 fjy3 = _mm_add_pd(fjy3,ty);
477 fjz3 = _mm_add_pd(fjz3,tz);
478
479 }
480
481 /**************************
482 * CALCULATE INTERACTIONS *
483 **************************/
484
485 if (gmx_mm_any_lt(rsq21,rcutoff2))
486 {
487
488 /* REACTION-FIELD ELECTROSTATICS */
489 velec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_add_pd(rinv21,_mm_mul_pd(krf,rsq21)),crf));
490 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
491
492 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
493
494 /* Update potential sum for this i atom from the interaction with this j atom. */
495 velec = _mm_and_pd(velec,cutoff_mask);
496 velecsum = _mm_add_pd(velecsum,velec);
497
498 fscal = felec;
499
500 fscal = _mm_and_pd(fscal,cutoff_mask);
501
502 /* Calculate temporary vectorial force */
503 tx = _mm_mul_pd(fscal,dx21);
504 ty = _mm_mul_pd(fscal,dy21);
505 tz = _mm_mul_pd(fscal,dz21);
506
507 /* Update vectorial force */
508 fix2 = _mm_add_pd(fix2,tx);
509 fiy2 = _mm_add_pd(fiy2,ty);
510 fiz2 = _mm_add_pd(fiz2,tz);
511
512 fjx1 = _mm_add_pd(fjx1,tx);
513 fjy1 = _mm_add_pd(fjy1,ty);
514 fjz1 = _mm_add_pd(fjz1,tz);
515
516 }
517
518 /**************************
519 * CALCULATE INTERACTIONS *
520 **************************/
521
522 if (gmx_mm_any_lt(rsq22,rcutoff2))
523 {
524
525 /* REACTION-FIELD ELECTROSTATICS */
526 velec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_add_pd(rinv22,_mm_mul_pd(krf,rsq22)),crf));
527 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
528
529 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
530
531 /* Update potential sum for this i atom from the interaction with this j atom. */
532 velec = _mm_and_pd(velec,cutoff_mask);
533 velecsum = _mm_add_pd(velecsum,velec);
534
535 fscal = felec;
536
537 fscal = _mm_and_pd(fscal,cutoff_mask);
538
539 /* Calculate temporary vectorial force */
540 tx = _mm_mul_pd(fscal,dx22);
541 ty = _mm_mul_pd(fscal,dy22);
542 tz = _mm_mul_pd(fscal,dz22);
543
544 /* Update vectorial force */
545 fix2 = _mm_add_pd(fix2,tx);
546 fiy2 = _mm_add_pd(fiy2,ty);
547 fiz2 = _mm_add_pd(fiz2,tz);
548
549 fjx2 = _mm_add_pd(fjx2,tx);
550 fjy2 = _mm_add_pd(fjy2,ty);
551 fjz2 = _mm_add_pd(fjz2,tz);
552
553 }
554
555 /**************************
556 * CALCULATE INTERACTIONS *
557 **************************/
558
559 if (gmx_mm_any_lt(rsq23,rcutoff2))
560 {
561
562 /* REACTION-FIELD ELECTROSTATICS */
563 velec = _mm_mul_pd(qq23,_mm_sub_pd(_mm_add_pd(rinv23,_mm_mul_pd(krf,rsq23)),crf));
564 felec = _mm_mul_pd(qq23,_mm_sub_pd(_mm_mul_pd(rinv23,rinvsq23),krf2));
565
566 cutoff_mask = _mm_cmplt_pd(rsq23,rcutoff2);
567
568 /* Update potential sum for this i atom from the interaction with this j atom. */
569 velec = _mm_and_pd(velec,cutoff_mask);
570 velecsum = _mm_add_pd(velecsum,velec);
571
572 fscal = felec;
573
574 fscal = _mm_and_pd(fscal,cutoff_mask);
575
576 /* Calculate temporary vectorial force */
577 tx = _mm_mul_pd(fscal,dx23);
578 ty = _mm_mul_pd(fscal,dy23);
579 tz = _mm_mul_pd(fscal,dz23);
580
581 /* Update vectorial force */
582 fix2 = _mm_add_pd(fix2,tx);
583 fiy2 = _mm_add_pd(fiy2,ty);
584 fiz2 = _mm_add_pd(fiz2,tz);
585
586 fjx3 = _mm_add_pd(fjx3,tx);
587 fjy3 = _mm_add_pd(fjy3,ty);
588 fjz3 = _mm_add_pd(fjz3,tz);
589
590 }
591
592 /**************************
593 * CALCULATE INTERACTIONS *
594 **************************/
595
596 if (gmx_mm_any_lt(rsq31,rcutoff2))
597 {
598
599 /* REACTION-FIELD ELECTROSTATICS */
600 velec = _mm_mul_pd(qq31,_mm_sub_pd(_mm_add_pd(rinv31,_mm_mul_pd(krf,rsq31)),crf));
601 felec = _mm_mul_pd(qq31,_mm_sub_pd(_mm_mul_pd(rinv31,rinvsq31),krf2));
602
603 cutoff_mask = _mm_cmplt_pd(rsq31,rcutoff2);
604
605 /* Update potential sum for this i atom from the interaction with this j atom. */
606 velec = _mm_and_pd(velec,cutoff_mask);
607 velecsum = _mm_add_pd(velecsum,velec);
608
609 fscal = felec;
610
611 fscal = _mm_and_pd(fscal,cutoff_mask);
612
613 /* Calculate temporary vectorial force */
614 tx = _mm_mul_pd(fscal,dx31);
615 ty = _mm_mul_pd(fscal,dy31);
616 tz = _mm_mul_pd(fscal,dz31);
617
618 /* Update vectorial force */
619 fix3 = _mm_add_pd(fix3,tx);
620 fiy3 = _mm_add_pd(fiy3,ty);
621 fiz3 = _mm_add_pd(fiz3,tz);
622
623 fjx1 = _mm_add_pd(fjx1,tx);
624 fjy1 = _mm_add_pd(fjy1,ty);
625 fjz1 = _mm_add_pd(fjz1,tz);
626
627 }
628
629 /**************************
630 * CALCULATE INTERACTIONS *
631 **************************/
632
633 if (gmx_mm_any_lt(rsq32,rcutoff2))
634 {
635
636 /* REACTION-FIELD ELECTROSTATICS */
637 velec = _mm_mul_pd(qq32,_mm_sub_pd(_mm_add_pd(rinv32,_mm_mul_pd(krf,rsq32)),crf));
638 felec = _mm_mul_pd(qq32,_mm_sub_pd(_mm_mul_pd(rinv32,rinvsq32),krf2));
639
640 cutoff_mask = _mm_cmplt_pd(rsq32,rcutoff2);
641
642 /* Update potential sum for this i atom from the interaction with this j atom. */
643 velec = _mm_and_pd(velec,cutoff_mask);
644 velecsum = _mm_add_pd(velecsum,velec);
645
646 fscal = felec;
647
648 fscal = _mm_and_pd(fscal,cutoff_mask);
649
650 /* Calculate temporary vectorial force */
651 tx = _mm_mul_pd(fscal,dx32);
652 ty = _mm_mul_pd(fscal,dy32);
653 tz = _mm_mul_pd(fscal,dz32);
654
655 /* Update vectorial force */
656 fix3 = _mm_add_pd(fix3,tx);
657 fiy3 = _mm_add_pd(fiy3,ty);
658 fiz3 = _mm_add_pd(fiz3,tz);
659
660 fjx2 = _mm_add_pd(fjx2,tx);
661 fjy2 = _mm_add_pd(fjy2,ty);
662 fjz2 = _mm_add_pd(fjz2,tz);
663
664 }
665
666 /**************************
667 * CALCULATE INTERACTIONS *
668 **************************/
669
670 if (gmx_mm_any_lt(rsq33,rcutoff2))
671 {
672
673 /* REACTION-FIELD ELECTROSTATICS */
674 velec = _mm_mul_pd(qq33,_mm_sub_pd(_mm_add_pd(rinv33,_mm_mul_pd(krf,rsq33)),crf));
675 felec = _mm_mul_pd(qq33,_mm_sub_pd(_mm_mul_pd(rinv33,rinvsq33),krf2));
676
677 cutoff_mask = _mm_cmplt_pd(rsq33,rcutoff2);
678
679 /* Update potential sum for this i atom from the interaction with this j atom. */
680 velec = _mm_and_pd(velec,cutoff_mask);
681 velecsum = _mm_add_pd(velecsum,velec);
682
683 fscal = felec;
684
685 fscal = _mm_and_pd(fscal,cutoff_mask);
686
687 /* Calculate temporary vectorial force */
688 tx = _mm_mul_pd(fscal,dx33);
689 ty = _mm_mul_pd(fscal,dy33);
690 tz = _mm_mul_pd(fscal,dz33);
691
692 /* Update vectorial force */
693 fix3 = _mm_add_pd(fix3,tx);
694 fiy3 = _mm_add_pd(fiy3,ty);
695 fiz3 = _mm_add_pd(fiz3,tz);
696
697 fjx3 = _mm_add_pd(fjx3,tx);
698 fjy3 = _mm_add_pd(fjy3,ty);
699 fjz3 = _mm_add_pd(fjz3,tz);
700
701 }
702
703 gmx_mm_decrement_4rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
704
705 /* Inner loop uses 383 flops */
706 }
707
708 if(jidx<j_index_end)
709 {
710
711 jnrA = jjnr[jidx];
712 j_coord_offsetA = DIM*jnrA;
713
714 /* load j atom coordinates */
715 gmx_mm_load_4rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
716 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
717 &jy2,&jz2,&jx3,&jy3,&jz3);
718
719 /* Calculate displacement vector */
720 dx00 = _mm_sub_pd(ix0,jx0);
721 dy00 = _mm_sub_pd(iy0,jy0);
722 dz00 = _mm_sub_pd(iz0,jz0);
723 dx11 = _mm_sub_pd(ix1,jx1);
724 dy11 = _mm_sub_pd(iy1,jy1);
725 dz11 = _mm_sub_pd(iz1,jz1);
726 dx12 = _mm_sub_pd(ix1,jx2);
727 dy12 = _mm_sub_pd(iy1,jy2);
728 dz12 = _mm_sub_pd(iz1,jz2);
729 dx13 = _mm_sub_pd(ix1,jx3);
730 dy13 = _mm_sub_pd(iy1,jy3);
731 dz13 = _mm_sub_pd(iz1,jz3);
732 dx21 = _mm_sub_pd(ix2,jx1);
733 dy21 = _mm_sub_pd(iy2,jy1);
734 dz21 = _mm_sub_pd(iz2,jz1);
735 dx22 = _mm_sub_pd(ix2,jx2);
736 dy22 = _mm_sub_pd(iy2,jy2);
737 dz22 = _mm_sub_pd(iz2,jz2);
738 dx23 = _mm_sub_pd(ix2,jx3);
739 dy23 = _mm_sub_pd(iy2,jy3);
740 dz23 = _mm_sub_pd(iz2,jz3);
741 dx31 = _mm_sub_pd(ix3,jx1);
742 dy31 = _mm_sub_pd(iy3,jy1);
743 dz31 = _mm_sub_pd(iz3,jz1);
744 dx32 = _mm_sub_pd(ix3,jx2);
745 dy32 = _mm_sub_pd(iy3,jy2);
746 dz32 = _mm_sub_pd(iz3,jz2);
747 dx33 = _mm_sub_pd(ix3,jx3);
748 dy33 = _mm_sub_pd(iy3,jy3);
749 dz33 = _mm_sub_pd(iz3,jz3);
750
751 /* Calculate squared distance and things based on it */
752 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
753 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
754 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
755 rsq13 = gmx_mm_calc_rsq_pd(dx13,dy13,dz13);
756 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
757 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
758 rsq23 = gmx_mm_calc_rsq_pd(dx23,dy23,dz23);
759 rsq31 = gmx_mm_calc_rsq_pd(dx31,dy31,dz31);
760 rsq32 = gmx_mm_calc_rsq_pd(dx32,dy32,dz32);
761 rsq33 = gmx_mm_calc_rsq_pd(dx33,dy33,dz33);
762
763 rinv00 = sse2_invsqrt_d(rsq00);
764 rinv11 = sse2_invsqrt_d(rsq11);
765 rinv12 = sse2_invsqrt_d(rsq12);
766 rinv13 = sse2_invsqrt_d(rsq13);
767 rinv21 = sse2_invsqrt_d(rsq21);
768 rinv22 = sse2_invsqrt_d(rsq22);
769 rinv23 = sse2_invsqrt_d(rsq23);
770 rinv31 = sse2_invsqrt_d(rsq31);
771 rinv32 = sse2_invsqrt_d(rsq32);
772 rinv33 = sse2_invsqrt_d(rsq33);
773
774 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
775 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
776 rinvsq13 = _mm_mul_pd(rinv13,rinv13);
777 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
778 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
779 rinvsq23 = _mm_mul_pd(rinv23,rinv23);
780 rinvsq31 = _mm_mul_pd(rinv31,rinv31);
781 rinvsq32 = _mm_mul_pd(rinv32,rinv32);
782 rinvsq33 = _mm_mul_pd(rinv33,rinv33);
783
784 fjx0 = _mm_setzero_pd();
785 fjy0 = _mm_setzero_pd();
786 fjz0 = _mm_setzero_pd();
787 fjx1 = _mm_setzero_pd();
788 fjy1 = _mm_setzero_pd();
789 fjz1 = _mm_setzero_pd();
790 fjx2 = _mm_setzero_pd();
791 fjy2 = _mm_setzero_pd();
792 fjz2 = _mm_setzero_pd();
793 fjx3 = _mm_setzero_pd();
794 fjy3 = _mm_setzero_pd();
795 fjz3 = _mm_setzero_pd();
796
797 /**************************
798 * CALCULATE INTERACTIONS *
799 **************************/
800
801 r00 = _mm_mul_pd(rsq00,rinv00);
802
803 /* Calculate table index by multiplying r with table scale and truncate to integer */
804 rt = _mm_mul_pd(r00,vftabscale);
805 vfitab = _mm_cvttpd_epi32(rt);
806 vfeps = _mm_sub_pd(rt,_mm_cvtepi32_pd(vfitab));
807 vfitab = _mm_slli_epi32(vfitab,3);
808
809 /* CUBIC SPLINE TABLE DISPERSION */
810 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
811 F = _mm_setzero_pd();
812 GMX_MM_TRANSPOSE2_PD(Y,F);
813 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
814 H = _mm_setzero_pd();
815 GMX_MM_TRANSPOSE2_PD(G,H);
816 Heps = _mm_mul_pd(vfeps,H);
817 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
818 VV = _mm_add_pd(Y,_mm_mul_pd(vfeps,Fp));
819 vvdw6 = _mm_mul_pd(c6_00,VV);
820 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
821 fvdw6 = _mm_mul_pd(c6_00,FF);
822
823 /* CUBIC SPLINE TABLE REPULSION */
824 vfitab = _mm_add_epi32(vfitab,ifour);
825 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
826 F = _mm_setzero_pd();
827 GMX_MM_TRANSPOSE2_PD(Y,F);
828 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
829 H = _mm_setzero_pd();
830 GMX_MM_TRANSPOSE2_PD(G,H);
831 Heps = _mm_mul_pd(vfeps,H);
832 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
833 VV = _mm_add_pd(Y,_mm_mul_pd(vfeps,Fp));
834 vvdw12 = _mm_mul_pd(c12_00,VV);
835 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
836 fvdw12 = _mm_mul_pd(c12_00,FF);
837 vvdw = _mm_add_pd(vvdw12,vvdw6);
838 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
839
840 /* Update potential sum for this i atom from the interaction with this j atom. */
841 vvdw = _mm_unpacklo_pd(vvdw,_mm_setzero_pd());
842 vvdwsum = _mm_add_pd(vvdwsum,vvdw);
843
844 fscal = fvdw;
845
846 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
847
848 /* Calculate temporary vectorial force */
849 tx = _mm_mul_pd(fscal,dx00);
850 ty = _mm_mul_pd(fscal,dy00);
851 tz = _mm_mul_pd(fscal,dz00);
852
853 /* Update vectorial force */
854 fix0 = _mm_add_pd(fix0,tx);
855 fiy0 = _mm_add_pd(fiy0,ty);
856 fiz0 = _mm_add_pd(fiz0,tz);
857
858 fjx0 = _mm_add_pd(fjx0,tx);
859 fjy0 = _mm_add_pd(fjy0,ty);
860 fjz0 = _mm_add_pd(fjz0,tz);
861
862 /**************************
863 * CALCULATE INTERACTIONS *
864 **************************/
865
866 if (gmx_mm_any_lt(rsq11,rcutoff2))
867 {
868
869 /* REACTION-FIELD ELECTROSTATICS */
870 velec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_add_pd(rinv11,_mm_mul_pd(krf,rsq11)),crf));
871 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
872
873 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
874
875 /* Update potential sum for this i atom from the interaction with this j atom. */
876 velec = _mm_and_pd(velec,cutoff_mask);
877 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
878 velecsum = _mm_add_pd(velecsum,velec);
879
880 fscal = felec;
881
882 fscal = _mm_and_pd(fscal,cutoff_mask);
883
884 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
885
886 /* Calculate temporary vectorial force */
887 tx = _mm_mul_pd(fscal,dx11);
888 ty = _mm_mul_pd(fscal,dy11);
889 tz = _mm_mul_pd(fscal,dz11);
890
891 /* Update vectorial force */
892 fix1 = _mm_add_pd(fix1,tx);
893 fiy1 = _mm_add_pd(fiy1,ty);
894 fiz1 = _mm_add_pd(fiz1,tz);
895
896 fjx1 = _mm_add_pd(fjx1,tx);
897 fjy1 = _mm_add_pd(fjy1,ty);
898 fjz1 = _mm_add_pd(fjz1,tz);
899
900 }
901
902 /**************************
903 * CALCULATE INTERACTIONS *
904 **************************/
905
906 if (gmx_mm_any_lt(rsq12,rcutoff2))
907 {
908
909 /* REACTION-FIELD ELECTROSTATICS */
910 velec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_add_pd(rinv12,_mm_mul_pd(krf,rsq12)),crf));
911 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
912
913 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
914
915 /* Update potential sum for this i atom from the interaction with this j atom. */
916 velec = _mm_and_pd(velec,cutoff_mask);
917 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
918 velecsum = _mm_add_pd(velecsum,velec);
919
920 fscal = felec;
921
922 fscal = _mm_and_pd(fscal,cutoff_mask);
923
924 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
925
926 /* Calculate temporary vectorial force */
927 tx = _mm_mul_pd(fscal,dx12);
928 ty = _mm_mul_pd(fscal,dy12);
929 tz = _mm_mul_pd(fscal,dz12);
930
931 /* Update vectorial force */
932 fix1 = _mm_add_pd(fix1,tx);
933 fiy1 = _mm_add_pd(fiy1,ty);
934 fiz1 = _mm_add_pd(fiz1,tz);
935
936 fjx2 = _mm_add_pd(fjx2,tx);
937 fjy2 = _mm_add_pd(fjy2,ty);
938 fjz2 = _mm_add_pd(fjz2,tz);
939
940 }
941
942 /**************************
943 * CALCULATE INTERACTIONS *
944 **************************/
945
946 if (gmx_mm_any_lt(rsq13,rcutoff2))
947 {
948
949 /* REACTION-FIELD ELECTROSTATICS */
950 velec = _mm_mul_pd(qq13,_mm_sub_pd(_mm_add_pd(rinv13,_mm_mul_pd(krf,rsq13)),crf));
951 felec = _mm_mul_pd(qq13,_mm_sub_pd(_mm_mul_pd(rinv13,rinvsq13),krf2));
952
953 cutoff_mask = _mm_cmplt_pd(rsq13,rcutoff2);
954
955 /* Update potential sum for this i atom from the interaction with this j atom. */
956 velec = _mm_and_pd(velec,cutoff_mask);
957 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
958 velecsum = _mm_add_pd(velecsum,velec);
959
960 fscal = felec;
961
962 fscal = _mm_and_pd(fscal,cutoff_mask);
963
964 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
965
966 /* Calculate temporary vectorial force */
967 tx = _mm_mul_pd(fscal,dx13);
968 ty = _mm_mul_pd(fscal,dy13);
969 tz = _mm_mul_pd(fscal,dz13);
970
971 /* Update vectorial force */
972 fix1 = _mm_add_pd(fix1,tx);
973 fiy1 = _mm_add_pd(fiy1,ty);
974 fiz1 = _mm_add_pd(fiz1,tz);
975
976 fjx3 = _mm_add_pd(fjx3,tx);
977 fjy3 = _mm_add_pd(fjy3,ty);
978 fjz3 = _mm_add_pd(fjz3,tz);
979
980 }
981
982 /**************************
983 * CALCULATE INTERACTIONS *
984 **************************/
985
986 if (gmx_mm_any_lt(rsq21,rcutoff2))
987 {
988
989 /* REACTION-FIELD ELECTROSTATICS */
990 velec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_add_pd(rinv21,_mm_mul_pd(krf,rsq21)),crf));
991 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
992
993 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
994
995 /* Update potential sum for this i atom from the interaction with this j atom. */
996 velec = _mm_and_pd(velec,cutoff_mask);
997 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
998 velecsum = _mm_add_pd(velecsum,velec);
999
1000 fscal = felec;
1001
1002 fscal = _mm_and_pd(fscal,cutoff_mask);
1003
1004 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1005
1006 /* Calculate temporary vectorial force */
1007 tx = _mm_mul_pd(fscal,dx21);
1008 ty = _mm_mul_pd(fscal,dy21);
1009 tz = _mm_mul_pd(fscal,dz21);
1010
1011 /* Update vectorial force */
1012 fix2 = _mm_add_pd(fix2,tx);
1013 fiy2 = _mm_add_pd(fiy2,ty);
1014 fiz2 = _mm_add_pd(fiz2,tz);
1015
1016 fjx1 = _mm_add_pd(fjx1,tx);
1017 fjy1 = _mm_add_pd(fjy1,ty);
1018 fjz1 = _mm_add_pd(fjz1,tz);
1019
1020 }
1021
1022 /**************************
1023 * CALCULATE INTERACTIONS *
1024 **************************/
1025
1026 if (gmx_mm_any_lt(rsq22,rcutoff2))
1027 {
1028
1029 /* REACTION-FIELD ELECTROSTATICS */
1030 velec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_add_pd(rinv22,_mm_mul_pd(krf,rsq22)),crf));
1031 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
1032
1033 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
1034
1035 /* Update potential sum for this i atom from the interaction with this j atom. */
1036 velec = _mm_and_pd(velec,cutoff_mask);
1037 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
1038 velecsum = _mm_add_pd(velecsum,velec);
1039
1040 fscal = felec;
1041
1042 fscal = _mm_and_pd(fscal,cutoff_mask);
1043
1044 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1045
1046 /* Calculate temporary vectorial force */
1047 tx = _mm_mul_pd(fscal,dx22);
1048 ty = _mm_mul_pd(fscal,dy22);
1049 tz = _mm_mul_pd(fscal,dz22);
1050
1051 /* Update vectorial force */
1052 fix2 = _mm_add_pd(fix2,tx);
1053 fiy2 = _mm_add_pd(fiy2,ty);
1054 fiz2 = _mm_add_pd(fiz2,tz);
1055
1056 fjx2 = _mm_add_pd(fjx2,tx);
1057 fjy2 = _mm_add_pd(fjy2,ty);
1058 fjz2 = _mm_add_pd(fjz2,tz);
1059
1060 }
1061
1062 /**************************
1063 * CALCULATE INTERACTIONS *
1064 **************************/
1065
1066 if (gmx_mm_any_lt(rsq23,rcutoff2))
1067 {
1068
1069 /* REACTION-FIELD ELECTROSTATICS */
1070 velec = _mm_mul_pd(qq23,_mm_sub_pd(_mm_add_pd(rinv23,_mm_mul_pd(krf,rsq23)),crf));
1071 felec = _mm_mul_pd(qq23,_mm_sub_pd(_mm_mul_pd(rinv23,rinvsq23),krf2));
1072
1073 cutoff_mask = _mm_cmplt_pd(rsq23,rcutoff2);
1074
1075 /* Update potential sum for this i atom from the interaction with this j atom. */
1076 velec = _mm_and_pd(velec,cutoff_mask);
1077 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
1078 velecsum = _mm_add_pd(velecsum,velec);
1079
1080 fscal = felec;
1081
1082 fscal = _mm_and_pd(fscal,cutoff_mask);
1083
1084 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1085
1086 /* Calculate temporary vectorial force */
1087 tx = _mm_mul_pd(fscal,dx23);
1088 ty = _mm_mul_pd(fscal,dy23);
1089 tz = _mm_mul_pd(fscal,dz23);
1090
1091 /* Update vectorial force */
1092 fix2 = _mm_add_pd(fix2,tx);
1093 fiy2 = _mm_add_pd(fiy2,ty);
1094 fiz2 = _mm_add_pd(fiz2,tz);
1095
1096 fjx3 = _mm_add_pd(fjx3,tx);
1097 fjy3 = _mm_add_pd(fjy3,ty);
1098 fjz3 = _mm_add_pd(fjz3,tz);
1099
1100 }
1101
1102 /**************************
1103 * CALCULATE INTERACTIONS *
1104 **************************/
1105
1106 if (gmx_mm_any_lt(rsq31,rcutoff2))
1107 {
1108
1109 /* REACTION-FIELD ELECTROSTATICS */
1110 velec = _mm_mul_pd(qq31,_mm_sub_pd(_mm_add_pd(rinv31,_mm_mul_pd(krf,rsq31)),crf));
1111 felec = _mm_mul_pd(qq31,_mm_sub_pd(_mm_mul_pd(rinv31,rinvsq31),krf2));
1112
1113 cutoff_mask = _mm_cmplt_pd(rsq31,rcutoff2);
1114
1115 /* Update potential sum for this i atom from the interaction with this j atom. */
1116 velec = _mm_and_pd(velec,cutoff_mask);
1117 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
1118 velecsum = _mm_add_pd(velecsum,velec);
1119
1120 fscal = felec;
1121
1122 fscal = _mm_and_pd(fscal,cutoff_mask);
1123
1124 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1125
1126 /* Calculate temporary vectorial force */
1127 tx = _mm_mul_pd(fscal,dx31);
1128 ty = _mm_mul_pd(fscal,dy31);
1129 tz = _mm_mul_pd(fscal,dz31);
1130
1131 /* Update vectorial force */
1132 fix3 = _mm_add_pd(fix3,tx);
1133 fiy3 = _mm_add_pd(fiy3,ty);
1134 fiz3 = _mm_add_pd(fiz3,tz);
1135
1136 fjx1 = _mm_add_pd(fjx1,tx);
1137 fjy1 = _mm_add_pd(fjy1,ty);
1138 fjz1 = _mm_add_pd(fjz1,tz);
1139
1140 }
1141
1142 /**************************
1143 * CALCULATE INTERACTIONS *
1144 **************************/
1145
1146 if (gmx_mm_any_lt(rsq32,rcutoff2))
1147 {
1148
1149 /* REACTION-FIELD ELECTROSTATICS */
1150 velec = _mm_mul_pd(qq32,_mm_sub_pd(_mm_add_pd(rinv32,_mm_mul_pd(krf,rsq32)),crf));
1151 felec = _mm_mul_pd(qq32,_mm_sub_pd(_mm_mul_pd(rinv32,rinvsq32),krf2));
1152
1153 cutoff_mask = _mm_cmplt_pd(rsq32,rcutoff2);
1154
1155 /* Update potential sum for this i atom from the interaction with this j atom. */
1156 velec = _mm_and_pd(velec,cutoff_mask);
1157 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
1158 velecsum = _mm_add_pd(velecsum,velec);
1159
1160 fscal = felec;
1161
1162 fscal = _mm_and_pd(fscal,cutoff_mask);
1163
1164 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1165
1166 /* Calculate temporary vectorial force */
1167 tx = _mm_mul_pd(fscal,dx32);
1168 ty = _mm_mul_pd(fscal,dy32);
1169 tz = _mm_mul_pd(fscal,dz32);
1170
1171 /* Update vectorial force */
1172 fix3 = _mm_add_pd(fix3,tx);
1173 fiy3 = _mm_add_pd(fiy3,ty);
1174 fiz3 = _mm_add_pd(fiz3,tz);
1175
1176 fjx2 = _mm_add_pd(fjx2,tx);
1177 fjy2 = _mm_add_pd(fjy2,ty);
1178 fjz2 = _mm_add_pd(fjz2,tz);
1179
1180 }
1181
1182 /**************************
1183 * CALCULATE INTERACTIONS *
1184 **************************/
1185
1186 if (gmx_mm_any_lt(rsq33,rcutoff2))
1187 {
1188
1189 /* REACTION-FIELD ELECTROSTATICS */
1190 velec = _mm_mul_pd(qq33,_mm_sub_pd(_mm_add_pd(rinv33,_mm_mul_pd(krf,rsq33)),crf));
1191 felec = _mm_mul_pd(qq33,_mm_sub_pd(_mm_mul_pd(rinv33,rinvsq33),krf2));
1192
1193 cutoff_mask = _mm_cmplt_pd(rsq33,rcutoff2);
1194
1195 /* Update potential sum for this i atom from the interaction with this j atom. */
1196 velec = _mm_and_pd(velec,cutoff_mask);
1197 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
1198 velecsum = _mm_add_pd(velecsum,velec);
1199
1200 fscal = felec;
1201
1202 fscal = _mm_and_pd(fscal,cutoff_mask);
1203
1204 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1205
1206 /* Calculate temporary vectorial force */
1207 tx = _mm_mul_pd(fscal,dx33);
1208 ty = _mm_mul_pd(fscal,dy33);
1209 tz = _mm_mul_pd(fscal,dz33);
1210
1211 /* Update vectorial force */
1212 fix3 = _mm_add_pd(fix3,tx);
1213 fiy3 = _mm_add_pd(fiy3,ty);
1214 fiz3 = _mm_add_pd(fiz3,tz);
1215
1216 fjx3 = _mm_add_pd(fjx3,tx);
1217 fjy3 = _mm_add_pd(fjy3,ty);
1218 fjz3 = _mm_add_pd(fjz3,tz);
1219
1220 }
1221
1222 gmx_mm_decrement_4rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
1223
1224 /* Inner loop uses 383 flops */
1225 }
1226
1227 /* End of innermost loop */
1228
1229 gmx_mm_update_iforce_4atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
1230 f+i_coord_offset,fshift+i_shift_offset);
1231
1232 ggid = gid[iidx];
1233 /* Update potential energies */
1234 gmx_mm_update_1pot_pd(velecsum,kernel_data->energygrp_elec+ggid);
1235 gmx_mm_update_1pot_pd(vvdwsum,kernel_data->energygrp_vdw+ggid);
1236
1237 /* Increment number of inner iterations */
1238 inneriter += j_index_end - j_index_start;
1239
1240 /* Outer loop uses 26 flops */
1241 }
1242
1243 /* Increment number of outer iterations */
1244 outeriter += nri;
1245
1246 /* Update outer/inner flops */
1247
1248 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4W4_VF,outeriter*26 + inneriter*383);
1249 }
1250 /*
1251 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwCSTab_GeomW4W4_F_sse2_double
1252 * Electrostatics interaction: ReactionField
1253 * VdW interaction: CubicSplineTable
1254 * Geometry: Water4-Water4
1255 * Calculate force/pot: Force
1256 */
1257 void
1258 nb_kernel_ElecRFCut_VdwCSTab_GeomW4W4_F_sse2_double
1259 (t_nblist * gmx_restrict nlist,
1260 rvec * gmx_restrict xx,
1261 rvec * gmx_restrict ff,
1262 struct t_forcerec * gmx_restrict fr,
1263 t_mdatoms * gmx_restrict mdatoms,
1264 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
1265 t_nrnb * gmx_restrict nrnb)
1266 {
1267 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
1268 * just 0 for non-waters.
1269 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
1270 * jnr indices corresponding to data put in the four positions in the SIMD register.
1271 */
1272 int i_shift_offset,i_coord_offset,outeriter,inneriter;
1273 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
1274 int jnrA,jnrB;
1275 int j_coord_offsetA,j_coord_offsetB;
1276 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
1277 real rcutoff_scalar;
1278 real *shiftvec,*fshift,*x,*f;
1279 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
1280 int vdwioffset0;
1281 __m128d ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
1282 int vdwioffset1;
1283 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
1284 int vdwioffset2;
1285 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
1286 int vdwioffset3;
1287 __m128d ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
1288 int vdwjidx0A,vdwjidx0B;
1289 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
1290 int vdwjidx1A,vdwjidx1B;
1291 __m128d jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
1292 int vdwjidx2A,vdwjidx2B;
1293 __m128d jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
1294 int vdwjidx3A,vdwjidx3B;
1295 __m128d jx3,jy3,jz3,fjx3,fjy3,fjz3,jq3,isaj3;
1296 __m128d dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
1297 __m128d dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
1298 __m128d dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
1299 __m128d dx13,dy13,dz13,rsq13,rinv13,rinvsq13,r13,qq13,c6_13,c12_13;
1300 __m128d dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
1301 __m128d dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
1302 __m128d dx23,dy23,dz23,rsq23,rinv23,rinvsq23,r23,qq23,c6_23,c12_23;
1303 __m128d dx31,dy31,dz31,rsq31,rinv31,rinvsq31,r31,qq31,c6_31,c12_31;
1304 __m128d dx32,dy32,dz32,rsq32,rinv32,rinvsq32,r32,qq32,c6_32,c12_32;
1305 __m128d dx33,dy33,dz33,rsq33,rinv33,rinvsq33,r33,qq33,c6_33,c12_33;
1306 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
1307 real *charge;
1308 int nvdwtype;
1309 __m128d rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
1310 int *vdwtype;
1311 real *vdwparam;
1312 __m128d one_sixth = _mm_set1_pd(1.0/6.0);
1313 __m128d one_twelfth = _mm_set1_pd(1.0/12.0);
1314 __m128i vfitab;
1315 __m128i ifour = _mm_set1_epi32(4);
1316 __m128d rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF;
1317 real *vftab;
1318 __m128d dummy_mask,cutoff_mask;
1319 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
1320 __m128d one = _mm_set1_pd(1.0);
1321 __m128d two = _mm_set1_pd(2.0);
1322 x = xx[0];
1323 f = ff[0];
1324
1325 nri = nlist->nri;
1326 iinr = nlist->iinr;
1327 jindex = nlist->jindex;
1328 jjnr = nlist->jjnr;
1329 shiftidx = nlist->shift;
1330 gid = nlist->gid;
1331 shiftvec = fr->shift_vec[0];
1332 fshift = fr->fshift[0];
1333 facel = _mm_set1_pd(fr->ic->epsfac);
1334 charge = mdatoms->chargeA;
1335 krf = _mm_set1_pd(fr->ic->k_rf);
1336 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
1337 crf = _mm_set1_pd(fr->ic->c_rf);
1338 nvdwtype = fr->ntype;
1339 vdwparam = fr->nbfp;
1340 vdwtype = mdatoms->typeA;
1341
1342 vftab = kernel_data->table_vdw->data;
1343 vftabscale = _mm_set1_pd(kernel_data->table_vdw->scale);
1344
1345 /* Setup water-specific parameters */
1346 inr = nlist->iinr[0];
1347 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
1348 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
1349 iq3 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+3]));
1350 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
1351
1352 jq1 = _mm_set1_pd(charge[inr+1]);
1353 jq2 = _mm_set1_pd(charge[inr+2]);
1354 jq3 = _mm_set1_pd(charge[inr+3]);
1355 vdwjidx0A = 2*vdwtype[inr+0];
1356 c6_00 = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A]);
1357 c12_00 = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A+1]);
1358 qq11 = _mm_mul_pd(iq1,jq1);
1359 qq12 = _mm_mul_pd(iq1,jq2);
1360 qq13 = _mm_mul_pd(iq1,jq3);
1361 qq21 = _mm_mul_pd(iq2,jq1);
1362 qq22 = _mm_mul_pd(iq2,jq2);
1363 qq23 = _mm_mul_pd(iq2,jq3);
1364 qq31 = _mm_mul_pd(iq3,jq1);
1365 qq32 = _mm_mul_pd(iq3,jq2);
1366 qq33 = _mm_mul_pd(iq3,jq3);
1367
1368 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
1369 rcutoff_scalar = fr->ic->rcoulomb;
1370 rcutoff = _mm_set1_pd(rcutoff_scalar);
1371 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
1372
1373 /* Avoid stupid compiler warnings */
1374 jnrA = jnrB = 0;
1375 j_coord_offsetA = 0;
1376 j_coord_offsetB = 0;
1377
1378 outeriter = 0;
1379 inneriter = 0;
1380
1381 /* Start outer loop over neighborlists */
1382 for(iidx=0; iidx<nri; iidx++)
1383 {
1384 /* Load shift vector for this list */
1385 i_shift_offset = DIM*shiftidx[iidx];
1386
1387 /* Load limits for loop over neighbors */
1388 j_index_start = jindex[iidx];
1389 j_index_end = jindex[iidx+1];
1390
1391 /* Get outer coordinate index */
1392 inr = iinr[iidx];
1393 i_coord_offset = DIM*inr;
1394
1395 /* Load i particle coords and add shift vector */
1396 gmx_mm_load_shift_and_4rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
1397 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
1398
1399 fix0 = _mm_setzero_pd();
1400 fiy0 = _mm_setzero_pd();
1401 fiz0 = _mm_setzero_pd();
1402 fix1 = _mm_setzero_pd();
1403 fiy1 = _mm_setzero_pd();
1404 fiz1 = _mm_setzero_pd();
1405 fix2 = _mm_setzero_pd();
1406 fiy2 = _mm_setzero_pd();
1407 fiz2 = _mm_setzero_pd();
1408 fix3 = _mm_setzero_pd();
1409 fiy3 = _mm_setzero_pd();
1410 fiz3 = _mm_setzero_pd();
1411
1412 /* Start inner kernel loop */
1413 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
1414 {
1415
1416 /* Get j neighbor index, and coordinate index */
1417 jnrA = jjnr[jidx];
1418 jnrB = jjnr[jidx+1];
1419 j_coord_offsetA = DIM*jnrA;
1420 j_coord_offsetB = DIM*jnrB;
1421
1422 /* load j atom coordinates */
1423 gmx_mm_load_4rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
1424 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
1425 &jy2,&jz2,&jx3,&jy3,&jz3);
1426
1427 /* Calculate displacement vector */
1428 dx00 = _mm_sub_pd(ix0,jx0);
1429 dy00 = _mm_sub_pd(iy0,jy0);
1430 dz00 = _mm_sub_pd(iz0,jz0);
1431 dx11 = _mm_sub_pd(ix1,jx1);
1432 dy11 = _mm_sub_pd(iy1,jy1);
1433 dz11 = _mm_sub_pd(iz1,jz1);
1434 dx12 = _mm_sub_pd(ix1,jx2);
1435 dy12 = _mm_sub_pd(iy1,jy2);
1436 dz12 = _mm_sub_pd(iz1,jz2);
1437 dx13 = _mm_sub_pd(ix1,jx3);
1438 dy13 = _mm_sub_pd(iy1,jy3);
1439 dz13 = _mm_sub_pd(iz1,jz3);
1440 dx21 = _mm_sub_pd(ix2,jx1);
1441 dy21 = _mm_sub_pd(iy2,jy1);
1442 dz21 = _mm_sub_pd(iz2,jz1);
1443 dx22 = _mm_sub_pd(ix2,jx2);
1444 dy22 = _mm_sub_pd(iy2,jy2);
1445 dz22 = _mm_sub_pd(iz2,jz2);
1446 dx23 = _mm_sub_pd(ix2,jx3);
1447 dy23 = _mm_sub_pd(iy2,jy3);
1448 dz23 = _mm_sub_pd(iz2,jz3);
1449 dx31 = _mm_sub_pd(ix3,jx1);
1450 dy31 = _mm_sub_pd(iy3,jy1);
1451 dz31 = _mm_sub_pd(iz3,jz1);
1452 dx32 = _mm_sub_pd(ix3,jx2);
1453 dy32 = _mm_sub_pd(iy3,jy2);
1454 dz32 = _mm_sub_pd(iz3,jz2);
1455 dx33 = _mm_sub_pd(ix3,jx3);
1456 dy33 = _mm_sub_pd(iy3,jy3);
1457 dz33 = _mm_sub_pd(iz3,jz3);
1458
1459 /* Calculate squared distance and things based on it */
1460 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
1461 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
1462 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
1463 rsq13 = gmx_mm_calc_rsq_pd(dx13,dy13,dz13);
1464 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
1465 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
1466 rsq23 = gmx_mm_calc_rsq_pd(dx23,dy23,dz23);
1467 rsq31 = gmx_mm_calc_rsq_pd(dx31,dy31,dz31);
1468 rsq32 = gmx_mm_calc_rsq_pd(dx32,dy32,dz32);
1469 rsq33 = gmx_mm_calc_rsq_pd(dx33,dy33,dz33);
1470
1471 rinv00 = sse2_invsqrt_d(rsq00);
1472 rinv11 = sse2_invsqrt_d(rsq11);
1473 rinv12 = sse2_invsqrt_d(rsq12);
1474 rinv13 = sse2_invsqrt_d(rsq13);
1475 rinv21 = sse2_invsqrt_d(rsq21);
1476 rinv22 = sse2_invsqrt_d(rsq22);
1477 rinv23 = sse2_invsqrt_d(rsq23);
1478 rinv31 = sse2_invsqrt_d(rsq31);
1479 rinv32 = sse2_invsqrt_d(rsq32);
1480 rinv33 = sse2_invsqrt_d(rsq33);
1481
1482 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
1483 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
1484 rinvsq13 = _mm_mul_pd(rinv13,rinv13);
1485 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
1486 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
1487 rinvsq23 = _mm_mul_pd(rinv23,rinv23);
1488 rinvsq31 = _mm_mul_pd(rinv31,rinv31);
1489 rinvsq32 = _mm_mul_pd(rinv32,rinv32);
1490 rinvsq33 = _mm_mul_pd(rinv33,rinv33);
1491
1492 fjx0 = _mm_setzero_pd();
1493 fjy0 = _mm_setzero_pd();
1494 fjz0 = _mm_setzero_pd();
1495 fjx1 = _mm_setzero_pd();
1496 fjy1 = _mm_setzero_pd();
1497 fjz1 = _mm_setzero_pd();
1498 fjx2 = _mm_setzero_pd();
1499 fjy2 = _mm_setzero_pd();
1500 fjz2 = _mm_setzero_pd();
1501 fjx3 = _mm_setzero_pd();
1502 fjy3 = _mm_setzero_pd();
1503 fjz3 = _mm_setzero_pd();
1504
1505 /**************************
1506 * CALCULATE INTERACTIONS *
1507 **************************/
1508
1509 r00 = _mm_mul_pd(rsq00,rinv00);
1510
1511 /* Calculate table index by multiplying r with table scale and truncate to integer */
1512 rt = _mm_mul_pd(r00,vftabscale);
1513 vfitab = _mm_cvttpd_epi32(rt);
1514 vfeps = _mm_sub_pd(rt,_mm_cvtepi32_pd(vfitab));
1515 vfitab = _mm_slli_epi32(vfitab,3);
1516
1517 /* CUBIC SPLINE TABLE DISPERSION */
1518 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
1519 F = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) );
1520 GMX_MM_TRANSPOSE2_PD(Y,F);
1521 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
1522 H = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) +2);
1523 GMX_MM_TRANSPOSE2_PD(G,H);
1524 Heps = _mm_mul_pd(vfeps,H);
1525 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
1526 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
1527 fvdw6 = _mm_mul_pd(c6_00,FF);
1528
1529 /* CUBIC SPLINE TABLE REPULSION */
1530 vfitab = _mm_add_epi32(vfitab,ifour);
1531 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
1532 F = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) );
1533 GMX_MM_TRANSPOSE2_PD(Y,F);
1534 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
1535 H = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,1) +2);
1536 GMX_MM_TRANSPOSE2_PD(G,H);
1537 Heps = _mm_mul_pd(vfeps,H);
1538 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
1539 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
1540 fvdw12 = _mm_mul_pd(c12_00,FF);
1541 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
1542
1543 fscal = fvdw;
1544
1545 /* Calculate temporary vectorial force */
1546 tx = _mm_mul_pd(fscal,dx00);
1547 ty = _mm_mul_pd(fscal,dy00);
1548 tz = _mm_mul_pd(fscal,dz00);
1549
1550 /* Update vectorial force */
1551 fix0 = _mm_add_pd(fix0,tx);
1552 fiy0 = _mm_add_pd(fiy0,ty);
1553 fiz0 = _mm_add_pd(fiz0,tz);
1554
1555 fjx0 = _mm_add_pd(fjx0,tx);
1556 fjy0 = _mm_add_pd(fjy0,ty);
1557 fjz0 = _mm_add_pd(fjz0,tz);
1558
1559 /**************************
1560 * CALCULATE INTERACTIONS *
1561 **************************/
1562
1563 if (gmx_mm_any_lt(rsq11,rcutoff2))
1564 {
1565
1566 /* REACTION-FIELD ELECTROSTATICS */
1567 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
1568
1569 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
1570
1571 fscal = felec;
1572
1573 fscal = _mm_and_pd(fscal,cutoff_mask);
1574
1575 /* Calculate temporary vectorial force */
1576 tx = _mm_mul_pd(fscal,dx11);
1577 ty = _mm_mul_pd(fscal,dy11);
1578 tz = _mm_mul_pd(fscal,dz11);
1579
1580 /* Update vectorial force */
1581 fix1 = _mm_add_pd(fix1,tx);
1582 fiy1 = _mm_add_pd(fiy1,ty);
1583 fiz1 = _mm_add_pd(fiz1,tz);
1584
1585 fjx1 = _mm_add_pd(fjx1,tx);
1586 fjy1 = _mm_add_pd(fjy1,ty);
1587 fjz1 = _mm_add_pd(fjz1,tz);
1588
1589 }
1590
1591 /**************************
1592 * CALCULATE INTERACTIONS *
1593 **************************/
1594
1595 if (gmx_mm_any_lt(rsq12,rcutoff2))
1596 {
1597
1598 /* REACTION-FIELD ELECTROSTATICS */
1599 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
1600
1601 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
1602
1603 fscal = felec;
1604
1605 fscal = _mm_and_pd(fscal,cutoff_mask);
1606
1607 /* Calculate temporary vectorial force */
1608 tx = _mm_mul_pd(fscal,dx12);
1609 ty = _mm_mul_pd(fscal,dy12);
1610 tz = _mm_mul_pd(fscal,dz12);
1611
1612 /* Update vectorial force */
1613 fix1 = _mm_add_pd(fix1,tx);
1614 fiy1 = _mm_add_pd(fiy1,ty);
1615 fiz1 = _mm_add_pd(fiz1,tz);
1616
1617 fjx2 = _mm_add_pd(fjx2,tx);
1618 fjy2 = _mm_add_pd(fjy2,ty);
1619 fjz2 = _mm_add_pd(fjz2,tz);
1620
1621 }
1622
1623 /**************************
1624 * CALCULATE INTERACTIONS *
1625 **************************/
1626
1627 if (gmx_mm_any_lt(rsq13,rcutoff2))
1628 {
1629
1630 /* REACTION-FIELD ELECTROSTATICS */
1631 felec = _mm_mul_pd(qq13,_mm_sub_pd(_mm_mul_pd(rinv13,rinvsq13),krf2));
1632
1633 cutoff_mask = _mm_cmplt_pd(rsq13,rcutoff2);
1634
1635 fscal = felec;
1636
1637 fscal = _mm_and_pd(fscal,cutoff_mask);
1638
1639 /* Calculate temporary vectorial force */
1640 tx = _mm_mul_pd(fscal,dx13);
1641 ty = _mm_mul_pd(fscal,dy13);
1642 tz = _mm_mul_pd(fscal,dz13);
1643
1644 /* Update vectorial force */
1645 fix1 = _mm_add_pd(fix1,tx);
1646 fiy1 = _mm_add_pd(fiy1,ty);
1647 fiz1 = _mm_add_pd(fiz1,tz);
1648
1649 fjx3 = _mm_add_pd(fjx3,tx);
1650 fjy3 = _mm_add_pd(fjy3,ty);
1651 fjz3 = _mm_add_pd(fjz3,tz);
1652
1653 }
1654
1655 /**************************
1656 * CALCULATE INTERACTIONS *
1657 **************************/
1658
1659 if (gmx_mm_any_lt(rsq21,rcutoff2))
1660 {
1661
1662 /* REACTION-FIELD ELECTROSTATICS */
1663 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
1664
1665 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
1666
1667 fscal = felec;
1668
1669 fscal = _mm_and_pd(fscal,cutoff_mask);
1670
1671 /* Calculate temporary vectorial force */
1672 tx = _mm_mul_pd(fscal,dx21);
1673 ty = _mm_mul_pd(fscal,dy21);
1674 tz = _mm_mul_pd(fscal,dz21);
1675
1676 /* Update vectorial force */
1677 fix2 = _mm_add_pd(fix2,tx);
1678 fiy2 = _mm_add_pd(fiy2,ty);
1679 fiz2 = _mm_add_pd(fiz2,tz);
1680
1681 fjx1 = _mm_add_pd(fjx1,tx);
1682 fjy1 = _mm_add_pd(fjy1,ty);
1683 fjz1 = _mm_add_pd(fjz1,tz);
1684
1685 }
1686
1687 /**************************
1688 * CALCULATE INTERACTIONS *
1689 **************************/
1690
1691 if (gmx_mm_any_lt(rsq22,rcutoff2))
1692 {
1693
1694 /* REACTION-FIELD ELECTROSTATICS */
1695 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
1696
1697 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
1698
1699 fscal = felec;
1700
1701 fscal = _mm_and_pd(fscal,cutoff_mask);
1702
1703 /* Calculate temporary vectorial force */
1704 tx = _mm_mul_pd(fscal,dx22);
1705 ty = _mm_mul_pd(fscal,dy22);
1706 tz = _mm_mul_pd(fscal,dz22);
1707
1708 /* Update vectorial force */
1709 fix2 = _mm_add_pd(fix2,tx);
1710 fiy2 = _mm_add_pd(fiy2,ty);
1711 fiz2 = _mm_add_pd(fiz2,tz);
1712
1713 fjx2 = _mm_add_pd(fjx2,tx);
1714 fjy2 = _mm_add_pd(fjy2,ty);
1715 fjz2 = _mm_add_pd(fjz2,tz);
1716
1717 }
1718
1719 /**************************
1720 * CALCULATE INTERACTIONS *
1721 **************************/
1722
1723 if (gmx_mm_any_lt(rsq23,rcutoff2))
1724 {
1725
1726 /* REACTION-FIELD ELECTROSTATICS */
1727 felec = _mm_mul_pd(qq23,_mm_sub_pd(_mm_mul_pd(rinv23,rinvsq23),krf2));
1728
1729 cutoff_mask = _mm_cmplt_pd(rsq23,rcutoff2);
1730
1731 fscal = felec;
1732
1733 fscal = _mm_and_pd(fscal,cutoff_mask);
1734
1735 /* Calculate temporary vectorial force */
1736 tx = _mm_mul_pd(fscal,dx23);
1737 ty = _mm_mul_pd(fscal,dy23);
1738 tz = _mm_mul_pd(fscal,dz23);
1739
1740 /* Update vectorial force */
1741 fix2 = _mm_add_pd(fix2,tx);
1742 fiy2 = _mm_add_pd(fiy2,ty);
1743 fiz2 = _mm_add_pd(fiz2,tz);
1744
1745 fjx3 = _mm_add_pd(fjx3,tx);
1746 fjy3 = _mm_add_pd(fjy3,ty);
1747 fjz3 = _mm_add_pd(fjz3,tz);
1748
1749 }
1750
1751 /**************************
1752 * CALCULATE INTERACTIONS *
1753 **************************/
1754
1755 if (gmx_mm_any_lt(rsq31,rcutoff2))
1756 {
1757
1758 /* REACTION-FIELD ELECTROSTATICS */
1759 felec = _mm_mul_pd(qq31,_mm_sub_pd(_mm_mul_pd(rinv31,rinvsq31),krf2));
1760
1761 cutoff_mask = _mm_cmplt_pd(rsq31,rcutoff2);
1762
1763 fscal = felec;
1764
1765 fscal = _mm_and_pd(fscal,cutoff_mask);
1766
1767 /* Calculate temporary vectorial force */
1768 tx = _mm_mul_pd(fscal,dx31);
1769 ty = _mm_mul_pd(fscal,dy31);
1770 tz = _mm_mul_pd(fscal,dz31);
1771
1772 /* Update vectorial force */
1773 fix3 = _mm_add_pd(fix3,tx);
1774 fiy3 = _mm_add_pd(fiy3,ty);
1775 fiz3 = _mm_add_pd(fiz3,tz);
1776
1777 fjx1 = _mm_add_pd(fjx1,tx);
1778 fjy1 = _mm_add_pd(fjy1,ty);
1779 fjz1 = _mm_add_pd(fjz1,tz);
1780
1781 }
1782
1783 /**************************
1784 * CALCULATE INTERACTIONS *
1785 **************************/
1786
1787 if (gmx_mm_any_lt(rsq32,rcutoff2))
1788 {
1789
1790 /* REACTION-FIELD ELECTROSTATICS */
1791 felec = _mm_mul_pd(qq32,_mm_sub_pd(_mm_mul_pd(rinv32,rinvsq32),krf2));
1792
1793 cutoff_mask = _mm_cmplt_pd(rsq32,rcutoff2);
1794
1795 fscal = felec;
1796
1797 fscal = _mm_and_pd(fscal,cutoff_mask);
1798
1799 /* Calculate temporary vectorial force */
1800 tx = _mm_mul_pd(fscal,dx32);
1801 ty = _mm_mul_pd(fscal,dy32);
1802 tz = _mm_mul_pd(fscal,dz32);
1803
1804 /* Update vectorial force */
1805 fix3 = _mm_add_pd(fix3,tx);
1806 fiy3 = _mm_add_pd(fiy3,ty);
1807 fiz3 = _mm_add_pd(fiz3,tz);
1808
1809 fjx2 = _mm_add_pd(fjx2,tx);
1810 fjy2 = _mm_add_pd(fjy2,ty);
1811 fjz2 = _mm_add_pd(fjz2,tz);
1812
1813 }
1814
1815 /**************************
1816 * CALCULATE INTERACTIONS *
1817 **************************/
1818
1819 if (gmx_mm_any_lt(rsq33,rcutoff2))
1820 {
1821
1822 /* REACTION-FIELD ELECTROSTATICS */
1823 felec = _mm_mul_pd(qq33,_mm_sub_pd(_mm_mul_pd(rinv33,rinvsq33),krf2));
1824
1825 cutoff_mask = _mm_cmplt_pd(rsq33,rcutoff2);
1826
1827 fscal = felec;
1828
1829 fscal = _mm_and_pd(fscal,cutoff_mask);
1830
1831 /* Calculate temporary vectorial force */
1832 tx = _mm_mul_pd(fscal,dx33);
1833 ty = _mm_mul_pd(fscal,dy33);
1834 tz = _mm_mul_pd(fscal,dz33);
1835
1836 /* Update vectorial force */
1837 fix3 = _mm_add_pd(fix3,tx);
1838 fiy3 = _mm_add_pd(fiy3,ty);
1839 fiz3 = _mm_add_pd(fiz3,tz);
1840
1841 fjx3 = _mm_add_pd(fjx3,tx);
1842 fjy3 = _mm_add_pd(fjy3,ty);
1843 fjz3 = _mm_add_pd(fjz3,tz);
1844
1845 }
1846
1847 gmx_mm_decrement_4rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
1848
1849 /* Inner loop uses 321 flops */
1850 }
1851
1852 if(jidx<j_index_end)
1853 {
1854
1855 jnrA = jjnr[jidx];
1856 j_coord_offsetA = DIM*jnrA;
1857
1858 /* load j atom coordinates */
1859 gmx_mm_load_4rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
1860 &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
1861 &jy2,&jz2,&jx3,&jy3,&jz3);
1862
1863 /* Calculate displacement vector */
1864 dx00 = _mm_sub_pd(ix0,jx0);
1865 dy00 = _mm_sub_pd(iy0,jy0);
1866 dz00 = _mm_sub_pd(iz0,jz0);
1867 dx11 = _mm_sub_pd(ix1,jx1);
1868 dy11 = _mm_sub_pd(iy1,jy1);
1869 dz11 = _mm_sub_pd(iz1,jz1);
1870 dx12 = _mm_sub_pd(ix1,jx2);
1871 dy12 = _mm_sub_pd(iy1,jy2);
1872 dz12 = _mm_sub_pd(iz1,jz2);
1873 dx13 = _mm_sub_pd(ix1,jx3);
1874 dy13 = _mm_sub_pd(iy1,jy3);
1875 dz13 = _mm_sub_pd(iz1,jz3);
1876 dx21 = _mm_sub_pd(ix2,jx1);
1877 dy21 = _mm_sub_pd(iy2,jy1);
1878 dz21 = _mm_sub_pd(iz2,jz1);
1879 dx22 = _mm_sub_pd(ix2,jx2);
1880 dy22 = _mm_sub_pd(iy2,jy2);
1881 dz22 = _mm_sub_pd(iz2,jz2);
1882 dx23 = _mm_sub_pd(ix2,jx3);
1883 dy23 = _mm_sub_pd(iy2,jy3);
1884 dz23 = _mm_sub_pd(iz2,jz3);
1885 dx31 = _mm_sub_pd(ix3,jx1);
1886 dy31 = _mm_sub_pd(iy3,jy1);
1887 dz31 = _mm_sub_pd(iz3,jz1);
1888 dx32 = _mm_sub_pd(ix3,jx2);
1889 dy32 = _mm_sub_pd(iy3,jy2);
1890 dz32 = _mm_sub_pd(iz3,jz2);
1891 dx33 = _mm_sub_pd(ix3,jx3);
1892 dy33 = _mm_sub_pd(iy3,jy3);
1893 dz33 = _mm_sub_pd(iz3,jz3);
1894
1895 /* Calculate squared distance and things based on it */
1896 rsq00 = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
1897 rsq11 = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
1898 rsq12 = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
1899 rsq13 = gmx_mm_calc_rsq_pd(dx13,dy13,dz13);
1900 rsq21 = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
1901 rsq22 = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
1902 rsq23 = gmx_mm_calc_rsq_pd(dx23,dy23,dz23);
1903 rsq31 = gmx_mm_calc_rsq_pd(dx31,dy31,dz31);
1904 rsq32 = gmx_mm_calc_rsq_pd(dx32,dy32,dz32);
1905 rsq33 = gmx_mm_calc_rsq_pd(dx33,dy33,dz33);
1906
1907 rinv00 = sse2_invsqrt_d(rsq00);
1908 rinv11 = sse2_invsqrt_d(rsq11);
1909 rinv12 = sse2_invsqrt_d(rsq12);
1910 rinv13 = sse2_invsqrt_d(rsq13);
1911 rinv21 = sse2_invsqrt_d(rsq21);
1912 rinv22 = sse2_invsqrt_d(rsq22);
1913 rinv23 = sse2_invsqrt_d(rsq23);
1914 rinv31 = sse2_invsqrt_d(rsq31);
1915 rinv32 = sse2_invsqrt_d(rsq32);
1916 rinv33 = sse2_invsqrt_d(rsq33);
1917
1918 rinvsq11 = _mm_mul_pd(rinv11,rinv11);
1919 rinvsq12 = _mm_mul_pd(rinv12,rinv12);
1920 rinvsq13 = _mm_mul_pd(rinv13,rinv13);
1921 rinvsq21 = _mm_mul_pd(rinv21,rinv21);
1922 rinvsq22 = _mm_mul_pd(rinv22,rinv22);
1923 rinvsq23 = _mm_mul_pd(rinv23,rinv23);
1924 rinvsq31 = _mm_mul_pd(rinv31,rinv31);
1925 rinvsq32 = _mm_mul_pd(rinv32,rinv32);
1926 rinvsq33 = _mm_mul_pd(rinv33,rinv33);
1927
1928 fjx0 = _mm_setzero_pd();
1929 fjy0 = _mm_setzero_pd();
1930 fjz0 = _mm_setzero_pd();
1931 fjx1 = _mm_setzero_pd();
1932 fjy1 = _mm_setzero_pd();
1933 fjz1 = _mm_setzero_pd();
1934 fjx2 = _mm_setzero_pd();
1935 fjy2 = _mm_setzero_pd();
1936 fjz2 = _mm_setzero_pd();
1937 fjx3 = _mm_setzero_pd();
1938 fjy3 = _mm_setzero_pd();
1939 fjz3 = _mm_setzero_pd();
1940
1941 /**************************
1942 * CALCULATE INTERACTIONS *
1943 **************************/
1944
1945 r00 = _mm_mul_pd(rsq00,rinv00);
1946
1947 /* Calculate table index by multiplying r with table scale and truncate to integer */
1948 rt = _mm_mul_pd(r00,vftabscale);
1949 vfitab = _mm_cvttpd_epi32(rt);
1950 vfeps = _mm_sub_pd(rt,_mm_cvtepi32_pd(vfitab));
1951 vfitab = _mm_slli_epi32(vfitab,3);
1952
1953 /* CUBIC SPLINE TABLE DISPERSION */
1954 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
1955 F = _mm_setzero_pd();
1956 GMX_MM_TRANSPOSE2_PD(Y,F);
1957 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
1958 H = _mm_setzero_pd();
1959 GMX_MM_TRANSPOSE2_PD(G,H);
1960 Heps = _mm_mul_pd(vfeps,H);
1961 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
1962 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
1963 fvdw6 = _mm_mul_pd(c6_00,FF);
1964
1965 /* CUBIC SPLINE TABLE REPULSION */
1966 vfitab = _mm_add_epi32(vfitab,ifour);
1967 Y = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) );
1968 F = _mm_setzero_pd();
1969 GMX_MM_TRANSPOSE2_PD(Y,F);
1970 G = _mm_load_pd( vftab + gmx_mm_extract_epi32(vfitab,0) +2);
1971 H = _mm_setzero_pd();
1972 GMX_MM_TRANSPOSE2_PD(G,H);
1973 Heps = _mm_mul_pd(vfeps,H);
1974 Fp = _mm_add_pd(F,_mm_mul_pd(vfeps,_mm_add_pd(G,Heps)));
1975 FF = _mm_add_pd(Fp,_mm_mul_pd(vfeps,_mm_add_pd(G,_mm_add_pd(Heps,Heps))));
1976 fvdw12 = _mm_mul_pd(c12_00,FF);
1977 fvdw = _mm_xor_pd(signbit,_mm_mul_pd(_mm_add_pd(fvdw6,fvdw12),_mm_mul_pd(vftabscale,rinv00)));
1978
1979 fscal = fvdw;
1980
1981 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1982
1983 /* Calculate temporary vectorial force */
1984 tx = _mm_mul_pd(fscal,dx00);
1985 ty = _mm_mul_pd(fscal,dy00);
1986 tz = _mm_mul_pd(fscal,dz00);
1987
1988 /* Update vectorial force */
1989 fix0 = _mm_add_pd(fix0,tx);
1990 fiy0 = _mm_add_pd(fiy0,ty);
1991 fiz0 = _mm_add_pd(fiz0,tz);
1992
1993 fjx0 = _mm_add_pd(fjx0,tx);
1994 fjy0 = _mm_add_pd(fjy0,ty);
1995 fjz0 = _mm_add_pd(fjz0,tz);
1996
1997 /**************************
1998 * CALCULATE INTERACTIONS *
1999 **************************/
2000
2001 if (gmx_mm_any_lt(rsq11,rcutoff2))
2002 {
2003
2004 /* REACTION-FIELD ELECTROSTATICS */
2005 felec = _mm_mul_pd(qq11,_mm_sub_pd(_mm_mul_pd(rinv11,rinvsq11),krf2));
2006
2007 cutoff_mask = _mm_cmplt_pd(rsq11,rcutoff2);
2008
2009 fscal = felec;
2010
2011 fscal = _mm_and_pd(fscal,cutoff_mask);
2012
2013 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2014
2015 /* Calculate temporary vectorial force */
2016 tx = _mm_mul_pd(fscal,dx11);
2017 ty = _mm_mul_pd(fscal,dy11);
2018 tz = _mm_mul_pd(fscal,dz11);
2019
2020 /* Update vectorial force */
2021 fix1 = _mm_add_pd(fix1,tx);
2022 fiy1 = _mm_add_pd(fiy1,ty);
2023 fiz1 = _mm_add_pd(fiz1,tz);
2024
2025 fjx1 = _mm_add_pd(fjx1,tx);
2026 fjy1 = _mm_add_pd(fjy1,ty);
2027 fjz1 = _mm_add_pd(fjz1,tz);
2028
2029 }
2030
2031 /**************************
2032 * CALCULATE INTERACTIONS *
2033 **************************/
2034
2035 if (gmx_mm_any_lt(rsq12,rcutoff2))
2036 {
2037
2038 /* REACTION-FIELD ELECTROSTATICS */
2039 felec = _mm_mul_pd(qq12,_mm_sub_pd(_mm_mul_pd(rinv12,rinvsq12),krf2));
2040
2041 cutoff_mask = _mm_cmplt_pd(rsq12,rcutoff2);
2042
2043 fscal = felec;
2044
2045 fscal = _mm_and_pd(fscal,cutoff_mask);
2046
2047 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2048
2049 /* Calculate temporary vectorial force */
2050 tx = _mm_mul_pd(fscal,dx12);
2051 ty = _mm_mul_pd(fscal,dy12);
2052 tz = _mm_mul_pd(fscal,dz12);
2053
2054 /* Update vectorial force */
2055 fix1 = _mm_add_pd(fix1,tx);
2056 fiy1 = _mm_add_pd(fiy1,ty);
2057 fiz1 = _mm_add_pd(fiz1,tz);
2058
2059 fjx2 = _mm_add_pd(fjx2,tx);
2060 fjy2 = _mm_add_pd(fjy2,ty);
2061 fjz2 = _mm_add_pd(fjz2,tz);
2062
2063 }
2064
2065 /**************************
2066 * CALCULATE INTERACTIONS *
2067 **************************/
2068
2069 if (gmx_mm_any_lt(rsq13,rcutoff2))
2070 {
2071
2072 /* REACTION-FIELD ELECTROSTATICS */
2073 felec = _mm_mul_pd(qq13,_mm_sub_pd(_mm_mul_pd(rinv13,rinvsq13),krf2));
2074
2075 cutoff_mask = _mm_cmplt_pd(rsq13,rcutoff2);
2076
2077 fscal = felec;
2078
2079 fscal = _mm_and_pd(fscal,cutoff_mask);
2080
2081 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2082
2083 /* Calculate temporary vectorial force */
2084 tx = _mm_mul_pd(fscal,dx13);
2085 ty = _mm_mul_pd(fscal,dy13);
2086 tz = _mm_mul_pd(fscal,dz13);
2087
2088 /* Update vectorial force */
2089 fix1 = _mm_add_pd(fix1,tx);
2090 fiy1 = _mm_add_pd(fiy1,ty);
2091 fiz1 = _mm_add_pd(fiz1,tz);
2092
2093 fjx3 = _mm_add_pd(fjx3,tx);
2094 fjy3 = _mm_add_pd(fjy3,ty);
2095 fjz3 = _mm_add_pd(fjz3,tz);
2096
2097 }
2098
2099 /**************************
2100 * CALCULATE INTERACTIONS *
2101 **************************/
2102
2103 if (gmx_mm_any_lt(rsq21,rcutoff2))
2104 {
2105
2106 /* REACTION-FIELD ELECTROSTATICS */
2107 felec = _mm_mul_pd(qq21,_mm_sub_pd(_mm_mul_pd(rinv21,rinvsq21),krf2));
2108
2109 cutoff_mask = _mm_cmplt_pd(rsq21,rcutoff2);
2110
2111 fscal = felec;
2112
2113 fscal = _mm_and_pd(fscal,cutoff_mask);
2114
2115 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2116
2117 /* Calculate temporary vectorial force */
2118 tx = _mm_mul_pd(fscal,dx21);
2119 ty = _mm_mul_pd(fscal,dy21);
2120 tz = _mm_mul_pd(fscal,dz21);
2121
2122 /* Update vectorial force */
2123 fix2 = _mm_add_pd(fix2,tx);
2124 fiy2 = _mm_add_pd(fiy2,ty);
2125 fiz2 = _mm_add_pd(fiz2,tz);
2126
2127 fjx1 = _mm_add_pd(fjx1,tx);
2128 fjy1 = _mm_add_pd(fjy1,ty);
2129 fjz1 = _mm_add_pd(fjz1,tz);
2130
2131 }
2132
2133 /**************************
2134 * CALCULATE INTERACTIONS *
2135 **************************/
2136
2137 if (gmx_mm_any_lt(rsq22,rcutoff2))
2138 {
2139
2140 /* REACTION-FIELD ELECTROSTATICS */
2141 felec = _mm_mul_pd(qq22,_mm_sub_pd(_mm_mul_pd(rinv22,rinvsq22),krf2));
2142
2143 cutoff_mask = _mm_cmplt_pd(rsq22,rcutoff2);
2144
2145 fscal = felec;
2146
2147 fscal = _mm_and_pd(fscal,cutoff_mask);
2148
2149 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2150
2151 /* Calculate temporary vectorial force */
2152 tx = _mm_mul_pd(fscal,dx22);
2153 ty = _mm_mul_pd(fscal,dy22);
2154 tz = _mm_mul_pd(fscal,dz22);
2155
2156 /* Update vectorial force */
2157 fix2 = _mm_add_pd(fix2,tx);
2158 fiy2 = _mm_add_pd(fiy2,ty);
2159 fiz2 = _mm_add_pd(fiz2,tz);
2160
2161 fjx2 = _mm_add_pd(fjx2,tx);
2162 fjy2 = _mm_add_pd(fjy2,ty);
2163 fjz2 = _mm_add_pd(fjz2,tz);
2164
2165 }
2166
2167 /**************************
2168 * CALCULATE INTERACTIONS *
2169 **************************/
2170
2171 if (gmx_mm_any_lt(rsq23,rcutoff2))
2172 {
2173
2174 /* REACTION-FIELD ELECTROSTATICS */
2175 felec = _mm_mul_pd(qq23,_mm_sub_pd(_mm_mul_pd(rinv23,rinvsq23),krf2));
2176
2177 cutoff_mask = _mm_cmplt_pd(rsq23,rcutoff2);
2178
2179 fscal = felec;
2180
2181 fscal = _mm_and_pd(fscal,cutoff_mask);
2182
2183 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2184
2185 /* Calculate temporary vectorial force */
2186 tx = _mm_mul_pd(fscal,dx23);
2187 ty = _mm_mul_pd(fscal,dy23);
2188 tz = _mm_mul_pd(fscal,dz23);
2189
2190 /* Update vectorial force */
2191 fix2 = _mm_add_pd(fix2,tx);
2192 fiy2 = _mm_add_pd(fiy2,ty);
2193 fiz2 = _mm_add_pd(fiz2,tz);
2194
2195 fjx3 = _mm_add_pd(fjx3,tx);
2196 fjy3 = _mm_add_pd(fjy3,ty);
2197 fjz3 = _mm_add_pd(fjz3,tz);
2198
2199 }
2200
2201 /**************************
2202 * CALCULATE INTERACTIONS *
2203 **************************/
2204
2205 if (gmx_mm_any_lt(rsq31,rcutoff2))
2206 {
2207
2208 /* REACTION-FIELD ELECTROSTATICS */
2209 felec = _mm_mul_pd(qq31,_mm_sub_pd(_mm_mul_pd(rinv31,rinvsq31),krf2));
2210
2211 cutoff_mask = _mm_cmplt_pd(rsq31,rcutoff2);
2212
2213 fscal = felec;
2214
2215 fscal = _mm_and_pd(fscal,cutoff_mask);
2216
2217 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2218
2219 /* Calculate temporary vectorial force */
2220 tx = _mm_mul_pd(fscal,dx31);
2221 ty = _mm_mul_pd(fscal,dy31);
2222 tz = _mm_mul_pd(fscal,dz31);
2223
2224 /* Update vectorial force */
2225 fix3 = _mm_add_pd(fix3,tx);
2226 fiy3 = _mm_add_pd(fiy3,ty);
2227 fiz3 = _mm_add_pd(fiz3,tz);
2228
2229 fjx1 = _mm_add_pd(fjx1,tx);
2230 fjy1 = _mm_add_pd(fjy1,ty);
2231 fjz1 = _mm_add_pd(fjz1,tz);
2232
2233 }
2234
2235 /**************************
2236 * CALCULATE INTERACTIONS *
2237 **************************/
2238
2239 if (gmx_mm_any_lt(rsq32,rcutoff2))
2240 {
2241
2242 /* REACTION-FIELD ELECTROSTATICS */
2243 felec = _mm_mul_pd(qq32,_mm_sub_pd(_mm_mul_pd(rinv32,rinvsq32),krf2));
2244
2245 cutoff_mask = _mm_cmplt_pd(rsq32,rcutoff2);
2246
2247 fscal = felec;
2248
2249 fscal = _mm_and_pd(fscal,cutoff_mask);
2250
2251 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2252
2253 /* Calculate temporary vectorial force */
2254 tx = _mm_mul_pd(fscal,dx32);
2255 ty = _mm_mul_pd(fscal,dy32);
2256 tz = _mm_mul_pd(fscal,dz32);
2257
2258 /* Update vectorial force */
2259 fix3 = _mm_add_pd(fix3,tx);
2260 fiy3 = _mm_add_pd(fiy3,ty);
2261 fiz3 = _mm_add_pd(fiz3,tz);
2262
2263 fjx2 = _mm_add_pd(fjx2,tx);
2264 fjy2 = _mm_add_pd(fjy2,ty);
2265 fjz2 = _mm_add_pd(fjz2,tz);
2266
2267 }
2268
2269 /**************************
2270 * CALCULATE INTERACTIONS *
2271 **************************/
2272
2273 if (gmx_mm_any_lt(rsq33,rcutoff2))
2274 {
2275
2276 /* REACTION-FIELD ELECTROSTATICS */
2277 felec = _mm_mul_pd(qq33,_mm_sub_pd(_mm_mul_pd(rinv33,rinvsq33),krf2));
2278
2279 cutoff_mask = _mm_cmplt_pd(rsq33,rcutoff2);
2280
2281 fscal = felec;
2282
2283 fscal = _mm_and_pd(fscal,cutoff_mask);
2284
2285 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
2286
2287 /* Calculate temporary vectorial force */
2288 tx = _mm_mul_pd(fscal,dx33);
2289 ty = _mm_mul_pd(fscal,dy33);
2290 tz = _mm_mul_pd(fscal,dz33);
2291
2292 /* Update vectorial force */
2293 fix3 = _mm_add_pd(fix3,tx);
2294 fiy3 = _mm_add_pd(fiy3,ty);
2295 fiz3 = _mm_add_pd(fiz3,tz);
2296
2297 fjx3 = _mm_add_pd(fjx3,tx);
2298 fjy3 = _mm_add_pd(fjy3,ty);
2299 fjz3 = _mm_add_pd(fjz3,tz);
2300
2301 }
2302
2303 gmx_mm_decrement_4rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
2304
2305 /* Inner loop uses 321 flops */
2306 }
2307
2308 /* End of innermost loop */
2309
2310 gmx_mm_update_iforce_4atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
2311 f+i_coord_offset,fshift+i_shift_offset);
2312
2313 /* Increment number of inner iterations */
2314 inneriter += j_index_end - j_index_start;
2315
2316 /* Outer loop uses 24 flops */
2317 }
2318
2319 /* Increment number of outer iterations */
2320 outeriter += nri;
2321
2322 /* Update outer/inner flops */
2323
2324 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4W4_F,outeriter*24 + inneriter*321);
2325 }
The ultimate molecular dynamics simulation package