src/gromacs/gmxlib/nonbonded/nb_kernel_avx_128_fma_double/nb_kernel_ElecCSTab_VdwLJ_GeomW4W4_avx_128_fma_double.c

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