?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Removed simple.h from nb_kernel_sse2_XX
[gromacs.git] / src / gromacs / gmxlib / nonbonded / nb_kernel_sse2_single / nb_kernel_ElecEwSh_VdwLJSh_GeomW4P1_sse2_single.c
blob0c3d07aa017a0f16fa11ef102fbec2001c2692ac
1 /*
2 * This file is part of the GROMACS molecular simulation package.
3 *
4 * Copyright (c) 2012,2013,2014,2015, 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_single 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/math/vec.h"
46 #include "gromacs/legacyheaders/nrnb.h"
47
48 #include "gromacs/simd/math_x86_sse2_single.h"
49 #include "kernelutil_x86_sse2_single.h"
50
51 /*
52 * Gromacs nonbonded kernel: nb_kernel_ElecEwSh_VdwLJSh_GeomW4P1_VF_sse2_single
53 * Electrostatics interaction: Ewald
54 * VdW interaction: LennardJones
55 * Geometry: Water4-Particle
56 * Calculate force/pot: PotentialAndForce
57 */
58 void
59 nb_kernel_ElecEwSh_VdwLJSh_GeomW4P1_VF_sse2_single
60 (t_nblist * gmx_restrict nlist,
61 rvec * gmx_restrict xx,
62 rvec * gmx_restrict ff,
63 t_forcerec * gmx_restrict fr,
64 t_mdatoms * gmx_restrict mdatoms,
65 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
66 t_nrnb * gmx_restrict nrnb)
67 {
68 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
69 * just 0 for non-waters.
70 * Suffixes A,B,C,D refer to j loop unrolling done with SSE, e.g. for the four different
71 * jnr indices corresponding to data put in the four positions in the SIMD register.
72 */
73 int i_shift_offset,i_coord_offset,outeriter,inneriter;
74 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
75 int jnrA,jnrB,jnrC,jnrD;
76 int jnrlistA,jnrlistB,jnrlistC,jnrlistD;
77 int j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
78 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
79 real rcutoff_scalar;
80 real *shiftvec,*fshift,*x,*f;
81 real *fjptrA,*fjptrB,*fjptrC,*fjptrD;
82 real scratch[4*DIM];
83 __m128 tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
84 int vdwioffset0;
85 __m128 ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
86 int vdwioffset1;
87 __m128 ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
88 int vdwioffset2;
89 __m128 ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
90 int vdwioffset3;
91 __m128 ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
92 int vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
93 __m128 jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
94 __m128 dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
95 __m128 dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
96 __m128 dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
97 __m128 dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
98 __m128 velec,felec,velecsum,facel,crf,krf,krf2;
99 real *charge;
100 int nvdwtype;
101 __m128 rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
102 int *vdwtype;
103 real *vdwparam;
104 __m128 one_sixth = _mm_set1_ps(1.0/6.0);
105 __m128 one_twelfth = _mm_set1_ps(1.0/12.0);
106 __m128i ewitab;
107 __m128 ewtabscale,eweps,sh_ewald,ewrt,ewtabhalfspace,ewtabF,ewtabFn,ewtabD,ewtabV;
108 real *ewtab;
109 __m128 dummy_mask,cutoff_mask;
110 __m128 signbit = _mm_castsi128_ps( _mm_set1_epi32(0x80000000) );
111 __m128 one = _mm_set1_ps(1.0);
112 __m128 two = _mm_set1_ps(2.0);
113 x = xx[0];
114 f = ff[0];
115
116 nri = nlist->nri;
117 iinr = nlist->iinr;
118 jindex = nlist->jindex;
119 jjnr = nlist->jjnr;
120 shiftidx = nlist->shift;
121 gid = nlist->gid;
122 shiftvec = fr->shift_vec[0];
123 fshift = fr->fshift[0];
124 facel = _mm_set1_ps(fr->epsfac);
125 charge = mdatoms->chargeA;
126 nvdwtype = fr->ntype;
127 vdwparam = fr->nbfp;
128 vdwtype = mdatoms->typeA;
129
130 sh_ewald = _mm_set1_ps(fr->ic->sh_ewald);
131 ewtab = fr->ic->tabq_coul_FDV0;
132 ewtabscale = _mm_set1_ps(fr->ic->tabq_scale);
133 ewtabhalfspace = _mm_set1_ps(0.5/fr->ic->tabq_scale);
134
135 /* Setup water-specific parameters */
136 inr = nlist->iinr[0];
137 iq1 = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+1]));
138 iq2 = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+2]));
139 iq3 = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+3]));
140 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
141
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->rcoulomb;
144 rcutoff = _mm_set1_ps(rcutoff_scalar);
145 rcutoff2 = _mm_mul_ps(rcutoff,rcutoff);
146
147 sh_vdw_invrcut6 = _mm_set1_ps(fr->ic->sh_invrc6);
148 rvdw = _mm_set1_ps(fr->rvdw);
149
150 /* Avoid stupid compiler warnings */
151 jnrA = jnrB = jnrC = jnrD = 0;
152 j_coord_offsetA = 0;
153 j_coord_offsetB = 0;
154 j_coord_offsetC = 0;
155 j_coord_offsetD = 0;
156
157 outeriter = 0;
158 inneriter = 0;
159
160 for(iidx=0;iidx<4*DIM;iidx++)
161 {
162 scratch[iidx] = 0.0;
163 }
164
165 /* Start outer loop over neighborlists */
166 for(iidx=0; iidx<nri; iidx++)
167 {
168 /* Load shift vector for this list */
169 i_shift_offset = DIM*shiftidx[iidx];
170
171 /* Load limits for loop over neighbors */
172 j_index_start = jindex[iidx];
173 j_index_end = jindex[iidx+1];
174
175 /* Get outer coordinate index */
176 inr = iinr[iidx];
177 i_coord_offset = DIM*inr;
178
179 /* Load i particle coords and add shift vector */
180 gmx_mm_load_shift_and_4rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,
181 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
182
183 fix0 = _mm_setzero_ps();
184 fiy0 = _mm_setzero_ps();
185 fiz0 = _mm_setzero_ps();
186 fix1 = _mm_setzero_ps();
187 fiy1 = _mm_setzero_ps();
188 fiz1 = _mm_setzero_ps();
189 fix2 = _mm_setzero_ps();
190 fiy2 = _mm_setzero_ps();
191 fiz2 = _mm_setzero_ps();
192 fix3 = _mm_setzero_ps();
193 fiy3 = _mm_setzero_ps();
194 fiz3 = _mm_setzero_ps();
195
196 /* Reset potential sums */
197 velecsum = _mm_setzero_ps();
198 vvdwsum = _mm_setzero_ps();
199
200 /* Start inner kernel loop */
201 for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
202 {
203
204 /* Get j neighbor index, and coordinate index */
205 jnrA = jjnr[jidx];
206 jnrB = jjnr[jidx+1];
207 jnrC = jjnr[jidx+2];
208 jnrD = jjnr[jidx+3];
209 j_coord_offsetA = DIM*jnrA;
210 j_coord_offsetB = DIM*jnrB;
211 j_coord_offsetC = DIM*jnrC;
212 j_coord_offsetD = DIM*jnrD;
213
214 /* load j atom coordinates */
215 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
216 x+j_coord_offsetC,x+j_coord_offsetD,
217 &jx0,&jy0,&jz0);
218
219 /* Calculate displacement vector */
220 dx00 = _mm_sub_ps(ix0,jx0);
221 dy00 = _mm_sub_ps(iy0,jy0);
222 dz00 = _mm_sub_ps(iz0,jz0);
223 dx10 = _mm_sub_ps(ix1,jx0);
224 dy10 = _mm_sub_ps(iy1,jy0);
225 dz10 = _mm_sub_ps(iz1,jz0);
226 dx20 = _mm_sub_ps(ix2,jx0);
227 dy20 = _mm_sub_ps(iy2,jy0);
228 dz20 = _mm_sub_ps(iz2,jz0);
229 dx30 = _mm_sub_ps(ix3,jx0);
230 dy30 = _mm_sub_ps(iy3,jy0);
231 dz30 = _mm_sub_ps(iz3,jz0);
232
233 /* Calculate squared distance and things based on it */
234 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
235 rsq10 = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
236 rsq20 = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
237 rsq30 = gmx_mm_calc_rsq_ps(dx30,dy30,dz30);
238
239 rinv10 = gmx_mm_invsqrt_ps(rsq10);
240 rinv20 = gmx_mm_invsqrt_ps(rsq20);
241 rinv30 = gmx_mm_invsqrt_ps(rsq30);
242
243 rinvsq00 = gmx_mm_inv_ps(rsq00);
244 rinvsq10 = _mm_mul_ps(rinv10,rinv10);
245 rinvsq20 = _mm_mul_ps(rinv20,rinv20);
246 rinvsq30 = _mm_mul_ps(rinv30,rinv30);
247
248 /* Load parameters for j particles */
249 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
250 charge+jnrC+0,charge+jnrD+0);
251 vdwjidx0A = 2*vdwtype[jnrA+0];
252 vdwjidx0B = 2*vdwtype[jnrB+0];
253 vdwjidx0C = 2*vdwtype[jnrC+0];
254 vdwjidx0D = 2*vdwtype[jnrD+0];
255
256 fjx0 = _mm_setzero_ps();
257 fjy0 = _mm_setzero_ps();
258 fjz0 = _mm_setzero_ps();
259
260 /**************************
261 * CALCULATE INTERACTIONS *
262 **************************/
263
264 if (gmx_mm_any_lt(rsq00,rcutoff2))
265 {
266
267 /* Compute parameters for interactions between i and j atoms */
268 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
269 vdwparam+vdwioffset0+vdwjidx0B,
270 vdwparam+vdwioffset0+vdwjidx0C,
271 vdwparam+vdwioffset0+vdwjidx0D,
272 &c6_00,&c12_00);
273
274 /* LENNARD-JONES DISPERSION/REPULSION */
275
276 rinvsix = _mm_mul_ps(_mm_mul_ps(rinvsq00,rinvsq00),rinvsq00);
277 vvdw6 = _mm_mul_ps(c6_00,rinvsix);
278 vvdw12 = _mm_mul_ps(c12_00,_mm_mul_ps(rinvsix,rinvsix));
279 vvdw = _mm_sub_ps(_mm_mul_ps( _mm_sub_ps(vvdw12 , _mm_mul_ps(c12_00,_mm_mul_ps(sh_vdw_invrcut6,sh_vdw_invrcut6))), one_twelfth) ,
280 _mm_mul_ps( _mm_sub_ps(vvdw6,_mm_mul_ps(c6_00,sh_vdw_invrcut6)),one_sixth));
281 fvdw = _mm_mul_ps(_mm_sub_ps(vvdw12,vvdw6),rinvsq00);
282
283 cutoff_mask = _mm_cmplt_ps(rsq00,rcutoff2);
284
285 /* Update potential sum for this i atom from the interaction with this j atom. */
286 vvdw = _mm_and_ps(vvdw,cutoff_mask);
287 vvdwsum = _mm_add_ps(vvdwsum,vvdw);
288
289 fscal = fvdw;
290
291 fscal = _mm_and_ps(fscal,cutoff_mask);
292
293 /* Calculate temporary vectorial force */
294 tx = _mm_mul_ps(fscal,dx00);
295 ty = _mm_mul_ps(fscal,dy00);
296 tz = _mm_mul_ps(fscal,dz00);
297
298 /* Update vectorial force */
299 fix0 = _mm_add_ps(fix0,tx);
300 fiy0 = _mm_add_ps(fiy0,ty);
301 fiz0 = _mm_add_ps(fiz0,tz);
302
303 fjx0 = _mm_add_ps(fjx0,tx);
304 fjy0 = _mm_add_ps(fjy0,ty);
305 fjz0 = _mm_add_ps(fjz0,tz);
306
307 }
308
309 /**************************
310 * CALCULATE INTERACTIONS *
311 **************************/
312
313 if (gmx_mm_any_lt(rsq10,rcutoff2))
314 {
315
316 r10 = _mm_mul_ps(rsq10,rinv10);
317
318 /* Compute parameters for interactions between i and j atoms */
319 qq10 = _mm_mul_ps(iq1,jq0);
320
321 /* EWALD ELECTROSTATICS */
322
323 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
324 ewrt = _mm_mul_ps(r10,ewtabscale);
325 ewitab = _mm_cvttps_epi32(ewrt);
326 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
327 ewitab = _mm_slli_epi32(ewitab,2);
328 ewtabF = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,0) );
329 ewtabD = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,1) );
330 ewtabV = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,2) );
331 ewtabFn = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,3) );
332 _MM_TRANSPOSE4_PS(ewtabF,ewtabD,ewtabV,ewtabFn);
333 felec = _mm_add_ps(ewtabF,_mm_mul_ps(eweps,ewtabD));
334 velec = _mm_sub_ps(ewtabV,_mm_mul_ps(_mm_mul_ps(ewtabhalfspace,eweps),_mm_add_ps(ewtabF,felec)));
335 velec = _mm_mul_ps(qq10,_mm_sub_ps(_mm_sub_ps(rinv10,sh_ewald),velec));
336 felec = _mm_mul_ps(_mm_mul_ps(qq10,rinv10),_mm_sub_ps(rinvsq10,felec));
337
338 cutoff_mask = _mm_cmplt_ps(rsq10,rcutoff2);
339
340 /* Update potential sum for this i atom from the interaction with this j atom. */
341 velec = _mm_and_ps(velec,cutoff_mask);
342 velecsum = _mm_add_ps(velecsum,velec);
343
344 fscal = felec;
345
346 fscal = _mm_and_ps(fscal,cutoff_mask);
347
348 /* Calculate temporary vectorial force */
349 tx = _mm_mul_ps(fscal,dx10);
350 ty = _mm_mul_ps(fscal,dy10);
351 tz = _mm_mul_ps(fscal,dz10);
352
353 /* Update vectorial force */
354 fix1 = _mm_add_ps(fix1,tx);
355 fiy1 = _mm_add_ps(fiy1,ty);
356 fiz1 = _mm_add_ps(fiz1,tz);
357
358 fjx0 = _mm_add_ps(fjx0,tx);
359 fjy0 = _mm_add_ps(fjy0,ty);
360 fjz0 = _mm_add_ps(fjz0,tz);
361
362 }
363
364 /**************************
365 * CALCULATE INTERACTIONS *
366 **************************/
367
368 if (gmx_mm_any_lt(rsq20,rcutoff2))
369 {
370
371 r20 = _mm_mul_ps(rsq20,rinv20);
372
373 /* Compute parameters for interactions between i and j atoms */
374 qq20 = _mm_mul_ps(iq2,jq0);
375
376 /* EWALD ELECTROSTATICS */
377
378 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
379 ewrt = _mm_mul_ps(r20,ewtabscale);
380 ewitab = _mm_cvttps_epi32(ewrt);
381 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
382 ewitab = _mm_slli_epi32(ewitab,2);
383 ewtabF = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,0) );
384 ewtabD = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,1) );
385 ewtabV = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,2) );
386 ewtabFn = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,3) );
387 _MM_TRANSPOSE4_PS(ewtabF,ewtabD,ewtabV,ewtabFn);
388 felec = _mm_add_ps(ewtabF,_mm_mul_ps(eweps,ewtabD));
389 velec = _mm_sub_ps(ewtabV,_mm_mul_ps(_mm_mul_ps(ewtabhalfspace,eweps),_mm_add_ps(ewtabF,felec)));
390 velec = _mm_mul_ps(qq20,_mm_sub_ps(_mm_sub_ps(rinv20,sh_ewald),velec));
391 felec = _mm_mul_ps(_mm_mul_ps(qq20,rinv20),_mm_sub_ps(rinvsq20,felec));
392
393 cutoff_mask = _mm_cmplt_ps(rsq20,rcutoff2);
394
395 /* Update potential sum for this i atom from the interaction with this j atom. */
396 velec = _mm_and_ps(velec,cutoff_mask);
397 velecsum = _mm_add_ps(velecsum,velec);
398
399 fscal = felec;
400
401 fscal = _mm_and_ps(fscal,cutoff_mask);
402
403 /* Calculate temporary vectorial force */
404 tx = _mm_mul_ps(fscal,dx20);
405 ty = _mm_mul_ps(fscal,dy20);
406 tz = _mm_mul_ps(fscal,dz20);
407
408 /* Update vectorial force */
409 fix2 = _mm_add_ps(fix2,tx);
410 fiy2 = _mm_add_ps(fiy2,ty);
411 fiz2 = _mm_add_ps(fiz2,tz);
412
413 fjx0 = _mm_add_ps(fjx0,tx);
414 fjy0 = _mm_add_ps(fjy0,ty);
415 fjz0 = _mm_add_ps(fjz0,tz);
416
417 }
418
419 /**************************
420 * CALCULATE INTERACTIONS *
421 **************************/
422
423 if (gmx_mm_any_lt(rsq30,rcutoff2))
424 {
425
426 r30 = _mm_mul_ps(rsq30,rinv30);
427
428 /* Compute parameters for interactions between i and j atoms */
429 qq30 = _mm_mul_ps(iq3,jq0);
430
431 /* EWALD ELECTROSTATICS */
432
433 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
434 ewrt = _mm_mul_ps(r30,ewtabscale);
435 ewitab = _mm_cvttps_epi32(ewrt);
436 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
437 ewitab = _mm_slli_epi32(ewitab,2);
438 ewtabF = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,0) );
439 ewtabD = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,1) );
440 ewtabV = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,2) );
441 ewtabFn = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,3) );
442 _MM_TRANSPOSE4_PS(ewtabF,ewtabD,ewtabV,ewtabFn);
443 felec = _mm_add_ps(ewtabF,_mm_mul_ps(eweps,ewtabD));
444 velec = _mm_sub_ps(ewtabV,_mm_mul_ps(_mm_mul_ps(ewtabhalfspace,eweps),_mm_add_ps(ewtabF,felec)));
445 velec = _mm_mul_ps(qq30,_mm_sub_ps(_mm_sub_ps(rinv30,sh_ewald),velec));
446 felec = _mm_mul_ps(_mm_mul_ps(qq30,rinv30),_mm_sub_ps(rinvsq30,felec));
447
448 cutoff_mask = _mm_cmplt_ps(rsq30,rcutoff2);
449
450 /* Update potential sum for this i atom from the interaction with this j atom. */
451 velec = _mm_and_ps(velec,cutoff_mask);
452 velecsum = _mm_add_ps(velecsum,velec);
453
454 fscal = felec;
455
456 fscal = _mm_and_ps(fscal,cutoff_mask);
457
458 /* Calculate temporary vectorial force */
459 tx = _mm_mul_ps(fscal,dx30);
460 ty = _mm_mul_ps(fscal,dy30);
461 tz = _mm_mul_ps(fscal,dz30);
462
463 /* Update vectorial force */
464 fix3 = _mm_add_ps(fix3,tx);
465 fiy3 = _mm_add_ps(fiy3,ty);
466 fiz3 = _mm_add_ps(fiz3,tz);
467
468 fjx0 = _mm_add_ps(fjx0,tx);
469 fjy0 = _mm_add_ps(fjy0,ty);
470 fjz0 = _mm_add_ps(fjz0,tz);
471
472 }
473
474 fjptrA = f+j_coord_offsetA;
475 fjptrB = f+j_coord_offsetB;
476 fjptrC = f+j_coord_offsetC;
477 fjptrD = f+j_coord_offsetD;
478
479 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
480
481 /* Inner loop uses 179 flops */
482 }
483
484 if(jidx<j_index_end)
485 {
486
487 /* Get j neighbor index, and coordinate index */
488 jnrlistA = jjnr[jidx];
489 jnrlistB = jjnr[jidx+1];
490 jnrlistC = jjnr[jidx+2];
491 jnrlistD = jjnr[jidx+3];
492 /* Sign of each element will be negative for non-real atoms.
493 * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
494 * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
495 */
496 dummy_mask = gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
497 jnrA = (jnrlistA>=0) ? jnrlistA : 0;
498 jnrB = (jnrlistB>=0) ? jnrlistB : 0;
499 jnrC = (jnrlistC>=0) ? jnrlistC : 0;
500 jnrD = (jnrlistD>=0) ? jnrlistD : 0;
501 j_coord_offsetA = DIM*jnrA;
502 j_coord_offsetB = DIM*jnrB;
503 j_coord_offsetC = DIM*jnrC;
504 j_coord_offsetD = DIM*jnrD;
505
506 /* load j atom coordinates */
507 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
508 x+j_coord_offsetC,x+j_coord_offsetD,
509 &jx0,&jy0,&jz0);
510
511 /* Calculate displacement vector */
512 dx00 = _mm_sub_ps(ix0,jx0);
513 dy00 = _mm_sub_ps(iy0,jy0);
514 dz00 = _mm_sub_ps(iz0,jz0);
515 dx10 = _mm_sub_ps(ix1,jx0);
516 dy10 = _mm_sub_ps(iy1,jy0);
517 dz10 = _mm_sub_ps(iz1,jz0);
518 dx20 = _mm_sub_ps(ix2,jx0);
519 dy20 = _mm_sub_ps(iy2,jy0);
520 dz20 = _mm_sub_ps(iz2,jz0);
521 dx30 = _mm_sub_ps(ix3,jx0);
522 dy30 = _mm_sub_ps(iy3,jy0);
523 dz30 = _mm_sub_ps(iz3,jz0);
524
525 /* Calculate squared distance and things based on it */
526 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
527 rsq10 = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
528 rsq20 = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
529 rsq30 = gmx_mm_calc_rsq_ps(dx30,dy30,dz30);
530
531 rinv10 = gmx_mm_invsqrt_ps(rsq10);
532 rinv20 = gmx_mm_invsqrt_ps(rsq20);
533 rinv30 = gmx_mm_invsqrt_ps(rsq30);
534
535 rinvsq00 = gmx_mm_inv_ps(rsq00);
536 rinvsq10 = _mm_mul_ps(rinv10,rinv10);
537 rinvsq20 = _mm_mul_ps(rinv20,rinv20);
538 rinvsq30 = _mm_mul_ps(rinv30,rinv30);
539
540 /* Load parameters for j particles */
541 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
542 charge+jnrC+0,charge+jnrD+0);
543 vdwjidx0A = 2*vdwtype[jnrA+0];
544 vdwjidx0B = 2*vdwtype[jnrB+0];
545 vdwjidx0C = 2*vdwtype[jnrC+0];
546 vdwjidx0D = 2*vdwtype[jnrD+0];
547
548 fjx0 = _mm_setzero_ps();
549 fjy0 = _mm_setzero_ps();
550 fjz0 = _mm_setzero_ps();
551
552 /**************************
553 * CALCULATE INTERACTIONS *
554 **************************/
555
556 if (gmx_mm_any_lt(rsq00,rcutoff2))
557 {
558
559 /* Compute parameters for interactions between i and j atoms */
560 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
561 vdwparam+vdwioffset0+vdwjidx0B,
562 vdwparam+vdwioffset0+vdwjidx0C,
563 vdwparam+vdwioffset0+vdwjidx0D,
564 &c6_00,&c12_00);
565
566 /* LENNARD-JONES DISPERSION/REPULSION */
567
568 rinvsix = _mm_mul_ps(_mm_mul_ps(rinvsq00,rinvsq00),rinvsq00);
569 vvdw6 = _mm_mul_ps(c6_00,rinvsix);
570 vvdw12 = _mm_mul_ps(c12_00,_mm_mul_ps(rinvsix,rinvsix));
571 vvdw = _mm_sub_ps(_mm_mul_ps( _mm_sub_ps(vvdw12 , _mm_mul_ps(c12_00,_mm_mul_ps(sh_vdw_invrcut6,sh_vdw_invrcut6))), one_twelfth) ,
572 _mm_mul_ps( _mm_sub_ps(vvdw6,_mm_mul_ps(c6_00,sh_vdw_invrcut6)),one_sixth));
573 fvdw = _mm_mul_ps(_mm_sub_ps(vvdw12,vvdw6),rinvsq00);
574
575 cutoff_mask = _mm_cmplt_ps(rsq00,rcutoff2);
576
577 /* Update potential sum for this i atom from the interaction with this j atom. */
578 vvdw = _mm_and_ps(vvdw,cutoff_mask);
579 vvdw = _mm_andnot_ps(dummy_mask,vvdw);
580 vvdwsum = _mm_add_ps(vvdwsum,vvdw);
581
582 fscal = fvdw;
583
584 fscal = _mm_and_ps(fscal,cutoff_mask);
585
586 fscal = _mm_andnot_ps(dummy_mask,fscal);
587
588 /* Calculate temporary vectorial force */
589 tx = _mm_mul_ps(fscal,dx00);
590 ty = _mm_mul_ps(fscal,dy00);
591 tz = _mm_mul_ps(fscal,dz00);
592
593 /* Update vectorial force */
594 fix0 = _mm_add_ps(fix0,tx);
595 fiy0 = _mm_add_ps(fiy0,ty);
596 fiz0 = _mm_add_ps(fiz0,tz);
597
598 fjx0 = _mm_add_ps(fjx0,tx);
599 fjy0 = _mm_add_ps(fjy0,ty);
600 fjz0 = _mm_add_ps(fjz0,tz);
601
602 }
603
604 /**************************
605 * CALCULATE INTERACTIONS *
606 **************************/
607
608 if (gmx_mm_any_lt(rsq10,rcutoff2))
609 {
610
611 r10 = _mm_mul_ps(rsq10,rinv10);
612 r10 = _mm_andnot_ps(dummy_mask,r10);
613
614 /* Compute parameters for interactions between i and j atoms */
615 qq10 = _mm_mul_ps(iq1,jq0);
616
617 /* EWALD ELECTROSTATICS */
618
619 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
620 ewrt = _mm_mul_ps(r10,ewtabscale);
621 ewitab = _mm_cvttps_epi32(ewrt);
622 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
623 ewitab = _mm_slli_epi32(ewitab,2);
624 ewtabF = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,0) );
625 ewtabD = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,1) );
626 ewtabV = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,2) );
627 ewtabFn = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,3) );
628 _MM_TRANSPOSE4_PS(ewtabF,ewtabD,ewtabV,ewtabFn);
629 felec = _mm_add_ps(ewtabF,_mm_mul_ps(eweps,ewtabD));
630 velec = _mm_sub_ps(ewtabV,_mm_mul_ps(_mm_mul_ps(ewtabhalfspace,eweps),_mm_add_ps(ewtabF,felec)));
631 velec = _mm_mul_ps(qq10,_mm_sub_ps(_mm_sub_ps(rinv10,sh_ewald),velec));
632 felec = _mm_mul_ps(_mm_mul_ps(qq10,rinv10),_mm_sub_ps(rinvsq10,felec));
633
634 cutoff_mask = _mm_cmplt_ps(rsq10,rcutoff2);
635
636 /* Update potential sum for this i atom from the interaction with this j atom. */
637 velec = _mm_and_ps(velec,cutoff_mask);
638 velec = _mm_andnot_ps(dummy_mask,velec);
639 velecsum = _mm_add_ps(velecsum,velec);
640
641 fscal = felec;
642
643 fscal = _mm_and_ps(fscal,cutoff_mask);
644
645 fscal = _mm_andnot_ps(dummy_mask,fscal);
646
647 /* Calculate temporary vectorial force */
648 tx = _mm_mul_ps(fscal,dx10);
649 ty = _mm_mul_ps(fscal,dy10);
650 tz = _mm_mul_ps(fscal,dz10);
651
652 /* Update vectorial force */
653 fix1 = _mm_add_ps(fix1,tx);
654 fiy1 = _mm_add_ps(fiy1,ty);
655 fiz1 = _mm_add_ps(fiz1,tz);
656
657 fjx0 = _mm_add_ps(fjx0,tx);
658 fjy0 = _mm_add_ps(fjy0,ty);
659 fjz0 = _mm_add_ps(fjz0,tz);
660
661 }
662
663 /**************************
664 * CALCULATE INTERACTIONS *
665 **************************/
666
667 if (gmx_mm_any_lt(rsq20,rcutoff2))
668 {
669
670 r20 = _mm_mul_ps(rsq20,rinv20);
671 r20 = _mm_andnot_ps(dummy_mask,r20);
672
673 /* Compute parameters for interactions between i and j atoms */
674 qq20 = _mm_mul_ps(iq2,jq0);
675
676 /* EWALD ELECTROSTATICS */
677
678 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
679 ewrt = _mm_mul_ps(r20,ewtabscale);
680 ewitab = _mm_cvttps_epi32(ewrt);
681 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
682 ewitab = _mm_slli_epi32(ewitab,2);
683 ewtabF = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,0) );
684 ewtabD = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,1) );
685 ewtabV = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,2) );
686 ewtabFn = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,3) );
687 _MM_TRANSPOSE4_PS(ewtabF,ewtabD,ewtabV,ewtabFn);
688 felec = _mm_add_ps(ewtabF,_mm_mul_ps(eweps,ewtabD));
689 velec = _mm_sub_ps(ewtabV,_mm_mul_ps(_mm_mul_ps(ewtabhalfspace,eweps),_mm_add_ps(ewtabF,felec)));
690 velec = _mm_mul_ps(qq20,_mm_sub_ps(_mm_sub_ps(rinv20,sh_ewald),velec));
691 felec = _mm_mul_ps(_mm_mul_ps(qq20,rinv20),_mm_sub_ps(rinvsq20,felec));
692
693 cutoff_mask = _mm_cmplt_ps(rsq20,rcutoff2);
694
695 /* Update potential sum for this i atom from the interaction with this j atom. */
696 velec = _mm_and_ps(velec,cutoff_mask);
697 velec = _mm_andnot_ps(dummy_mask,velec);
698 velecsum = _mm_add_ps(velecsum,velec);
699
700 fscal = felec;
701
702 fscal = _mm_and_ps(fscal,cutoff_mask);
703
704 fscal = _mm_andnot_ps(dummy_mask,fscal);
705
706 /* Calculate temporary vectorial force */
707 tx = _mm_mul_ps(fscal,dx20);
708 ty = _mm_mul_ps(fscal,dy20);
709 tz = _mm_mul_ps(fscal,dz20);
710
711 /* Update vectorial force */
712 fix2 = _mm_add_ps(fix2,tx);
713 fiy2 = _mm_add_ps(fiy2,ty);
714 fiz2 = _mm_add_ps(fiz2,tz);
715
716 fjx0 = _mm_add_ps(fjx0,tx);
717 fjy0 = _mm_add_ps(fjy0,ty);
718 fjz0 = _mm_add_ps(fjz0,tz);
719
720 }
721
722 /**************************
723 * CALCULATE INTERACTIONS *
724 **************************/
725
726 if (gmx_mm_any_lt(rsq30,rcutoff2))
727 {
728
729 r30 = _mm_mul_ps(rsq30,rinv30);
730 r30 = _mm_andnot_ps(dummy_mask,r30);
731
732 /* Compute parameters for interactions between i and j atoms */
733 qq30 = _mm_mul_ps(iq3,jq0);
734
735 /* EWALD ELECTROSTATICS */
736
737 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
738 ewrt = _mm_mul_ps(r30,ewtabscale);
739 ewitab = _mm_cvttps_epi32(ewrt);
740 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
741 ewitab = _mm_slli_epi32(ewitab,2);
742 ewtabF = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,0) );
743 ewtabD = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,1) );
744 ewtabV = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,2) );
745 ewtabFn = _mm_load_ps( ewtab + gmx_mm_extract_epi32(ewitab,3) );
746 _MM_TRANSPOSE4_PS(ewtabF,ewtabD,ewtabV,ewtabFn);
747 felec = _mm_add_ps(ewtabF,_mm_mul_ps(eweps,ewtabD));
748 velec = _mm_sub_ps(ewtabV,_mm_mul_ps(_mm_mul_ps(ewtabhalfspace,eweps),_mm_add_ps(ewtabF,felec)));
749 velec = _mm_mul_ps(qq30,_mm_sub_ps(_mm_sub_ps(rinv30,sh_ewald),velec));
750 felec = _mm_mul_ps(_mm_mul_ps(qq30,rinv30),_mm_sub_ps(rinvsq30,felec));
751
752 cutoff_mask = _mm_cmplt_ps(rsq30,rcutoff2);
753
754 /* Update potential sum for this i atom from the interaction with this j atom. */
755 velec = _mm_and_ps(velec,cutoff_mask);
756 velec = _mm_andnot_ps(dummy_mask,velec);
757 velecsum = _mm_add_ps(velecsum,velec);
758
759 fscal = felec;
760
761 fscal = _mm_and_ps(fscal,cutoff_mask);
762
763 fscal = _mm_andnot_ps(dummy_mask,fscal);
764
765 /* Calculate temporary vectorial force */
766 tx = _mm_mul_ps(fscal,dx30);
767 ty = _mm_mul_ps(fscal,dy30);
768 tz = _mm_mul_ps(fscal,dz30);
769
770 /* Update vectorial force */
771 fix3 = _mm_add_ps(fix3,tx);
772 fiy3 = _mm_add_ps(fiy3,ty);
773 fiz3 = _mm_add_ps(fiz3,tz);
774
775 fjx0 = _mm_add_ps(fjx0,tx);
776 fjy0 = _mm_add_ps(fjy0,ty);
777 fjz0 = _mm_add_ps(fjz0,tz);
778
779 }
780
781 fjptrA = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
782 fjptrB = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
783 fjptrC = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
784 fjptrD = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
785
786 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
787
788 /* Inner loop uses 182 flops */
789 }
790
791 /* End of innermost loop */
792
793 gmx_mm_update_iforce_4atom_swizzle_ps(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
794 f+i_coord_offset,fshift+i_shift_offset);
795
796 ggid = gid[iidx];
797 /* Update potential energies */
798 gmx_mm_update_1pot_ps(velecsum,kernel_data->energygrp_elec+ggid);
799 gmx_mm_update_1pot_ps(vvdwsum,kernel_data->energygrp_vdw+ggid);
800
801 /* Increment number of inner iterations */
802 inneriter += j_index_end - j_index_start;
803
804 /* Outer loop uses 26 flops */
805 }
806
807 /* Increment number of outer iterations */
808 outeriter += nri;
809
810 /* Update outer/inner flops */
811
812 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4_VF,outeriter*26 + inneriter*182);
813 }
814 /*
815 * Gromacs nonbonded kernel: nb_kernel_ElecEwSh_VdwLJSh_GeomW4P1_F_sse2_single
816 * Electrostatics interaction: Ewald
817 * VdW interaction: LennardJones
818 * Geometry: Water4-Particle
819 * Calculate force/pot: Force
820 */
821 void
822 nb_kernel_ElecEwSh_VdwLJSh_GeomW4P1_F_sse2_single
823 (t_nblist * gmx_restrict nlist,
824 rvec * gmx_restrict xx,
825 rvec * gmx_restrict ff,
826 t_forcerec * gmx_restrict fr,
827 t_mdatoms * gmx_restrict mdatoms,
828 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
829 t_nrnb * gmx_restrict nrnb)
830 {
831 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
832 * just 0 for non-waters.
833 * Suffixes A,B,C,D refer to j loop unrolling done with SSE, e.g. for the four different
834 * jnr indices corresponding to data put in the four positions in the SIMD register.
835 */
836 int i_shift_offset,i_coord_offset,outeriter,inneriter;
837 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
838 int jnrA,jnrB,jnrC,jnrD;
839 int jnrlistA,jnrlistB,jnrlistC,jnrlistD;
840 int j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
841 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
842 real rcutoff_scalar;
843 real *shiftvec,*fshift,*x,*f;
844 real *fjptrA,*fjptrB,*fjptrC,*fjptrD;
845 real scratch[4*DIM];
846 __m128 tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
847 int vdwioffset0;
848 __m128 ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
849 int vdwioffset1;
850 __m128 ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
851 int vdwioffset2;
852 __m128 ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
853 int vdwioffset3;
854 __m128 ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
855 int vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
856 __m128 jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
857 __m128 dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
858 __m128 dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
859 __m128 dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
860 __m128 dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
861 __m128 velec,felec,velecsum,facel,crf,krf,krf2;
862 real *charge;
863 int nvdwtype;
864 __m128 rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
865 int *vdwtype;
866 real *vdwparam;
867 __m128 one_sixth = _mm_set1_ps(1.0/6.0);
868 __m128 one_twelfth = _mm_set1_ps(1.0/12.0);
869 __m128i ewitab;
870 __m128 ewtabscale,eweps,sh_ewald,ewrt,ewtabhalfspace,ewtabF,ewtabFn,ewtabD,ewtabV;
871 real *ewtab;
872 __m128 dummy_mask,cutoff_mask;
873 __m128 signbit = _mm_castsi128_ps( _mm_set1_epi32(0x80000000) );
874 __m128 one = _mm_set1_ps(1.0);
875 __m128 two = _mm_set1_ps(2.0);
876 x = xx[0];
877 f = ff[0];
878
879 nri = nlist->nri;
880 iinr = nlist->iinr;
881 jindex = nlist->jindex;
882 jjnr = nlist->jjnr;
883 shiftidx = nlist->shift;
884 gid = nlist->gid;
885 shiftvec = fr->shift_vec[0];
886 fshift = fr->fshift[0];
887 facel = _mm_set1_ps(fr->epsfac);
888 charge = mdatoms->chargeA;
889 nvdwtype = fr->ntype;
890 vdwparam = fr->nbfp;
891 vdwtype = mdatoms->typeA;
892
893 sh_ewald = _mm_set1_ps(fr->ic->sh_ewald);
894 ewtab = fr->ic->tabq_coul_F;
895 ewtabscale = _mm_set1_ps(fr->ic->tabq_scale);
896 ewtabhalfspace = _mm_set1_ps(0.5/fr->ic->tabq_scale);
897
898 /* Setup water-specific parameters */
899 inr = nlist->iinr[0];
900 iq1 = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+1]));
901 iq2 = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+2]));
902 iq3 = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+3]));
903 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
904
905 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
906 rcutoff_scalar = fr->rcoulomb;
907 rcutoff = _mm_set1_ps(rcutoff_scalar);
908 rcutoff2 = _mm_mul_ps(rcutoff,rcutoff);
909
910 sh_vdw_invrcut6 = _mm_set1_ps(fr->ic->sh_invrc6);
911 rvdw = _mm_set1_ps(fr->rvdw);
912
913 /* Avoid stupid compiler warnings */
914 jnrA = jnrB = jnrC = jnrD = 0;
915 j_coord_offsetA = 0;
916 j_coord_offsetB = 0;
917 j_coord_offsetC = 0;
918 j_coord_offsetD = 0;
919
920 outeriter = 0;
921 inneriter = 0;
922
923 for(iidx=0;iidx<4*DIM;iidx++)
924 {
925 scratch[iidx] = 0.0;
926 }
927
928 /* Start outer loop over neighborlists */
929 for(iidx=0; iidx<nri; iidx++)
930 {
931 /* Load shift vector for this list */
932 i_shift_offset = DIM*shiftidx[iidx];
933
934 /* Load limits for loop over neighbors */
935 j_index_start = jindex[iidx];
936 j_index_end = jindex[iidx+1];
937
938 /* Get outer coordinate index */
939 inr = iinr[iidx];
940 i_coord_offset = DIM*inr;
941
942 /* Load i particle coords and add shift vector */
943 gmx_mm_load_shift_and_4rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,
944 &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
945
946 fix0 = _mm_setzero_ps();
947 fiy0 = _mm_setzero_ps();
948 fiz0 = _mm_setzero_ps();
949 fix1 = _mm_setzero_ps();
950 fiy1 = _mm_setzero_ps();
951 fiz1 = _mm_setzero_ps();
952 fix2 = _mm_setzero_ps();
953 fiy2 = _mm_setzero_ps();
954 fiz2 = _mm_setzero_ps();
955 fix3 = _mm_setzero_ps();
956 fiy3 = _mm_setzero_ps();
957 fiz3 = _mm_setzero_ps();
958
959 /* Start inner kernel loop */
960 for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
961 {
962
963 /* Get j neighbor index, and coordinate index */
964 jnrA = jjnr[jidx];
965 jnrB = jjnr[jidx+1];
966 jnrC = jjnr[jidx+2];
967 jnrD = jjnr[jidx+3];
968 j_coord_offsetA = DIM*jnrA;
969 j_coord_offsetB = DIM*jnrB;
970 j_coord_offsetC = DIM*jnrC;
971 j_coord_offsetD = DIM*jnrD;
972
973 /* load j atom coordinates */
974 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
975 x+j_coord_offsetC,x+j_coord_offsetD,
976 &jx0,&jy0,&jz0);
977
978 /* Calculate displacement vector */
979 dx00 = _mm_sub_ps(ix0,jx0);
980 dy00 = _mm_sub_ps(iy0,jy0);
981 dz00 = _mm_sub_ps(iz0,jz0);
982 dx10 = _mm_sub_ps(ix1,jx0);
983 dy10 = _mm_sub_ps(iy1,jy0);
984 dz10 = _mm_sub_ps(iz1,jz0);
985 dx20 = _mm_sub_ps(ix2,jx0);
986 dy20 = _mm_sub_ps(iy2,jy0);
987 dz20 = _mm_sub_ps(iz2,jz0);
988 dx30 = _mm_sub_ps(ix3,jx0);
989 dy30 = _mm_sub_ps(iy3,jy0);
990 dz30 = _mm_sub_ps(iz3,jz0);
991
992 /* Calculate squared distance and things based on it */
993 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
994 rsq10 = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
995 rsq20 = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
996 rsq30 = gmx_mm_calc_rsq_ps(dx30,dy30,dz30);
997
998 rinv10 = gmx_mm_invsqrt_ps(rsq10);
999 rinv20 = gmx_mm_invsqrt_ps(rsq20);
1000 rinv30 = gmx_mm_invsqrt_ps(rsq30);
1001
1002 rinvsq00 = gmx_mm_inv_ps(rsq00);
1003 rinvsq10 = _mm_mul_ps(rinv10,rinv10);
1004 rinvsq20 = _mm_mul_ps(rinv20,rinv20);
1005 rinvsq30 = _mm_mul_ps(rinv30,rinv30);
1006
1007 /* Load parameters for j particles */
1008 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
1009 charge+jnrC+0,charge+jnrD+0);
1010 vdwjidx0A = 2*vdwtype[jnrA+0];
1011 vdwjidx0B = 2*vdwtype[jnrB+0];
1012 vdwjidx0C = 2*vdwtype[jnrC+0];
1013 vdwjidx0D = 2*vdwtype[jnrD+0];
1014
1015 fjx0 = _mm_setzero_ps();
1016 fjy0 = _mm_setzero_ps();
1017 fjz0 = _mm_setzero_ps();
1018
1019 /**************************
1020 * CALCULATE INTERACTIONS *
1021 **************************/
1022
1023 if (gmx_mm_any_lt(rsq00,rcutoff2))
1024 {
1025
1026 /* Compute parameters for interactions between i and j atoms */
1027 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
1028 vdwparam+vdwioffset0+vdwjidx0B,
1029 vdwparam+vdwioffset0+vdwjidx0C,
1030 vdwparam+vdwioffset0+vdwjidx0D,
1031 &c6_00,&c12_00);
1032
1033 /* LENNARD-JONES DISPERSION/REPULSION */
1034
1035 rinvsix = _mm_mul_ps(_mm_mul_ps(rinvsq00,rinvsq00),rinvsq00);
1036 fvdw = _mm_mul_ps(_mm_sub_ps(_mm_mul_ps(c12_00,rinvsix),c6_00),_mm_mul_ps(rinvsix,rinvsq00));
1037
1038 cutoff_mask = _mm_cmplt_ps(rsq00,rcutoff2);
1039
1040 fscal = fvdw;
1041
1042 fscal = _mm_and_ps(fscal,cutoff_mask);
1043
1044 /* Calculate temporary vectorial force */
1045 tx = _mm_mul_ps(fscal,dx00);
1046 ty = _mm_mul_ps(fscal,dy00);
1047 tz = _mm_mul_ps(fscal,dz00);
1048
1049 /* Update vectorial force */
1050 fix0 = _mm_add_ps(fix0,tx);
1051 fiy0 = _mm_add_ps(fiy0,ty);
1052 fiz0 = _mm_add_ps(fiz0,tz);
1053
1054 fjx0 = _mm_add_ps(fjx0,tx);
1055 fjy0 = _mm_add_ps(fjy0,ty);
1056 fjz0 = _mm_add_ps(fjz0,tz);
1057
1058 }
1059
1060 /**************************
1061 * CALCULATE INTERACTIONS *
1062 **************************/
1063
1064 if (gmx_mm_any_lt(rsq10,rcutoff2))
1065 {
1066
1067 r10 = _mm_mul_ps(rsq10,rinv10);
1068
1069 /* Compute parameters for interactions between i and j atoms */
1070 qq10 = _mm_mul_ps(iq1,jq0);
1071
1072 /* EWALD ELECTROSTATICS */
1073
1074 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1075 ewrt = _mm_mul_ps(r10,ewtabscale);
1076 ewitab = _mm_cvttps_epi32(ewrt);
1077 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
1078 gmx_mm_load_4pair_swizzle_ps(ewtab+gmx_mm_extract_epi32(ewitab,0),ewtab+gmx_mm_extract_epi32(ewitab,1),
1079 ewtab+gmx_mm_extract_epi32(ewitab,2),ewtab+gmx_mm_extract_epi32(ewitab,3),
1080 &ewtabF,&ewtabFn);
1081 felec = _mm_add_ps(_mm_mul_ps( _mm_sub_ps(one,eweps),ewtabF),_mm_mul_ps(eweps,ewtabFn));
1082 felec = _mm_mul_ps(_mm_mul_ps(qq10,rinv10),_mm_sub_ps(rinvsq10,felec));
1083
1084 cutoff_mask = _mm_cmplt_ps(rsq10,rcutoff2);
1085
1086 fscal = felec;
1087
1088 fscal = _mm_and_ps(fscal,cutoff_mask);
1089
1090 /* Calculate temporary vectorial force */
1091 tx = _mm_mul_ps(fscal,dx10);
1092 ty = _mm_mul_ps(fscal,dy10);
1093 tz = _mm_mul_ps(fscal,dz10);
1094
1095 /* Update vectorial force */
1096 fix1 = _mm_add_ps(fix1,tx);
1097 fiy1 = _mm_add_ps(fiy1,ty);
1098 fiz1 = _mm_add_ps(fiz1,tz);
1099
1100 fjx0 = _mm_add_ps(fjx0,tx);
1101 fjy0 = _mm_add_ps(fjy0,ty);
1102 fjz0 = _mm_add_ps(fjz0,tz);
1103
1104 }
1105
1106 /**************************
1107 * CALCULATE INTERACTIONS *
1108 **************************/
1109
1110 if (gmx_mm_any_lt(rsq20,rcutoff2))
1111 {
1112
1113 r20 = _mm_mul_ps(rsq20,rinv20);
1114
1115 /* Compute parameters for interactions between i and j atoms */
1116 qq20 = _mm_mul_ps(iq2,jq0);
1117
1118 /* EWALD ELECTROSTATICS */
1119
1120 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1121 ewrt = _mm_mul_ps(r20,ewtabscale);
1122 ewitab = _mm_cvttps_epi32(ewrt);
1123 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
1124 gmx_mm_load_4pair_swizzle_ps(ewtab+gmx_mm_extract_epi32(ewitab,0),ewtab+gmx_mm_extract_epi32(ewitab,1),
1125 ewtab+gmx_mm_extract_epi32(ewitab,2),ewtab+gmx_mm_extract_epi32(ewitab,3),
1126 &ewtabF,&ewtabFn);
1127 felec = _mm_add_ps(_mm_mul_ps( _mm_sub_ps(one,eweps),ewtabF),_mm_mul_ps(eweps,ewtabFn));
1128 felec = _mm_mul_ps(_mm_mul_ps(qq20,rinv20),_mm_sub_ps(rinvsq20,felec));
1129
1130 cutoff_mask = _mm_cmplt_ps(rsq20,rcutoff2);
1131
1132 fscal = felec;
1133
1134 fscal = _mm_and_ps(fscal,cutoff_mask);
1135
1136 /* Calculate temporary vectorial force */
1137 tx = _mm_mul_ps(fscal,dx20);
1138 ty = _mm_mul_ps(fscal,dy20);
1139 tz = _mm_mul_ps(fscal,dz20);
1140
1141 /* Update vectorial force */
1142 fix2 = _mm_add_ps(fix2,tx);
1143 fiy2 = _mm_add_ps(fiy2,ty);
1144 fiz2 = _mm_add_ps(fiz2,tz);
1145
1146 fjx0 = _mm_add_ps(fjx0,tx);
1147 fjy0 = _mm_add_ps(fjy0,ty);
1148 fjz0 = _mm_add_ps(fjz0,tz);
1149
1150 }
1151
1152 /**************************
1153 * CALCULATE INTERACTIONS *
1154 **************************/
1155
1156 if (gmx_mm_any_lt(rsq30,rcutoff2))
1157 {
1158
1159 r30 = _mm_mul_ps(rsq30,rinv30);
1160
1161 /* Compute parameters for interactions between i and j atoms */
1162 qq30 = _mm_mul_ps(iq3,jq0);
1163
1164 /* EWALD ELECTROSTATICS */
1165
1166 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1167 ewrt = _mm_mul_ps(r30,ewtabscale);
1168 ewitab = _mm_cvttps_epi32(ewrt);
1169 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
1170 gmx_mm_load_4pair_swizzle_ps(ewtab+gmx_mm_extract_epi32(ewitab,0),ewtab+gmx_mm_extract_epi32(ewitab,1),
1171 ewtab+gmx_mm_extract_epi32(ewitab,2),ewtab+gmx_mm_extract_epi32(ewitab,3),
1172 &ewtabF,&ewtabFn);
1173 felec = _mm_add_ps(_mm_mul_ps( _mm_sub_ps(one,eweps),ewtabF),_mm_mul_ps(eweps,ewtabFn));
1174 felec = _mm_mul_ps(_mm_mul_ps(qq30,rinv30),_mm_sub_ps(rinvsq30,felec));
1175
1176 cutoff_mask = _mm_cmplt_ps(rsq30,rcutoff2);
1177
1178 fscal = felec;
1179
1180 fscal = _mm_and_ps(fscal,cutoff_mask);
1181
1182 /* Calculate temporary vectorial force */
1183 tx = _mm_mul_ps(fscal,dx30);
1184 ty = _mm_mul_ps(fscal,dy30);
1185 tz = _mm_mul_ps(fscal,dz30);
1186
1187 /* Update vectorial force */
1188 fix3 = _mm_add_ps(fix3,tx);
1189 fiy3 = _mm_add_ps(fiy3,ty);
1190 fiz3 = _mm_add_ps(fiz3,tz);
1191
1192 fjx0 = _mm_add_ps(fjx0,tx);
1193 fjy0 = _mm_add_ps(fjy0,ty);
1194 fjz0 = _mm_add_ps(fjz0,tz);
1195
1196 }
1197
1198 fjptrA = f+j_coord_offsetA;
1199 fjptrB = f+j_coord_offsetB;
1200 fjptrC = f+j_coord_offsetC;
1201 fjptrD = f+j_coord_offsetD;
1202
1203 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
1204
1205 /* Inner loop uses 147 flops */
1206 }
1207
1208 if(jidx<j_index_end)
1209 {
1210
1211 /* Get j neighbor index, and coordinate index */
1212 jnrlistA = jjnr[jidx];
1213 jnrlistB = jjnr[jidx+1];
1214 jnrlistC = jjnr[jidx+2];
1215 jnrlistD = jjnr[jidx+3];
1216 /* Sign of each element will be negative for non-real atoms.
1217 * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
1218 * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
1219 */
1220 dummy_mask = gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
1221 jnrA = (jnrlistA>=0) ? jnrlistA : 0;
1222 jnrB = (jnrlistB>=0) ? jnrlistB : 0;
1223 jnrC = (jnrlistC>=0) ? jnrlistC : 0;
1224 jnrD = (jnrlistD>=0) ? jnrlistD : 0;
1225 j_coord_offsetA = DIM*jnrA;
1226 j_coord_offsetB = DIM*jnrB;
1227 j_coord_offsetC = DIM*jnrC;
1228 j_coord_offsetD = DIM*jnrD;
1229
1230 /* load j atom coordinates */
1231 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
1232 x+j_coord_offsetC,x+j_coord_offsetD,
1233 &jx0,&jy0,&jz0);
1234
1235 /* Calculate displacement vector */
1236 dx00 = _mm_sub_ps(ix0,jx0);
1237 dy00 = _mm_sub_ps(iy0,jy0);
1238 dz00 = _mm_sub_ps(iz0,jz0);
1239 dx10 = _mm_sub_ps(ix1,jx0);
1240 dy10 = _mm_sub_ps(iy1,jy0);
1241 dz10 = _mm_sub_ps(iz1,jz0);
1242 dx20 = _mm_sub_ps(ix2,jx0);
1243 dy20 = _mm_sub_ps(iy2,jy0);
1244 dz20 = _mm_sub_ps(iz2,jz0);
1245 dx30 = _mm_sub_ps(ix3,jx0);
1246 dy30 = _mm_sub_ps(iy3,jy0);
1247 dz30 = _mm_sub_ps(iz3,jz0);
1248
1249 /* Calculate squared distance and things based on it */
1250 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
1251 rsq10 = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
1252 rsq20 = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
1253 rsq30 = gmx_mm_calc_rsq_ps(dx30,dy30,dz30);
1254
1255 rinv10 = gmx_mm_invsqrt_ps(rsq10);
1256 rinv20 = gmx_mm_invsqrt_ps(rsq20);
1257 rinv30 = gmx_mm_invsqrt_ps(rsq30);
1258
1259 rinvsq00 = gmx_mm_inv_ps(rsq00);
1260 rinvsq10 = _mm_mul_ps(rinv10,rinv10);
1261 rinvsq20 = _mm_mul_ps(rinv20,rinv20);
1262 rinvsq30 = _mm_mul_ps(rinv30,rinv30);
1263
1264 /* Load parameters for j particles */
1265 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
1266 charge+jnrC+0,charge+jnrD+0);
1267 vdwjidx0A = 2*vdwtype[jnrA+0];
1268 vdwjidx0B = 2*vdwtype[jnrB+0];
1269 vdwjidx0C = 2*vdwtype[jnrC+0];
1270 vdwjidx0D = 2*vdwtype[jnrD+0];
1271
1272 fjx0 = _mm_setzero_ps();
1273 fjy0 = _mm_setzero_ps();
1274 fjz0 = _mm_setzero_ps();
1275
1276 /**************************
1277 * CALCULATE INTERACTIONS *
1278 **************************/
1279
1280 if (gmx_mm_any_lt(rsq00,rcutoff2))
1281 {
1282
1283 /* Compute parameters for interactions between i and j atoms */
1284 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
1285 vdwparam+vdwioffset0+vdwjidx0B,
1286 vdwparam+vdwioffset0+vdwjidx0C,
1287 vdwparam+vdwioffset0+vdwjidx0D,
1288 &c6_00,&c12_00);
1289
1290 /* LENNARD-JONES DISPERSION/REPULSION */
1291
1292 rinvsix = _mm_mul_ps(_mm_mul_ps(rinvsq00,rinvsq00),rinvsq00);
1293 fvdw = _mm_mul_ps(_mm_sub_ps(_mm_mul_ps(c12_00,rinvsix),c6_00),_mm_mul_ps(rinvsix,rinvsq00));
1294
1295 cutoff_mask = _mm_cmplt_ps(rsq00,rcutoff2);
1296
1297 fscal = fvdw;
1298
1299 fscal = _mm_and_ps(fscal,cutoff_mask);
1300
1301 fscal = _mm_andnot_ps(dummy_mask,fscal);
1302
1303 /* Calculate temporary vectorial force */
1304 tx = _mm_mul_ps(fscal,dx00);
1305 ty = _mm_mul_ps(fscal,dy00);
1306 tz = _mm_mul_ps(fscal,dz00);
1307
1308 /* Update vectorial force */
1309 fix0 = _mm_add_ps(fix0,tx);
1310 fiy0 = _mm_add_ps(fiy0,ty);
1311 fiz0 = _mm_add_ps(fiz0,tz);
1312
1313 fjx0 = _mm_add_ps(fjx0,tx);
1314 fjy0 = _mm_add_ps(fjy0,ty);
1315 fjz0 = _mm_add_ps(fjz0,tz);
1316
1317 }
1318
1319 /**************************
1320 * CALCULATE INTERACTIONS *
1321 **************************/
1322
1323 if (gmx_mm_any_lt(rsq10,rcutoff2))
1324 {
1325
1326 r10 = _mm_mul_ps(rsq10,rinv10);
1327 r10 = _mm_andnot_ps(dummy_mask,r10);
1328
1329 /* Compute parameters for interactions between i and j atoms */
1330 qq10 = _mm_mul_ps(iq1,jq0);
1331
1332 /* EWALD ELECTROSTATICS */
1333
1334 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1335 ewrt = _mm_mul_ps(r10,ewtabscale);
1336 ewitab = _mm_cvttps_epi32(ewrt);
1337 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
1338 gmx_mm_load_4pair_swizzle_ps(ewtab+gmx_mm_extract_epi32(ewitab,0),ewtab+gmx_mm_extract_epi32(ewitab,1),
1339 ewtab+gmx_mm_extract_epi32(ewitab,2),ewtab+gmx_mm_extract_epi32(ewitab,3),
1340 &ewtabF,&ewtabFn);
1341 felec = _mm_add_ps(_mm_mul_ps( _mm_sub_ps(one,eweps),ewtabF),_mm_mul_ps(eweps,ewtabFn));
1342 felec = _mm_mul_ps(_mm_mul_ps(qq10,rinv10),_mm_sub_ps(rinvsq10,felec));
1343
1344 cutoff_mask = _mm_cmplt_ps(rsq10,rcutoff2);
1345
1346 fscal = felec;
1347
1348 fscal = _mm_and_ps(fscal,cutoff_mask);
1349
1350 fscal = _mm_andnot_ps(dummy_mask,fscal);
1351
1352 /* Calculate temporary vectorial force */
1353 tx = _mm_mul_ps(fscal,dx10);
1354 ty = _mm_mul_ps(fscal,dy10);
1355 tz = _mm_mul_ps(fscal,dz10);
1356
1357 /* Update vectorial force */
1358 fix1 = _mm_add_ps(fix1,tx);
1359 fiy1 = _mm_add_ps(fiy1,ty);
1360 fiz1 = _mm_add_ps(fiz1,tz);
1361
1362 fjx0 = _mm_add_ps(fjx0,tx);
1363 fjy0 = _mm_add_ps(fjy0,ty);
1364 fjz0 = _mm_add_ps(fjz0,tz);
1365
1366 }
1367
1368 /**************************
1369 * CALCULATE INTERACTIONS *
1370 **************************/
1371
1372 if (gmx_mm_any_lt(rsq20,rcutoff2))
1373 {
1374
1375 r20 = _mm_mul_ps(rsq20,rinv20);
1376 r20 = _mm_andnot_ps(dummy_mask,r20);
1377
1378 /* Compute parameters for interactions between i and j atoms */
1379 qq20 = _mm_mul_ps(iq2,jq0);
1380
1381 /* EWALD ELECTROSTATICS */
1382
1383 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1384 ewrt = _mm_mul_ps(r20,ewtabscale);
1385 ewitab = _mm_cvttps_epi32(ewrt);
1386 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
1387 gmx_mm_load_4pair_swizzle_ps(ewtab+gmx_mm_extract_epi32(ewitab,0),ewtab+gmx_mm_extract_epi32(ewitab,1),
1388 ewtab+gmx_mm_extract_epi32(ewitab,2),ewtab+gmx_mm_extract_epi32(ewitab,3),
1389 &ewtabF,&ewtabFn);
1390 felec = _mm_add_ps(_mm_mul_ps( _mm_sub_ps(one,eweps),ewtabF),_mm_mul_ps(eweps,ewtabFn));
1391 felec = _mm_mul_ps(_mm_mul_ps(qq20,rinv20),_mm_sub_ps(rinvsq20,felec));
1392
1393 cutoff_mask = _mm_cmplt_ps(rsq20,rcutoff2);
1394
1395 fscal = felec;
1396
1397 fscal = _mm_and_ps(fscal,cutoff_mask);
1398
1399 fscal = _mm_andnot_ps(dummy_mask,fscal);
1400
1401 /* Calculate temporary vectorial force */
1402 tx = _mm_mul_ps(fscal,dx20);
1403 ty = _mm_mul_ps(fscal,dy20);
1404 tz = _mm_mul_ps(fscal,dz20);
1405
1406 /* Update vectorial force */
1407 fix2 = _mm_add_ps(fix2,tx);
1408 fiy2 = _mm_add_ps(fiy2,ty);
1409 fiz2 = _mm_add_ps(fiz2,tz);
1410
1411 fjx0 = _mm_add_ps(fjx0,tx);
1412 fjy0 = _mm_add_ps(fjy0,ty);
1413 fjz0 = _mm_add_ps(fjz0,tz);
1414
1415 }
1416
1417 /**************************
1418 * CALCULATE INTERACTIONS *
1419 **************************/
1420
1421 if (gmx_mm_any_lt(rsq30,rcutoff2))
1422 {
1423
1424 r30 = _mm_mul_ps(rsq30,rinv30);
1425 r30 = _mm_andnot_ps(dummy_mask,r30);
1426
1427 /* Compute parameters for interactions between i and j atoms */
1428 qq30 = _mm_mul_ps(iq3,jq0);
1429
1430 /* EWALD ELECTROSTATICS */
1431
1432 /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1433 ewrt = _mm_mul_ps(r30,ewtabscale);
1434 ewitab = _mm_cvttps_epi32(ewrt);
1435 eweps = _mm_sub_ps(ewrt,_mm_cvtepi32_ps(ewitab));
1436 gmx_mm_load_4pair_swizzle_ps(ewtab+gmx_mm_extract_epi32(ewitab,0),ewtab+gmx_mm_extract_epi32(ewitab,1),
1437 ewtab+gmx_mm_extract_epi32(ewitab,2),ewtab+gmx_mm_extract_epi32(ewitab,3),
1438 &ewtabF,&ewtabFn);
1439 felec = _mm_add_ps(_mm_mul_ps( _mm_sub_ps(one,eweps),ewtabF),_mm_mul_ps(eweps,ewtabFn));
1440 felec = _mm_mul_ps(_mm_mul_ps(qq30,rinv30),_mm_sub_ps(rinvsq30,felec));
1441
1442 cutoff_mask = _mm_cmplt_ps(rsq30,rcutoff2);
1443
1444 fscal = felec;
1445
1446 fscal = _mm_and_ps(fscal,cutoff_mask);
1447
1448 fscal = _mm_andnot_ps(dummy_mask,fscal);
1449
1450 /* Calculate temporary vectorial force */
1451 tx = _mm_mul_ps(fscal,dx30);
1452 ty = _mm_mul_ps(fscal,dy30);
1453 tz = _mm_mul_ps(fscal,dz30);
1454
1455 /* Update vectorial force */
1456 fix3 = _mm_add_ps(fix3,tx);
1457 fiy3 = _mm_add_ps(fiy3,ty);
1458 fiz3 = _mm_add_ps(fiz3,tz);
1459
1460 fjx0 = _mm_add_ps(fjx0,tx);
1461 fjy0 = _mm_add_ps(fjy0,ty);
1462 fjz0 = _mm_add_ps(fjz0,tz);
1463
1464 }
1465
1466 fjptrA = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
1467 fjptrB = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
1468 fjptrC = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
1469 fjptrD = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
1470
1471 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
1472
1473 /* Inner loop uses 150 flops */
1474 }
1475
1476 /* End of innermost loop */
1477
1478 gmx_mm_update_iforce_4atom_swizzle_ps(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
1479 f+i_coord_offset,fshift+i_shift_offset);
1480
1481 /* Increment number of inner iterations */
1482 inneriter += j_index_end - j_index_start;
1483
1484 /* Outer loop uses 24 flops */
1485 }
1486
1487 /* Increment number of outer iterations */
1488 outeriter += nri;
1489
1490 /* Update outer/inner flops */
1491
1492 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4_F,outeriter*24 + inneriter*150);
1493 }
The ultimate molecular dynamics simulation package