src/gromacs/gmxlib/nonbonded/nb_kernel_avx_128_fma_single/nb_kernel_ElecCSTab_VdwNone_GeomW4W4_avx_128_fma_single.c

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