27f5f11fcc238a17ddc338c53d2fb05c90d832a0
[alexxy/gromacs.git] / src / gmxlib / nonbonded / nb_kernel_avx_128_fma_double / nb_kernel_ElecEwSh_VdwLJSh_GeomP1P1_avx_128_fma_double.c
1 /*
2  * Note: this file was generated by the Gromacs avx_128_fma_double kernel generator.
3  *
4  *                This source code is part of
5  *
6  *                 G   R   O   M   A   C   S
7  *
8  * Copyright (c) 2001-2012, The GROMACS Development Team
9  *
10  * Gromacs is a library for molecular simulation and trajectory analysis,
11  * written by Erik Lindahl, David van der Spoel, Berk Hess, and others - for
12  * a full list of developers and information, check out http://www.gromacs.org
13  *
14  * This program is free software; you can redistribute it and/or modify it under
15  * the terms of the GNU Lesser General Public License as published by the Free
16  * Software Foundation; either version 2 of the License, or (at your option) any
17  * later version.
18  *
19  * To help fund GROMACS development, we humbly ask that you cite
20  * the papers people have written on it - you can find them on the website.
21  */
22 #ifdef HAVE_CONFIG_H
23 #include <config.h>
24 #endif
25
26 #include <math.h>
27
28 #include "../nb_kernel.h"
29 #include "types/simple.h"
30 #include "vec.h"
31 #include "nrnb.h"
32
33 #include "gmx_math_x86_avx_128_fma_double.h"
34 #include "kernelutil_x86_avx_128_fma_double.h"
35
36 /*
37  * Gromacs nonbonded kernel:   nb_kernel_ElecEwSh_VdwLJSh_GeomP1P1_VF_avx_128_fma_double
38  * Electrostatics interaction: Ewald
39  * VdW interaction:            LennardJones
40  * Geometry:                   Particle-Particle
41  * Calculate force/pot:        PotentialAndForce
42  */
43 void
44 nb_kernel_ElecEwSh_VdwLJSh_GeomP1P1_VF_avx_128_fma_double
45                     (t_nblist * gmx_restrict                nlist,
46                      rvec * gmx_restrict                    xx,
47                      rvec * gmx_restrict                    ff,
48                      t_forcerec * gmx_restrict              fr,
49                      t_mdatoms * gmx_restrict               mdatoms,
50                      nb_kernel_data_t * gmx_restrict        kernel_data,
51                      t_nrnb * gmx_restrict                  nrnb)
52 {
53     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
54      * just 0 for non-waters.
55      * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
56      * jnr indices corresponding to data put in the four positions in the SIMD register.
57      */
58     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
59     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
60     int              jnrA,jnrB;
61     int              j_coord_offsetA,j_coord_offsetB;
62     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
63     real             rcutoff_scalar;
64     real             *shiftvec,*fshift,*x,*f;
65     __m128d          tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
66     int              vdwioffset0;
67     __m128d          ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
68     int              vdwjidx0A,vdwjidx0B;
69     __m128d          jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
70     __m128d          dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
71     __m128d          velec,felec,velecsum,facel,crf,krf,krf2;
72     real             *charge;
73     int              nvdwtype;
74     __m128d          rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
75     int              *vdwtype;
76     real             *vdwparam;
77     __m128d          one_sixth   = _mm_set1_pd(1.0/6.0);
78     __m128d          one_twelfth = _mm_set1_pd(1.0/12.0);
79     __m128i          ewitab;
80     __m128d          ewtabscale,eweps,twoeweps,sh_ewald,ewrt,ewtabhalfspace,ewtabF,ewtabFn,ewtabD,ewtabV;
81     real             *ewtab;
82     __m128d          dummy_mask,cutoff_mask;
83     __m128d          signbit   = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
84     __m128d          one     = _mm_set1_pd(1.0);
85     __m128d          two     = _mm_set1_pd(2.0);
86     x                = xx[0];
87     f                = ff[0];
88
89     nri              = nlist->nri;
90     iinr             = nlist->iinr;
91     jindex           = nlist->jindex;
92     jjnr             = nlist->jjnr;
93     shiftidx         = nlist->shift;
94     gid              = nlist->gid;
95     shiftvec         = fr->shift_vec[0];
96     fshift           = fr->fshift[0];
97     facel            = _mm_set1_pd(fr->epsfac);
98     charge           = mdatoms->chargeA;
99     nvdwtype         = fr->ntype;
100     vdwparam         = fr->nbfp;
101     vdwtype          = mdatoms->typeA;
102
103     sh_ewald         = _mm_set1_pd(fr->ic->sh_ewald);
104     ewtab            = fr->ic->tabq_coul_FDV0;
105     ewtabscale       = _mm_set1_pd(fr->ic->tabq_scale);
106     ewtabhalfspace   = _mm_set1_pd(0.5/fr->ic->tabq_scale);
107
108     /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
109     rcutoff_scalar   = fr->rcoulomb;
110     rcutoff          = _mm_set1_pd(rcutoff_scalar);
111     rcutoff2         = _mm_mul_pd(rcutoff,rcutoff);
112
113     sh_vdw_invrcut6  = _mm_set1_pd(fr->ic->sh_invrc6);
114     rvdw             = _mm_set1_pd(fr->rvdw);
115
116     /* Avoid stupid compiler warnings */
117     jnrA = jnrB = 0;
118     j_coord_offsetA = 0;
119     j_coord_offsetB = 0;
120
121     outeriter        = 0;
122     inneriter        = 0;
123
124     /* Start outer loop over neighborlists */
125     for(iidx=0; iidx<nri; iidx++)
126     {
127         /* Load shift vector for this list */
128         i_shift_offset   = DIM*shiftidx[iidx];
129
130         /* Load limits for loop over neighbors */
131         j_index_start    = jindex[iidx];
132         j_index_end      = jindex[iidx+1];
133
134         /* Get outer coordinate index */
135         inr              = iinr[iidx];
136         i_coord_offset   = DIM*inr;
137
138         /* Load i particle coords and add shift vector */
139         gmx_mm_load_shift_and_1rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,&ix0,&iy0,&iz0);
140
141         fix0             = _mm_setzero_pd();
142         fiy0             = _mm_setzero_pd();
143         fiz0             = _mm_setzero_pd();
144
145         /* Load parameters for i particles */
146         iq0              = _mm_mul_pd(facel,_mm_load1_pd(charge+inr+0));
147         vdwioffset0      = 2*nvdwtype*vdwtype[inr+0];
148
149         /* Reset potential sums */
150         velecsum         = _mm_setzero_pd();
151         vvdwsum          = _mm_setzero_pd();
152
153         /* Start inner kernel loop */
154         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
155         {
156
157             /* Get j neighbor index, and coordinate index */
158             jnrA             = jjnr[jidx];
159             jnrB             = jjnr[jidx+1];
160             j_coord_offsetA  = DIM*jnrA;
161             j_coord_offsetB  = DIM*jnrB;
162
163             /* load j atom coordinates */
164             gmx_mm_load_1rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
165                                               &jx0,&jy0,&jz0);
166
167             /* Calculate displacement vector */
168             dx00             = _mm_sub_pd(ix0,jx0);
169             dy00             = _mm_sub_pd(iy0,jy0);
170             dz00             = _mm_sub_pd(iz0,jz0);
171
172             /* Calculate squared distance and things based on it */
173             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
174
175             rinv00           = gmx_mm_invsqrt_pd(rsq00);
176
177             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
178
179             /* Load parameters for j particles */
180             jq0              = gmx_mm_load_2real_swizzle_pd(charge+jnrA+0,charge+jnrB+0);
181             vdwjidx0A        = 2*vdwtype[jnrA+0];
182             vdwjidx0B        = 2*vdwtype[jnrB+0];
183
184             /**************************
185              * CALCULATE INTERACTIONS *
186              **************************/
187
188             if (gmx_mm_any_lt(rsq00,rcutoff2))
189             {
190
191             r00              = _mm_mul_pd(rsq00,rinv00);
192
193             /* Compute parameters for interactions between i and j atoms */
194             qq00             = _mm_mul_pd(iq0,jq0);
195             gmx_mm_load_2pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,
196                                          vdwparam+vdwioffset0+vdwjidx0B,&c6_00,&c12_00);
197
198             /* EWALD ELECTROSTATICS */
199
200             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
201             ewrt             = _mm_mul_pd(r00,ewtabscale);
202             ewitab           = _mm_cvttpd_epi32(ewrt);
203 #ifdef __XOP__
204             eweps            = _mm_frcz_pd(ewrt);
205 #else
206             eweps            = _mm_sub_pd(ewrt,_mm_round_pd(ewrt, _MM_FROUND_FLOOR));
207 #endif
208             twoeweps         = _mm_add_pd(eweps,eweps);
209             ewitab           = _mm_slli_epi32(ewitab,2);
210             ewtabF           = _mm_load_pd( ewtab + _mm_extract_epi32(ewitab,0) );
211             ewtabD           = _mm_load_pd( ewtab + _mm_extract_epi32(ewitab,1) );
212             GMX_MM_TRANSPOSE2_PD(ewtabF,ewtabD);
213             ewtabV           = _mm_load_sd( ewtab + _mm_extract_epi32(ewitab,0) +2);
214             ewtabFn          = _mm_load_sd( ewtab + _mm_extract_epi32(ewitab,1) +2);
215             GMX_MM_TRANSPOSE2_PD(ewtabV,ewtabFn);
216             felec            = _mm_macc_pd(eweps,ewtabD,ewtabF);
217             velec            = _mm_nmacc_pd(_mm_mul_pd(ewtabhalfspace,eweps) ,_mm_add_pd(ewtabF,felec), ewtabV);
218             velec            = _mm_mul_pd(qq00,_mm_sub_pd(_mm_sub_pd(rinv00,sh_ewald),velec));
219             felec            = _mm_mul_pd(_mm_mul_pd(qq00,rinv00),_mm_sub_pd(rinvsq00,felec));
220
221             /* LENNARD-JONES DISPERSION/REPULSION */
222
223             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
224             vvdw6            = _mm_mul_pd(c6_00,rinvsix);
225             vvdw12           = _mm_mul_pd(c12_00,_mm_mul_pd(rinvsix,rinvsix));
226             vvdw             = _mm_msub_pd(_mm_nmacc_pd(c12_00,_mm_mul_pd(sh_vdw_invrcut6,sh_vdw_invrcut6),vvdw12),one_twelfth,
227                                            _mm_mul_pd(_mm_nmacc_pd( c6_00,sh_vdw_invrcut6,vvdw6),one_sixth));
228             fvdw             = _mm_mul_pd(_mm_sub_pd(vvdw12,vvdw6),rinvsq00);
229
230             cutoff_mask      = _mm_cmplt_pd(rsq00,rcutoff2);
231
232             /* Update potential sum for this i atom from the interaction with this j atom. */
233             velec            = _mm_and_pd(velec,cutoff_mask);
234             velecsum         = _mm_add_pd(velecsum,velec);
235             vvdw             = _mm_and_pd(vvdw,cutoff_mask);
236             vvdwsum          = _mm_add_pd(vvdwsum,vvdw);
237
238             fscal            = _mm_add_pd(felec,fvdw);
239
240             fscal            = _mm_and_pd(fscal,cutoff_mask);
241
242             /* Update vectorial force */
243             fix0             = _mm_macc_pd(dx00,fscal,fix0);
244             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
245             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
246             
247             gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,
248                                                    _mm_mul_pd(dx00,fscal),
249                                                    _mm_mul_pd(dy00,fscal),
250                                                    _mm_mul_pd(dz00,fscal));
251
252             }
253
254             /* Inner loop uses 67 flops */
255         }
256
257         if(jidx<j_index_end)
258         {
259
260             jnrA             = jjnr[jidx];
261             j_coord_offsetA  = DIM*jnrA;
262
263             /* load j atom coordinates */
264             gmx_mm_load_1rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
265                                               &jx0,&jy0,&jz0);
266
267             /* Calculate displacement vector */
268             dx00             = _mm_sub_pd(ix0,jx0);
269             dy00             = _mm_sub_pd(iy0,jy0);
270             dz00             = _mm_sub_pd(iz0,jz0);
271
272             /* Calculate squared distance and things based on it */
273             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
274
275             rinv00           = gmx_mm_invsqrt_pd(rsq00);
276
277             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
278
279             /* Load parameters for j particles */
280             jq0              = _mm_load_sd(charge+jnrA+0);
281             vdwjidx0A        = 2*vdwtype[jnrA+0];
282
283             /**************************
284              * CALCULATE INTERACTIONS *
285              **************************/
286
287             if (gmx_mm_any_lt(rsq00,rcutoff2))
288             {
289
290             r00              = _mm_mul_pd(rsq00,rinv00);
291
292             /* Compute parameters for interactions between i and j atoms */
293             qq00             = _mm_mul_pd(iq0,jq0);
294             gmx_mm_load_1pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,&c6_00,&c12_00);
295
296             /* EWALD ELECTROSTATICS */
297
298             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
299             ewrt             = _mm_mul_pd(r00,ewtabscale);
300             ewitab           = _mm_cvttpd_epi32(ewrt);
301 #ifdef __XOP__
302             eweps            = _mm_frcz_pd(ewrt);
303 #else
304             eweps            = _mm_sub_pd(ewrt,_mm_round_pd(ewrt, _MM_FROUND_FLOOR));
305 #endif
306             twoeweps         = _mm_add_pd(eweps,eweps);
307             ewitab           = _mm_slli_epi32(ewitab,2);
308             ewtabF           = _mm_load_pd( ewtab + _mm_extract_epi32(ewitab,0) );
309             ewtabD           = _mm_setzero_pd();
310             GMX_MM_TRANSPOSE2_PD(ewtabF,ewtabD);
311             ewtabV           = _mm_load_sd( ewtab + _mm_extract_epi32(ewitab,0) +2);
312             ewtabFn          = _mm_setzero_pd();
313             GMX_MM_TRANSPOSE2_PD(ewtabV,ewtabFn);
314             felec            = _mm_macc_pd(eweps,ewtabD,ewtabF);
315             velec            = _mm_nmacc_pd(_mm_mul_pd(ewtabhalfspace,eweps) ,_mm_add_pd(ewtabF,felec), ewtabV);
316             velec            = _mm_mul_pd(qq00,_mm_sub_pd(_mm_sub_pd(rinv00,sh_ewald),velec));
317             felec            = _mm_mul_pd(_mm_mul_pd(qq00,rinv00),_mm_sub_pd(rinvsq00,felec));
318
319             /* LENNARD-JONES DISPERSION/REPULSION */
320
321             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
322             vvdw6            = _mm_mul_pd(c6_00,rinvsix);
323             vvdw12           = _mm_mul_pd(c12_00,_mm_mul_pd(rinvsix,rinvsix));
324             vvdw             = _mm_msub_pd(_mm_nmacc_pd(c12_00,_mm_mul_pd(sh_vdw_invrcut6,sh_vdw_invrcut6),vvdw12),one_twelfth,
325                                            _mm_mul_pd(_mm_nmacc_pd( c6_00,sh_vdw_invrcut6,vvdw6),one_sixth));
326             fvdw             = _mm_mul_pd(_mm_sub_pd(vvdw12,vvdw6),rinvsq00);
327
328             cutoff_mask      = _mm_cmplt_pd(rsq00,rcutoff2);
329
330             /* Update potential sum for this i atom from the interaction with this j atom. */
331             velec            = _mm_and_pd(velec,cutoff_mask);
332             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
333             velecsum         = _mm_add_pd(velecsum,velec);
334             vvdw             = _mm_and_pd(vvdw,cutoff_mask);
335             vvdw             = _mm_unpacklo_pd(vvdw,_mm_setzero_pd());
336             vvdwsum          = _mm_add_pd(vvdwsum,vvdw);
337
338             fscal            = _mm_add_pd(felec,fvdw);
339
340             fscal            = _mm_and_pd(fscal,cutoff_mask);
341
342             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
343
344             /* Update vectorial force */
345             fix0             = _mm_macc_pd(dx00,fscal,fix0);
346             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
347             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
348             
349             gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f+j_coord_offsetA,
350                                                    _mm_mul_pd(dx00,fscal),
351                                                    _mm_mul_pd(dy00,fscal),
352                                                    _mm_mul_pd(dz00,fscal));
353
354             }
355
356             /* Inner loop uses 67 flops */
357         }
358
359         /* End of innermost loop */
360
361         gmx_mm_update_iforce_1atom_swizzle_pd(fix0,fiy0,fiz0,
362                                               f+i_coord_offset,fshift+i_shift_offset);
363
364         ggid                        = gid[iidx];
365         /* Update potential energies */
366         gmx_mm_update_1pot_pd(velecsum,kernel_data->energygrp_elec+ggid);
367         gmx_mm_update_1pot_pd(vvdwsum,kernel_data->energygrp_vdw+ggid);
368
369         /* Increment number of inner iterations */
370         inneriter                  += j_index_end - j_index_start;
371
372         /* Outer loop uses 9 flops */
373     }
374
375     /* Increment number of outer iterations */
376     outeriter        += nri;
377
378     /* Update outer/inner flops */
379
380     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_VF,outeriter*9 + inneriter*67);
381 }
382 /*
383  * Gromacs nonbonded kernel:   nb_kernel_ElecEwSh_VdwLJSh_GeomP1P1_F_avx_128_fma_double
384  * Electrostatics interaction: Ewald
385  * VdW interaction:            LennardJones
386  * Geometry:                   Particle-Particle
387  * Calculate force/pot:        Force
388  */
389 void
390 nb_kernel_ElecEwSh_VdwLJSh_GeomP1P1_F_avx_128_fma_double
391                     (t_nblist * gmx_restrict                nlist,
392                      rvec * gmx_restrict                    xx,
393                      rvec * gmx_restrict                    ff,
394                      t_forcerec * gmx_restrict              fr,
395                      t_mdatoms * gmx_restrict               mdatoms,
396                      nb_kernel_data_t * gmx_restrict        kernel_data,
397                      t_nrnb * gmx_restrict                  nrnb)
398 {
399     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
400      * just 0 for non-waters.
401      * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
402      * jnr indices corresponding to data put in the four positions in the SIMD register.
403      */
404     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
405     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
406     int              jnrA,jnrB;
407     int              j_coord_offsetA,j_coord_offsetB;
408     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
409     real             rcutoff_scalar;
410     real             *shiftvec,*fshift,*x,*f;
411     __m128d          tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
412     int              vdwioffset0;
413     __m128d          ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
414     int              vdwjidx0A,vdwjidx0B;
415     __m128d          jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
416     __m128d          dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
417     __m128d          velec,felec,velecsum,facel,crf,krf,krf2;
418     real             *charge;
419     int              nvdwtype;
420     __m128d          rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
421     int              *vdwtype;
422     real             *vdwparam;
423     __m128d          one_sixth   = _mm_set1_pd(1.0/6.0);
424     __m128d          one_twelfth = _mm_set1_pd(1.0/12.0);
425     __m128i          ewitab;
426     __m128d          ewtabscale,eweps,twoeweps,sh_ewald,ewrt,ewtabhalfspace,ewtabF,ewtabFn,ewtabD,ewtabV;
427     real             *ewtab;
428     __m128d          dummy_mask,cutoff_mask;
429     __m128d          signbit   = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
430     __m128d          one     = _mm_set1_pd(1.0);
431     __m128d          two     = _mm_set1_pd(2.0);
432     x                = xx[0];
433     f                = ff[0];
434
435     nri              = nlist->nri;
436     iinr             = nlist->iinr;
437     jindex           = nlist->jindex;
438     jjnr             = nlist->jjnr;
439     shiftidx         = nlist->shift;
440     gid              = nlist->gid;
441     shiftvec         = fr->shift_vec[0];
442     fshift           = fr->fshift[0];
443     facel            = _mm_set1_pd(fr->epsfac);
444     charge           = mdatoms->chargeA;
445     nvdwtype         = fr->ntype;
446     vdwparam         = fr->nbfp;
447     vdwtype          = mdatoms->typeA;
448
449     sh_ewald         = _mm_set1_pd(fr->ic->sh_ewald);
450     ewtab            = fr->ic->tabq_coul_F;
451     ewtabscale       = _mm_set1_pd(fr->ic->tabq_scale);
452     ewtabhalfspace   = _mm_set1_pd(0.5/fr->ic->tabq_scale);
453
454     /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
455     rcutoff_scalar   = fr->rcoulomb;
456     rcutoff          = _mm_set1_pd(rcutoff_scalar);
457     rcutoff2         = _mm_mul_pd(rcutoff,rcutoff);
458
459     sh_vdw_invrcut6  = _mm_set1_pd(fr->ic->sh_invrc6);
460     rvdw             = _mm_set1_pd(fr->rvdw);
461
462     /* Avoid stupid compiler warnings */
463     jnrA = jnrB = 0;
464     j_coord_offsetA = 0;
465     j_coord_offsetB = 0;
466
467     outeriter        = 0;
468     inneriter        = 0;
469
470     /* Start outer loop over neighborlists */
471     for(iidx=0; iidx<nri; iidx++)
472     {
473         /* Load shift vector for this list */
474         i_shift_offset   = DIM*shiftidx[iidx];
475
476         /* Load limits for loop over neighbors */
477         j_index_start    = jindex[iidx];
478         j_index_end      = jindex[iidx+1];
479
480         /* Get outer coordinate index */
481         inr              = iinr[iidx];
482         i_coord_offset   = DIM*inr;
483
484         /* Load i particle coords and add shift vector */
485         gmx_mm_load_shift_and_1rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,&ix0,&iy0,&iz0);
486
487         fix0             = _mm_setzero_pd();
488         fiy0             = _mm_setzero_pd();
489         fiz0             = _mm_setzero_pd();
490
491         /* Load parameters for i particles */
492         iq0              = _mm_mul_pd(facel,_mm_load1_pd(charge+inr+0));
493         vdwioffset0      = 2*nvdwtype*vdwtype[inr+0];
494
495         /* Start inner kernel loop */
496         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
497         {
498
499             /* Get j neighbor index, and coordinate index */
500             jnrA             = jjnr[jidx];
501             jnrB             = jjnr[jidx+1];
502             j_coord_offsetA  = DIM*jnrA;
503             j_coord_offsetB  = DIM*jnrB;
504
505             /* load j atom coordinates */
506             gmx_mm_load_1rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
507                                               &jx0,&jy0,&jz0);
508
509             /* Calculate displacement vector */
510             dx00             = _mm_sub_pd(ix0,jx0);
511             dy00             = _mm_sub_pd(iy0,jy0);
512             dz00             = _mm_sub_pd(iz0,jz0);
513
514             /* Calculate squared distance and things based on it */
515             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
516
517             rinv00           = gmx_mm_invsqrt_pd(rsq00);
518
519             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
520
521             /* Load parameters for j particles */
522             jq0              = gmx_mm_load_2real_swizzle_pd(charge+jnrA+0,charge+jnrB+0);
523             vdwjidx0A        = 2*vdwtype[jnrA+0];
524             vdwjidx0B        = 2*vdwtype[jnrB+0];
525
526             /**************************
527              * CALCULATE INTERACTIONS *
528              **************************/
529
530             if (gmx_mm_any_lt(rsq00,rcutoff2))
531             {
532
533             r00              = _mm_mul_pd(rsq00,rinv00);
534
535             /* Compute parameters for interactions between i and j atoms */
536             qq00             = _mm_mul_pd(iq0,jq0);
537             gmx_mm_load_2pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,
538                                          vdwparam+vdwioffset0+vdwjidx0B,&c6_00,&c12_00);
539
540             /* EWALD ELECTROSTATICS */
541
542             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
543             ewrt             = _mm_mul_pd(r00,ewtabscale);
544             ewitab           = _mm_cvttpd_epi32(ewrt);
545 #ifdef __XOP__
546             eweps            = _mm_frcz_pd(ewrt);
547 #else
548             eweps            = _mm_sub_pd(ewrt,_mm_round_pd(ewrt, _MM_FROUND_FLOOR));
549 #endif
550             twoeweps         = _mm_add_pd(eweps,eweps);
551             gmx_mm_load_2pair_swizzle_pd(ewtab+_mm_extract_epi32(ewitab,0),ewtab+_mm_extract_epi32(ewitab,1),
552                                          &ewtabF,&ewtabFn);
553             felec            = _mm_macc_pd(eweps,ewtabFn,_mm_mul_pd( _mm_sub_pd(one,eweps),ewtabF));
554             felec            = _mm_mul_pd(_mm_mul_pd(qq00,rinv00),_mm_sub_pd(rinvsq00,felec));
555
556             /* LENNARD-JONES DISPERSION/REPULSION */
557
558             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
559             fvdw             = _mm_mul_pd(_mm_msub_pd(c12_00,rinvsix,c6_00),_mm_mul_pd(rinvsix,rinvsq00));
560
561             cutoff_mask      = _mm_cmplt_pd(rsq00,rcutoff2);
562
563             fscal            = _mm_add_pd(felec,fvdw);
564
565             fscal            = _mm_and_pd(fscal,cutoff_mask);
566
567             /* Update vectorial force */
568             fix0             = _mm_macc_pd(dx00,fscal,fix0);
569             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
570             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
571             
572             gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,
573                                                    _mm_mul_pd(dx00,fscal),
574                                                    _mm_mul_pd(dy00,fscal),
575                                                    _mm_mul_pd(dz00,fscal));
576
577             }
578
579             /* Inner loop uses 49 flops */
580         }
581
582         if(jidx<j_index_end)
583         {
584
585             jnrA             = jjnr[jidx];
586             j_coord_offsetA  = DIM*jnrA;
587
588             /* load j atom coordinates */
589             gmx_mm_load_1rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
590                                               &jx0,&jy0,&jz0);
591
592             /* Calculate displacement vector */
593             dx00             = _mm_sub_pd(ix0,jx0);
594             dy00             = _mm_sub_pd(iy0,jy0);
595             dz00             = _mm_sub_pd(iz0,jz0);
596
597             /* Calculate squared distance and things based on it */
598             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
599
600             rinv00           = gmx_mm_invsqrt_pd(rsq00);
601
602             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
603
604             /* Load parameters for j particles */
605             jq0              = _mm_load_sd(charge+jnrA+0);
606             vdwjidx0A        = 2*vdwtype[jnrA+0];
607
608             /**************************
609              * CALCULATE INTERACTIONS *
610              **************************/
611
612             if (gmx_mm_any_lt(rsq00,rcutoff2))
613             {
614
615             r00              = _mm_mul_pd(rsq00,rinv00);
616
617             /* Compute parameters for interactions between i and j atoms */
618             qq00             = _mm_mul_pd(iq0,jq0);
619             gmx_mm_load_1pair_swizzle_pd(vdwparam+vdwioffset0+vdwjidx0A,&c6_00,&c12_00);
620
621             /* EWALD ELECTROSTATICS */
622
623             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
624             ewrt             = _mm_mul_pd(r00,ewtabscale);
625             ewitab           = _mm_cvttpd_epi32(ewrt);
626 #ifdef __XOP__
627             eweps            = _mm_frcz_pd(ewrt);
628 #else
629             eweps            = _mm_sub_pd(ewrt,_mm_round_pd(ewrt, _MM_FROUND_FLOOR));
630 #endif
631             twoeweps         = _mm_add_pd(eweps,eweps);
632             gmx_mm_load_1pair_swizzle_pd(ewtab+_mm_extract_epi32(ewitab,0),&ewtabF,&ewtabFn);
633             felec            = _mm_macc_pd(eweps,ewtabFn,_mm_mul_pd( _mm_sub_pd(one,eweps),ewtabF));
634             felec            = _mm_mul_pd(_mm_mul_pd(qq00,rinv00),_mm_sub_pd(rinvsq00,felec));
635
636             /* LENNARD-JONES DISPERSION/REPULSION */
637
638             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
639             fvdw             = _mm_mul_pd(_mm_msub_pd(c12_00,rinvsix,c6_00),_mm_mul_pd(rinvsix,rinvsq00));
640
641             cutoff_mask      = _mm_cmplt_pd(rsq00,rcutoff2);
642
643             fscal            = _mm_add_pd(felec,fvdw);
644
645             fscal            = _mm_and_pd(fscal,cutoff_mask);
646
647             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
648
649             /* Update vectorial force */
650             fix0             = _mm_macc_pd(dx00,fscal,fix0);
651             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
652             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
653             
654             gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f+j_coord_offsetA,
655                                                    _mm_mul_pd(dx00,fscal),
656                                                    _mm_mul_pd(dy00,fscal),
657                                                    _mm_mul_pd(dz00,fscal));
658
659             }
660
661             /* Inner loop uses 49 flops */
662         }
663
664         /* End of innermost loop */
665
666         gmx_mm_update_iforce_1atom_swizzle_pd(fix0,fiy0,fiz0,
667                                               f+i_coord_offset,fshift+i_shift_offset);
668
669         /* Increment number of inner iterations */
670         inneriter                  += j_index_end - j_index_start;
671
672         /* Outer loop uses 7 flops */
673     }
674
675     /* Increment number of outer iterations */
676     outeriter        += nri;
677
678     /* Update outer/inner flops */
679
680     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_F,outeriter*7 + inneriter*49);
681 }