2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014,2015,2017, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
36 * Note: this file was generated by the GROMACS sse4_1_single kernel generator.
44 #include "../nb_kernel.h"
45 #include "gromacs/gmxlib/nrnb.h"
47 #include "kernelutil_x86_sse4_1_single.h"
50 * Gromacs nonbonded kernel: nb_kernel_ElecCoul_VdwCSTab_GeomP1P1_VF_sse4_1_single
51 * Electrostatics interaction: Coulomb
52 * VdW interaction: CubicSplineTable
53 * Geometry: Particle-Particle
54 * Calculate force/pot: PotentialAndForce
57 nb_kernel_ElecCoul_VdwCSTab_GeomP1P1_VF_sse4_1_single
58 (t_nblist * gmx_restrict nlist,
59 rvec * gmx_restrict xx,
60 rvec * gmx_restrict ff,
61 struct t_forcerec * gmx_restrict fr,
62 t_mdatoms * gmx_restrict mdatoms,
63 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
64 t_nrnb * gmx_restrict nrnb)
66 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
67 * just 0 for non-waters.
68 * Suffixes A,B,C,D refer to j loop unrolling done with SSE, e.g. for the four different
69 * jnr indices corresponding to data put in the four positions in the SIMD register.
71 int i_shift_offset,i_coord_offset,outeriter,inneriter;
72 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
73 int jnrA,jnrB,jnrC,jnrD;
74 int jnrlistA,jnrlistB,jnrlistC,jnrlistD;
75 int j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
76 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
78 real *shiftvec,*fshift,*x,*f;
79 real *fjptrA,*fjptrB,*fjptrC,*fjptrD;
81 __m128 tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
83 __m128 ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
84 int vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
85 __m128 jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
86 __m128 dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
87 __m128 velec,felec,velecsum,facel,crf,krf,krf2;
90 __m128 rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
93 __m128 one_sixth = _mm_set1_ps(1.0/6.0);
94 __m128 one_twelfth = _mm_set1_ps(1.0/12.0);
96 __m128i ifour = _mm_set1_epi32(4);
97 __m128 rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF;
99 __m128 dummy_mask,cutoff_mask;
100 __m128 signbit = _mm_castsi128_ps( _mm_set1_epi32(0x80000000) );
101 __m128 one = _mm_set1_ps(1.0);
102 __m128 two = _mm_set1_ps(2.0);
108 jindex = nlist->jindex;
110 shiftidx = nlist->shift;
112 shiftvec = fr->shift_vec[0];
113 fshift = fr->fshift[0];
114 facel = _mm_set1_ps(fr->ic->epsfac);
115 charge = mdatoms->chargeA;
116 nvdwtype = fr->ntype;
118 vdwtype = mdatoms->typeA;
120 vftab = kernel_data->table_vdw->data;
121 vftabscale = _mm_set1_ps(kernel_data->table_vdw->scale);
123 /* Avoid stupid compiler warnings */
124 jnrA = jnrB = jnrC = jnrD = 0;
133 for(iidx=0;iidx<4*DIM;iidx++)
138 /* Start outer loop over neighborlists */
139 for(iidx=0; iidx<nri; iidx++)
141 /* Load shift vector for this list */
142 i_shift_offset = DIM*shiftidx[iidx];
144 /* Load limits for loop over neighbors */
145 j_index_start = jindex[iidx];
146 j_index_end = jindex[iidx+1];
148 /* Get outer coordinate index */
150 i_coord_offset = DIM*inr;
152 /* Load i particle coords and add shift vector */
153 gmx_mm_load_shift_and_1rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,&ix0,&iy0,&iz0);
155 fix0 = _mm_setzero_ps();
156 fiy0 = _mm_setzero_ps();
157 fiz0 = _mm_setzero_ps();
159 /* Load parameters for i particles */
160 iq0 = _mm_mul_ps(facel,_mm_load1_ps(charge+inr+0));
161 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
163 /* Reset potential sums */
164 velecsum = _mm_setzero_ps();
165 vvdwsum = _mm_setzero_ps();
167 /* Start inner kernel loop */
168 for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
171 /* Get j neighbor index, and coordinate index */
176 j_coord_offsetA = DIM*jnrA;
177 j_coord_offsetB = DIM*jnrB;
178 j_coord_offsetC = DIM*jnrC;
179 j_coord_offsetD = DIM*jnrD;
181 /* load j atom coordinates */
182 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
183 x+j_coord_offsetC,x+j_coord_offsetD,
186 /* Calculate displacement vector */
187 dx00 = _mm_sub_ps(ix0,jx0);
188 dy00 = _mm_sub_ps(iy0,jy0);
189 dz00 = _mm_sub_ps(iz0,jz0);
191 /* Calculate squared distance and things based on it */
192 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
194 rinv00 = sse41_invsqrt_f(rsq00);
196 rinvsq00 = _mm_mul_ps(rinv00,rinv00);
198 /* Load parameters for j particles */
199 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
200 charge+jnrC+0,charge+jnrD+0);
201 vdwjidx0A = 2*vdwtype[jnrA+0];
202 vdwjidx0B = 2*vdwtype[jnrB+0];
203 vdwjidx0C = 2*vdwtype[jnrC+0];
204 vdwjidx0D = 2*vdwtype[jnrD+0];
206 /**************************
207 * CALCULATE INTERACTIONS *
208 **************************/
210 r00 = _mm_mul_ps(rsq00,rinv00);
212 /* Compute parameters for interactions between i and j atoms */
213 qq00 = _mm_mul_ps(iq0,jq0);
214 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
215 vdwparam+vdwioffset0+vdwjidx0B,
216 vdwparam+vdwioffset0+vdwjidx0C,
217 vdwparam+vdwioffset0+vdwjidx0D,
220 /* Calculate table index by multiplying r with table scale and truncate to integer */
221 rt = _mm_mul_ps(r00,vftabscale);
222 vfitab = _mm_cvttps_epi32(rt);
223 vfeps = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
224 vfitab = _mm_slli_epi32(vfitab,3);
226 /* COULOMB ELECTROSTATICS */
227 velec = _mm_mul_ps(qq00,rinv00);
228 felec = _mm_mul_ps(velec,rinvsq00);
230 /* CUBIC SPLINE TABLE DISPERSION */
231 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
232 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
233 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
234 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
235 _MM_TRANSPOSE4_PS(Y,F,G,H);
236 Heps = _mm_mul_ps(vfeps,H);
237 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
238 VV = _mm_add_ps(Y,_mm_mul_ps(vfeps,Fp));
239 vvdw6 = _mm_mul_ps(c6_00,VV);
240 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
241 fvdw6 = _mm_mul_ps(c6_00,FF);
243 /* CUBIC SPLINE TABLE REPULSION */
244 vfitab = _mm_add_epi32(vfitab,ifour);
245 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
246 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
247 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
248 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
249 _MM_TRANSPOSE4_PS(Y,F,G,H);
250 Heps = _mm_mul_ps(vfeps,H);
251 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
252 VV = _mm_add_ps(Y,_mm_mul_ps(vfeps,Fp));
253 vvdw12 = _mm_mul_ps(c12_00,VV);
254 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
255 fvdw12 = _mm_mul_ps(c12_00,FF);
256 vvdw = _mm_add_ps(vvdw12,vvdw6);
257 fvdw = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
259 /* Update potential sum for this i atom from the interaction with this j atom. */
260 velecsum = _mm_add_ps(velecsum,velec);
261 vvdwsum = _mm_add_ps(vvdwsum,vvdw);
263 fscal = _mm_add_ps(felec,fvdw);
265 /* Calculate temporary vectorial force */
266 tx = _mm_mul_ps(fscal,dx00);
267 ty = _mm_mul_ps(fscal,dy00);
268 tz = _mm_mul_ps(fscal,dz00);
270 /* Update vectorial force */
271 fix0 = _mm_add_ps(fix0,tx);
272 fiy0 = _mm_add_ps(fiy0,ty);
273 fiz0 = _mm_add_ps(fiz0,tz);
275 fjptrA = f+j_coord_offsetA;
276 fjptrB = f+j_coord_offsetB;
277 fjptrC = f+j_coord_offsetC;
278 fjptrD = f+j_coord_offsetD;
279 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,tx,ty,tz);
281 /* Inner loop uses 63 flops */
287 /* Get j neighbor index, and coordinate index */
288 jnrlistA = jjnr[jidx];
289 jnrlistB = jjnr[jidx+1];
290 jnrlistC = jjnr[jidx+2];
291 jnrlistD = jjnr[jidx+3];
292 /* Sign of each element will be negative for non-real atoms.
293 * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
294 * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
296 dummy_mask = gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
297 jnrA = (jnrlistA>=0) ? jnrlistA : 0;
298 jnrB = (jnrlistB>=0) ? jnrlistB : 0;
299 jnrC = (jnrlistC>=0) ? jnrlistC : 0;
300 jnrD = (jnrlistD>=0) ? jnrlistD : 0;
301 j_coord_offsetA = DIM*jnrA;
302 j_coord_offsetB = DIM*jnrB;
303 j_coord_offsetC = DIM*jnrC;
304 j_coord_offsetD = DIM*jnrD;
306 /* load j atom coordinates */
307 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
308 x+j_coord_offsetC,x+j_coord_offsetD,
311 /* Calculate displacement vector */
312 dx00 = _mm_sub_ps(ix0,jx0);
313 dy00 = _mm_sub_ps(iy0,jy0);
314 dz00 = _mm_sub_ps(iz0,jz0);
316 /* Calculate squared distance and things based on it */
317 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
319 rinv00 = sse41_invsqrt_f(rsq00);
321 rinvsq00 = _mm_mul_ps(rinv00,rinv00);
323 /* Load parameters for j particles */
324 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
325 charge+jnrC+0,charge+jnrD+0);
326 vdwjidx0A = 2*vdwtype[jnrA+0];
327 vdwjidx0B = 2*vdwtype[jnrB+0];
328 vdwjidx0C = 2*vdwtype[jnrC+0];
329 vdwjidx0D = 2*vdwtype[jnrD+0];
331 /**************************
332 * CALCULATE INTERACTIONS *
333 **************************/
335 r00 = _mm_mul_ps(rsq00,rinv00);
336 r00 = _mm_andnot_ps(dummy_mask,r00);
338 /* Compute parameters for interactions between i and j atoms */
339 qq00 = _mm_mul_ps(iq0,jq0);
340 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
341 vdwparam+vdwioffset0+vdwjidx0B,
342 vdwparam+vdwioffset0+vdwjidx0C,
343 vdwparam+vdwioffset0+vdwjidx0D,
346 /* Calculate table index by multiplying r with table scale and truncate to integer */
347 rt = _mm_mul_ps(r00,vftabscale);
348 vfitab = _mm_cvttps_epi32(rt);
349 vfeps = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
350 vfitab = _mm_slli_epi32(vfitab,3);
352 /* COULOMB ELECTROSTATICS */
353 velec = _mm_mul_ps(qq00,rinv00);
354 felec = _mm_mul_ps(velec,rinvsq00);
356 /* CUBIC SPLINE TABLE DISPERSION */
357 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
358 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
359 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
360 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
361 _MM_TRANSPOSE4_PS(Y,F,G,H);
362 Heps = _mm_mul_ps(vfeps,H);
363 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
364 VV = _mm_add_ps(Y,_mm_mul_ps(vfeps,Fp));
365 vvdw6 = _mm_mul_ps(c6_00,VV);
366 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
367 fvdw6 = _mm_mul_ps(c6_00,FF);
369 /* CUBIC SPLINE TABLE REPULSION */
370 vfitab = _mm_add_epi32(vfitab,ifour);
371 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
372 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
373 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
374 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
375 _MM_TRANSPOSE4_PS(Y,F,G,H);
376 Heps = _mm_mul_ps(vfeps,H);
377 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
378 VV = _mm_add_ps(Y,_mm_mul_ps(vfeps,Fp));
379 vvdw12 = _mm_mul_ps(c12_00,VV);
380 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
381 fvdw12 = _mm_mul_ps(c12_00,FF);
382 vvdw = _mm_add_ps(vvdw12,vvdw6);
383 fvdw = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
385 /* Update potential sum for this i atom from the interaction with this j atom. */
386 velec = _mm_andnot_ps(dummy_mask,velec);
387 velecsum = _mm_add_ps(velecsum,velec);
388 vvdw = _mm_andnot_ps(dummy_mask,vvdw);
389 vvdwsum = _mm_add_ps(vvdwsum,vvdw);
391 fscal = _mm_add_ps(felec,fvdw);
393 fscal = _mm_andnot_ps(dummy_mask,fscal);
395 /* Calculate temporary vectorial force */
396 tx = _mm_mul_ps(fscal,dx00);
397 ty = _mm_mul_ps(fscal,dy00);
398 tz = _mm_mul_ps(fscal,dz00);
400 /* Update vectorial force */
401 fix0 = _mm_add_ps(fix0,tx);
402 fiy0 = _mm_add_ps(fiy0,ty);
403 fiz0 = _mm_add_ps(fiz0,tz);
405 fjptrA = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
406 fjptrB = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
407 fjptrC = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
408 fjptrD = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
409 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,tx,ty,tz);
411 /* Inner loop uses 64 flops */
414 /* End of innermost loop */
416 gmx_mm_update_iforce_1atom_swizzle_ps(fix0,fiy0,fiz0,
417 f+i_coord_offset,fshift+i_shift_offset);
420 /* Update potential energies */
421 gmx_mm_update_1pot_ps(velecsum,kernel_data->energygrp_elec+ggid);
422 gmx_mm_update_1pot_ps(vvdwsum,kernel_data->energygrp_vdw+ggid);
424 /* Increment number of inner iterations */
425 inneriter += j_index_end - j_index_start;
427 /* Outer loop uses 9 flops */
430 /* Increment number of outer iterations */
433 /* Update outer/inner flops */
435 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_VF,outeriter*9 + inneriter*64);
438 * Gromacs nonbonded kernel: nb_kernel_ElecCoul_VdwCSTab_GeomP1P1_F_sse4_1_single
439 * Electrostatics interaction: Coulomb
440 * VdW interaction: CubicSplineTable
441 * Geometry: Particle-Particle
442 * Calculate force/pot: Force
445 nb_kernel_ElecCoul_VdwCSTab_GeomP1P1_F_sse4_1_single
446 (t_nblist * gmx_restrict nlist,
447 rvec * gmx_restrict xx,
448 rvec * gmx_restrict ff,
449 struct t_forcerec * gmx_restrict fr,
450 t_mdatoms * gmx_restrict mdatoms,
451 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
452 t_nrnb * gmx_restrict nrnb)
454 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
455 * just 0 for non-waters.
456 * Suffixes A,B,C,D refer to j loop unrolling done with SSE, e.g. for the four different
457 * jnr indices corresponding to data put in the four positions in the SIMD register.
459 int i_shift_offset,i_coord_offset,outeriter,inneriter;
460 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
461 int jnrA,jnrB,jnrC,jnrD;
462 int jnrlistA,jnrlistB,jnrlistC,jnrlistD;
463 int j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
464 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
466 real *shiftvec,*fshift,*x,*f;
467 real *fjptrA,*fjptrB,*fjptrC,*fjptrD;
469 __m128 tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
471 __m128 ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
472 int vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
473 __m128 jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
474 __m128 dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
475 __m128 velec,felec,velecsum,facel,crf,krf,krf2;
478 __m128 rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
481 __m128 one_sixth = _mm_set1_ps(1.0/6.0);
482 __m128 one_twelfth = _mm_set1_ps(1.0/12.0);
484 __m128i ifour = _mm_set1_epi32(4);
485 __m128 rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF;
487 __m128 dummy_mask,cutoff_mask;
488 __m128 signbit = _mm_castsi128_ps( _mm_set1_epi32(0x80000000) );
489 __m128 one = _mm_set1_ps(1.0);
490 __m128 two = _mm_set1_ps(2.0);
496 jindex = nlist->jindex;
498 shiftidx = nlist->shift;
500 shiftvec = fr->shift_vec[0];
501 fshift = fr->fshift[0];
502 facel = _mm_set1_ps(fr->ic->epsfac);
503 charge = mdatoms->chargeA;
504 nvdwtype = fr->ntype;
506 vdwtype = mdatoms->typeA;
508 vftab = kernel_data->table_vdw->data;
509 vftabscale = _mm_set1_ps(kernel_data->table_vdw->scale);
511 /* Avoid stupid compiler warnings */
512 jnrA = jnrB = jnrC = jnrD = 0;
521 for(iidx=0;iidx<4*DIM;iidx++)
526 /* Start outer loop over neighborlists */
527 for(iidx=0; iidx<nri; iidx++)
529 /* Load shift vector for this list */
530 i_shift_offset = DIM*shiftidx[iidx];
532 /* Load limits for loop over neighbors */
533 j_index_start = jindex[iidx];
534 j_index_end = jindex[iidx+1];
536 /* Get outer coordinate index */
538 i_coord_offset = DIM*inr;
540 /* Load i particle coords and add shift vector */
541 gmx_mm_load_shift_and_1rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,&ix0,&iy0,&iz0);
543 fix0 = _mm_setzero_ps();
544 fiy0 = _mm_setzero_ps();
545 fiz0 = _mm_setzero_ps();
547 /* Load parameters for i particles */
548 iq0 = _mm_mul_ps(facel,_mm_load1_ps(charge+inr+0));
549 vdwioffset0 = 2*nvdwtype*vdwtype[inr+0];
551 /* Start inner kernel loop */
552 for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
555 /* Get j neighbor index, and coordinate index */
560 j_coord_offsetA = DIM*jnrA;
561 j_coord_offsetB = DIM*jnrB;
562 j_coord_offsetC = DIM*jnrC;
563 j_coord_offsetD = DIM*jnrD;
565 /* load j atom coordinates */
566 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
567 x+j_coord_offsetC,x+j_coord_offsetD,
570 /* Calculate displacement vector */
571 dx00 = _mm_sub_ps(ix0,jx0);
572 dy00 = _mm_sub_ps(iy0,jy0);
573 dz00 = _mm_sub_ps(iz0,jz0);
575 /* Calculate squared distance and things based on it */
576 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
578 rinv00 = sse41_invsqrt_f(rsq00);
580 rinvsq00 = _mm_mul_ps(rinv00,rinv00);
582 /* Load parameters for j particles */
583 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
584 charge+jnrC+0,charge+jnrD+0);
585 vdwjidx0A = 2*vdwtype[jnrA+0];
586 vdwjidx0B = 2*vdwtype[jnrB+0];
587 vdwjidx0C = 2*vdwtype[jnrC+0];
588 vdwjidx0D = 2*vdwtype[jnrD+0];
590 /**************************
591 * CALCULATE INTERACTIONS *
592 **************************/
594 r00 = _mm_mul_ps(rsq00,rinv00);
596 /* Compute parameters for interactions between i and j atoms */
597 qq00 = _mm_mul_ps(iq0,jq0);
598 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
599 vdwparam+vdwioffset0+vdwjidx0B,
600 vdwparam+vdwioffset0+vdwjidx0C,
601 vdwparam+vdwioffset0+vdwjidx0D,
604 /* Calculate table index by multiplying r with table scale and truncate to integer */
605 rt = _mm_mul_ps(r00,vftabscale);
606 vfitab = _mm_cvttps_epi32(rt);
607 vfeps = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
608 vfitab = _mm_slli_epi32(vfitab,3);
610 /* COULOMB ELECTROSTATICS */
611 velec = _mm_mul_ps(qq00,rinv00);
612 felec = _mm_mul_ps(velec,rinvsq00);
614 /* CUBIC SPLINE TABLE DISPERSION */
615 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
616 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
617 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
618 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
619 _MM_TRANSPOSE4_PS(Y,F,G,H);
620 Heps = _mm_mul_ps(vfeps,H);
621 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
622 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
623 fvdw6 = _mm_mul_ps(c6_00,FF);
625 /* CUBIC SPLINE TABLE REPULSION */
626 vfitab = _mm_add_epi32(vfitab,ifour);
627 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
628 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
629 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
630 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
631 _MM_TRANSPOSE4_PS(Y,F,G,H);
632 Heps = _mm_mul_ps(vfeps,H);
633 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
634 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
635 fvdw12 = _mm_mul_ps(c12_00,FF);
636 fvdw = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
638 fscal = _mm_add_ps(felec,fvdw);
640 /* Calculate temporary vectorial force */
641 tx = _mm_mul_ps(fscal,dx00);
642 ty = _mm_mul_ps(fscal,dy00);
643 tz = _mm_mul_ps(fscal,dz00);
645 /* Update vectorial force */
646 fix0 = _mm_add_ps(fix0,tx);
647 fiy0 = _mm_add_ps(fiy0,ty);
648 fiz0 = _mm_add_ps(fiz0,tz);
650 fjptrA = f+j_coord_offsetA;
651 fjptrB = f+j_coord_offsetB;
652 fjptrC = f+j_coord_offsetC;
653 fjptrD = f+j_coord_offsetD;
654 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,tx,ty,tz);
656 /* Inner loop uses 54 flops */
662 /* Get j neighbor index, and coordinate index */
663 jnrlistA = jjnr[jidx];
664 jnrlistB = jjnr[jidx+1];
665 jnrlistC = jjnr[jidx+2];
666 jnrlistD = jjnr[jidx+3];
667 /* Sign of each element will be negative for non-real atoms.
668 * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
669 * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
671 dummy_mask = gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
672 jnrA = (jnrlistA>=0) ? jnrlistA : 0;
673 jnrB = (jnrlistB>=0) ? jnrlistB : 0;
674 jnrC = (jnrlistC>=0) ? jnrlistC : 0;
675 jnrD = (jnrlistD>=0) ? jnrlistD : 0;
676 j_coord_offsetA = DIM*jnrA;
677 j_coord_offsetB = DIM*jnrB;
678 j_coord_offsetC = DIM*jnrC;
679 j_coord_offsetD = DIM*jnrD;
681 /* load j atom coordinates */
682 gmx_mm_load_1rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
683 x+j_coord_offsetC,x+j_coord_offsetD,
686 /* Calculate displacement vector */
687 dx00 = _mm_sub_ps(ix0,jx0);
688 dy00 = _mm_sub_ps(iy0,jy0);
689 dz00 = _mm_sub_ps(iz0,jz0);
691 /* Calculate squared distance and things based on it */
692 rsq00 = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
694 rinv00 = sse41_invsqrt_f(rsq00);
696 rinvsq00 = _mm_mul_ps(rinv00,rinv00);
698 /* Load parameters for j particles */
699 jq0 = gmx_mm_load_4real_swizzle_ps(charge+jnrA+0,charge+jnrB+0,
700 charge+jnrC+0,charge+jnrD+0);
701 vdwjidx0A = 2*vdwtype[jnrA+0];
702 vdwjidx0B = 2*vdwtype[jnrB+0];
703 vdwjidx0C = 2*vdwtype[jnrC+0];
704 vdwjidx0D = 2*vdwtype[jnrD+0];
706 /**************************
707 * CALCULATE INTERACTIONS *
708 **************************/
710 r00 = _mm_mul_ps(rsq00,rinv00);
711 r00 = _mm_andnot_ps(dummy_mask,r00);
713 /* Compute parameters for interactions between i and j atoms */
714 qq00 = _mm_mul_ps(iq0,jq0);
715 gmx_mm_load_4pair_swizzle_ps(vdwparam+vdwioffset0+vdwjidx0A,
716 vdwparam+vdwioffset0+vdwjidx0B,
717 vdwparam+vdwioffset0+vdwjidx0C,
718 vdwparam+vdwioffset0+vdwjidx0D,
721 /* Calculate table index by multiplying r with table scale and truncate to integer */
722 rt = _mm_mul_ps(r00,vftabscale);
723 vfitab = _mm_cvttps_epi32(rt);
724 vfeps = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
725 vfitab = _mm_slli_epi32(vfitab,3);
727 /* COULOMB ELECTROSTATICS */
728 velec = _mm_mul_ps(qq00,rinv00);
729 felec = _mm_mul_ps(velec,rinvsq00);
731 /* CUBIC SPLINE TABLE DISPERSION */
732 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
733 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
734 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
735 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
736 _MM_TRANSPOSE4_PS(Y,F,G,H);
737 Heps = _mm_mul_ps(vfeps,H);
738 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
739 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
740 fvdw6 = _mm_mul_ps(c6_00,FF);
742 /* CUBIC SPLINE TABLE REPULSION */
743 vfitab = _mm_add_epi32(vfitab,ifour);
744 Y = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,0) );
745 F = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,1) );
746 G = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,2) );
747 H = _mm_load_ps( vftab + gmx_mm_extract_epi32(vfitab,3) );
748 _MM_TRANSPOSE4_PS(Y,F,G,H);
749 Heps = _mm_mul_ps(vfeps,H);
750 Fp = _mm_add_ps(F,_mm_mul_ps(vfeps,_mm_add_ps(G,Heps)));
751 FF = _mm_add_ps(Fp,_mm_mul_ps(vfeps,_mm_add_ps(G,_mm_add_ps(Heps,Heps))));
752 fvdw12 = _mm_mul_ps(c12_00,FF);
753 fvdw = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
755 fscal = _mm_add_ps(felec,fvdw);
757 fscal = _mm_andnot_ps(dummy_mask,fscal);
759 /* Calculate temporary vectorial force */
760 tx = _mm_mul_ps(fscal,dx00);
761 ty = _mm_mul_ps(fscal,dy00);
762 tz = _mm_mul_ps(fscal,dz00);
764 /* Update vectorial force */
765 fix0 = _mm_add_ps(fix0,tx);
766 fiy0 = _mm_add_ps(fiy0,ty);
767 fiz0 = _mm_add_ps(fiz0,tz);
769 fjptrA = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
770 fjptrB = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
771 fjptrC = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
772 fjptrD = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
773 gmx_mm_decrement_1rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,tx,ty,tz);
775 /* Inner loop uses 55 flops */
778 /* End of innermost loop */
780 gmx_mm_update_iforce_1atom_swizzle_ps(fix0,fiy0,fiz0,
781 f+i_coord_offset,fshift+i_shift_offset);
783 /* Increment number of inner iterations */
784 inneriter += j_index_end - j_index_start;
786 /* Outer loop uses 7 flops */
789 /* Increment number of outer iterations */
792 /* Update outer/inner flops */
794 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_F,outeriter*7 + inneriter*55);