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