src/gromacs/gmxlib/nonbonded/nb_kernel_avx_128_fma_single/nb_kernel_ElecRF_VdwCSTab_GeomW3W3_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_ElecRF_VdwCSTab_GeomW3W3_VF_avx_128_fma_single
  52  * Electrostatics interaction: ReactionField
  53  * VdW interaction:            CubicSplineTable
  54  * Geometry:                   Water3-Water3
  55  * Calculate force/pot:        PotentialAndForce
  56  */
  57 void
  58 nb_kernel_ElecRF_VdwCSTab_GeomW3W3_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              vdwioffset0;
  84     __m128           ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
  85     int              vdwioffset1;
  86     __m128           ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
  87     int              vdwioffset2;
  88     __m128           ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
  89     int              vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
  90     __m128           jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
  91     int              vdwjidx1A,vdwjidx1B,vdwjidx1C,vdwjidx1D;
  92     __m128           jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
  93     int              vdwjidx2A,vdwjidx2B,vdwjidx2C,vdwjidx2D;
  94     __m128           jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
  95     __m128           dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
  96     __m128           dx01,dy01,dz01,rsq01,rinv01,rinvsq01,r01,qq01,c6_01,c12_01;
  97     __m128           dx02,dy02,dz02,rsq02,rinv02,rinvsq02,r02,qq02,c6_02,c12_02;
  98     __m128           dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
  99     __m128           dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
 100     __m128           dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
 101     __m128           dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
 102     __m128           dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
 103     __m128           dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
 104     __m128           velec,felec,velecsum,facel,crf,krf,krf2;
 105     real             *charge;
 106     int              nvdwtype;
 107     __m128           rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
 108     int              *vdwtype;
 109     real             *vdwparam;
 110     __m128           one_sixth   = _mm_set1_ps(1.0/6.0);
 111     __m128           one_twelfth = _mm_set1_ps(1.0/12.0);
 112     __m128i          vfitab;
 113     __m128i          ifour       = _mm_set1_epi32(4);
 114     __m128           rt,vfeps,twovfeps,vftabscale,Y,F,G,H,Fp,VV,FF;
 115     real             *vftab;
 116     __m128           dummy_mask,cutoff_mask;
 117     __m128           signbit = _mm_castsi128_ps( _mm_set1_epi32(0x80000000) );
 118     __m128           one     = _mm_set1_ps(1.0);
 119     __m128           two     = _mm_set1_ps(2.0);
 120     x                = xx[0];
 121     f                = ff[0];
 122
 123     nri              = nlist->nri;
 124     iinr             = nlist->iinr;
 125     jindex           = nlist->jindex;
 126     jjnr             = nlist->jjnr;
 127     shiftidx         = nlist->shift;
 128     gid              = nlist->gid;
 129     shiftvec         = fr->shift_vec[0];
 130     fshift           = fr->fshift[0];
 131     facel            = _mm_set1_ps(fr->epsfac);
 132     charge           = mdatoms->chargeA;
 133     krf              = _mm_set1_ps(fr->ic->k_rf);
 134     krf2             = _mm_set1_ps(fr->ic->k_rf*2.0);
 135     crf              = _mm_set1_ps(fr->ic->c_rf);
 136     nvdwtype         = fr->ntype;
 137     vdwparam         = fr->nbfp;
 138     vdwtype          = mdatoms->typeA;
 139
 140     vftab            = kernel_data->table_vdw->data;
 141     vftabscale       = _mm_set1_ps(kernel_data->table_vdw->scale);
 142
 143     /* Setup water-specific parameters */
 144     inr              = nlist->iinr[0];
 145     iq0              = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+0]));
 146     iq1              = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+1]));
 147     iq2              = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+2]));
 148     vdwioffset0      = 2*nvdwtype*vdwtype[inr+0];
 149
 150     jq0              = _mm_set1_ps(charge[inr+0]);
 151     jq1              = _mm_set1_ps(charge[inr+1]);
 152     jq2              = _mm_set1_ps(charge[inr+2]);
 153     vdwjidx0A        = 2*vdwtype[inr+0];
 154     qq00             = _mm_mul_ps(iq0,jq0);
 155     c6_00            = _mm_set1_ps(vdwparam[vdwioffset0+vdwjidx0A]);
 156     c12_00           = _mm_set1_ps(vdwparam[vdwioffset0+vdwjidx0A+1]);
 157     qq01             = _mm_mul_ps(iq0,jq1);
 158     qq02             = _mm_mul_ps(iq0,jq2);
 159     qq10             = _mm_mul_ps(iq1,jq0);
 160     qq11             = _mm_mul_ps(iq1,jq1);
 161     qq12             = _mm_mul_ps(iq1,jq2);
 162     qq20             = _mm_mul_ps(iq2,jq0);
 163     qq21             = _mm_mul_ps(iq2,jq1);
 164     qq22             = _mm_mul_ps(iq2,jq2);
 165
 166     /* Avoid stupid compiler warnings */
 167     jnrA = jnrB = jnrC = jnrD = 0;
 168     j_coord_offsetA = 0;
 169     j_coord_offsetB = 0;
 170     j_coord_offsetC = 0;
 171     j_coord_offsetD = 0;
 172
 173     outeriter        = 0;
 174     inneriter        = 0;
 175
 176     for(iidx=0;iidx<4*DIM;iidx++)
 177     {
 178         scratch[iidx] = 0.0;
 179     }
 180
 181     /* Start outer loop over neighborlists */
 182     for(iidx=0; iidx<nri; iidx++)
 183     {
 184         /* Load shift vector for this list */
 185         i_shift_offset   = DIM*shiftidx[iidx];
 186
 187         /* Load limits for loop over neighbors */
 188         j_index_start    = jindex[iidx];
 189         j_index_end      = jindex[iidx+1];
 190
 191         /* Get outer coordinate index */
 192         inr              = iinr[iidx];
 193         i_coord_offset   = DIM*inr;
 194
 195         /* Load i particle coords and add shift vector */
 196         gmx_mm_load_shift_and_3rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,
 197                                                  &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2);
 198
 199         fix0             = _mm_setzero_ps();
 200         fiy0             = _mm_setzero_ps();
 201         fiz0             = _mm_setzero_ps();
 202         fix1             = _mm_setzero_ps();
 203         fiy1             = _mm_setzero_ps();
 204         fiz1             = _mm_setzero_ps();
 205         fix2             = _mm_setzero_ps();
 206         fiy2             = _mm_setzero_ps();
 207         fiz2             = _mm_setzero_ps();
 208
 209         /* Reset potential sums */
 210         velecsum         = _mm_setzero_ps();
 211         vvdwsum          = _mm_setzero_ps();
 212
 213         /* Start inner kernel loop */
 214         for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
 215         {
 216
 217             /* Get j neighbor index, and coordinate index */
 218             jnrA             = jjnr[jidx];
 219             jnrB             = jjnr[jidx+1];
 220             jnrC             = jjnr[jidx+2];
 221             jnrD             = jjnr[jidx+3];
 222             j_coord_offsetA  = DIM*jnrA;
 223             j_coord_offsetB  = DIM*jnrB;
 224             j_coord_offsetC  = DIM*jnrC;
 225             j_coord_offsetD  = DIM*jnrD;
 226
 227             /* load j atom coordinates */
 228             gmx_mm_load_3rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
 229                                               x+j_coord_offsetC,x+j_coord_offsetD,
 230                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
 231
 232             /* Calculate displacement vector */
 233             dx00             = _mm_sub_ps(ix0,jx0);
 234             dy00             = _mm_sub_ps(iy0,jy0);
 235             dz00             = _mm_sub_ps(iz0,jz0);
 236             dx01             = _mm_sub_ps(ix0,jx1);
 237             dy01             = _mm_sub_ps(iy0,jy1);
 238             dz01             = _mm_sub_ps(iz0,jz1);
 239             dx02             = _mm_sub_ps(ix0,jx2);
 240             dy02             = _mm_sub_ps(iy0,jy2);
 241             dz02             = _mm_sub_ps(iz0,jz2);
 242             dx10             = _mm_sub_ps(ix1,jx0);
 243             dy10             = _mm_sub_ps(iy1,jy0);
 244             dz10             = _mm_sub_ps(iz1,jz0);
 245             dx11             = _mm_sub_ps(ix1,jx1);
 246             dy11             = _mm_sub_ps(iy1,jy1);
 247             dz11             = _mm_sub_ps(iz1,jz1);
 248             dx12             = _mm_sub_ps(ix1,jx2);
 249             dy12             = _mm_sub_ps(iy1,jy2);
 250             dz12             = _mm_sub_ps(iz1,jz2);
 251             dx20             = _mm_sub_ps(ix2,jx0);
 252             dy20             = _mm_sub_ps(iy2,jy0);
 253             dz20             = _mm_sub_ps(iz2,jz0);
 254             dx21             = _mm_sub_ps(ix2,jx1);
 255             dy21             = _mm_sub_ps(iy2,jy1);
 256             dz21             = _mm_sub_ps(iz2,jz1);
 257             dx22             = _mm_sub_ps(ix2,jx2);
 258             dy22             = _mm_sub_ps(iy2,jy2);
 259             dz22             = _mm_sub_ps(iz2,jz2);
 260
 261             /* Calculate squared distance and things based on it */
 262             rsq00            = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
 263             rsq01            = gmx_mm_calc_rsq_ps(dx01,dy01,dz01);
 264             rsq02            = gmx_mm_calc_rsq_ps(dx02,dy02,dz02);
 265             rsq10            = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
 266             rsq11            = gmx_mm_calc_rsq_ps(dx11,dy11,dz11);
 267             rsq12            = gmx_mm_calc_rsq_ps(dx12,dy12,dz12);
 268             rsq20            = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
 269             rsq21            = gmx_mm_calc_rsq_ps(dx21,dy21,dz21);
 270             rsq22            = gmx_mm_calc_rsq_ps(dx22,dy22,dz22);
 271
 272             rinv00           = gmx_mm_invsqrt_ps(rsq00);
 273             rinv01           = gmx_mm_invsqrt_ps(rsq01);
 274             rinv02           = gmx_mm_invsqrt_ps(rsq02);
 275             rinv10           = gmx_mm_invsqrt_ps(rsq10);
 276             rinv11           = gmx_mm_invsqrt_ps(rsq11);
 277             rinv12           = gmx_mm_invsqrt_ps(rsq12);
 278             rinv20           = gmx_mm_invsqrt_ps(rsq20);
 279             rinv21           = gmx_mm_invsqrt_ps(rsq21);
 280             rinv22           = gmx_mm_invsqrt_ps(rsq22);
 281
 282             rinvsq00         = _mm_mul_ps(rinv00,rinv00);
 283             rinvsq01         = _mm_mul_ps(rinv01,rinv01);
 284             rinvsq02         = _mm_mul_ps(rinv02,rinv02);
 285             rinvsq10         = _mm_mul_ps(rinv10,rinv10);
 286             rinvsq11         = _mm_mul_ps(rinv11,rinv11);
 287             rinvsq12         = _mm_mul_ps(rinv12,rinv12);
 288             rinvsq20         = _mm_mul_ps(rinv20,rinv20);
 289             rinvsq21         = _mm_mul_ps(rinv21,rinv21);
 290             rinvsq22         = _mm_mul_ps(rinv22,rinv22);
 291
 292             fjx0             = _mm_setzero_ps();
 293             fjy0             = _mm_setzero_ps();
 294             fjz0             = _mm_setzero_ps();
 295             fjx1             = _mm_setzero_ps();
 296             fjy1             = _mm_setzero_ps();
 297             fjz1             = _mm_setzero_ps();
 298             fjx2             = _mm_setzero_ps();
 299             fjy2             = _mm_setzero_ps();
 300             fjz2             = _mm_setzero_ps();
 301
 302             /**************************
 303              * CALCULATE INTERACTIONS *
 304              **************************/
 305
 306             r00              = _mm_mul_ps(rsq00,rinv00);
 307
 308             /* Calculate table index by multiplying r with table scale and truncate to integer */
 309             rt               = _mm_mul_ps(r00,vftabscale);
 310             vfitab           = _mm_cvttps_epi32(rt);
 311 #ifdef __XOP__
 312             vfeps            = _mm_frcz_ps(rt);
 313 #else
 314             vfeps            = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
 315 #endif
 316             twovfeps         = _mm_add_ps(vfeps,vfeps);
 317             vfitab           = _mm_slli_epi32(vfitab,3);
 318
 319             /* REACTION-FIELD ELECTROSTATICS */
 320             velec            = _mm_mul_ps(qq00,_mm_sub_ps(_mm_macc_ps(krf,rsq00,rinv00),crf));
 321             felec            = _mm_mul_ps(qq00,_mm_msub_ps(rinv00,rinvsq00,krf2));
 322
 323             /* CUBIC SPLINE TABLE DISPERSION */
 324             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
 325             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
 326             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
 327             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
 328             _MM_TRANSPOSE4_PS(Y,F,G,H);
 329             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
 330             VV               = _mm_macc_ps(vfeps,Fp,Y);
 331             vvdw6            = _mm_mul_ps(c6_00,VV);
 332             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
 333             fvdw6            = _mm_mul_ps(c6_00,FF);
 334
 335             /* CUBIC SPLINE TABLE REPULSION */
 336             vfitab           = _mm_add_epi32(vfitab,ifour);
 337             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
 338             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
 339             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
 340             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
 341             _MM_TRANSPOSE4_PS(Y,F,G,H);
 342             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
 343             VV               = _mm_macc_ps(vfeps,Fp,Y);
 344             vvdw12           = _mm_mul_ps(c12_00,VV);
 345             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
 346             fvdw12           = _mm_mul_ps(c12_00,FF);
 347             vvdw             = _mm_add_ps(vvdw12,vvdw6);
 348             fvdw             = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
 349
 350             /* Update potential sum for this i atom from the interaction with this j atom. */
 351             velecsum         = _mm_add_ps(velecsum,velec);
 352             vvdwsum          = _mm_add_ps(vvdwsum,vvdw);
 353
 354             fscal            = _mm_add_ps(felec,fvdw);
 355
 356              /* Update vectorial force */
 357             fix0             = _mm_macc_ps(dx00,fscal,fix0);
 358             fiy0             = _mm_macc_ps(dy00,fscal,fiy0);
 359             fiz0             = _mm_macc_ps(dz00,fscal,fiz0);
 360
 361             fjx0             = _mm_macc_ps(dx00,fscal,fjx0);
 362             fjy0             = _mm_macc_ps(dy00,fscal,fjy0);
 363             fjz0             = _mm_macc_ps(dz00,fscal,fjz0);
 364
 365             /**************************
 366              * CALCULATE INTERACTIONS *
 367              **************************/
 368
 369             /* REACTION-FIELD ELECTROSTATICS */
 370             velec            = _mm_mul_ps(qq01,_mm_sub_ps(_mm_macc_ps(krf,rsq01,rinv01),crf));
 371             felec            = _mm_mul_ps(qq01,_mm_msub_ps(rinv01,rinvsq01,krf2));
 372
 373             /* Update potential sum for this i atom from the interaction with this j atom. */
 374             velecsum         = _mm_add_ps(velecsum,velec);
 375
 376             fscal            = felec;
 377
 378              /* Update vectorial force */
 379             fix0             = _mm_macc_ps(dx01,fscal,fix0);
 380             fiy0             = _mm_macc_ps(dy01,fscal,fiy0);
 381             fiz0             = _mm_macc_ps(dz01,fscal,fiz0);
 382
 383             fjx1             = _mm_macc_ps(dx01,fscal,fjx1);
 384             fjy1             = _mm_macc_ps(dy01,fscal,fjy1);
 385             fjz1             = _mm_macc_ps(dz01,fscal,fjz1);
 386
 387             /**************************
 388              * CALCULATE INTERACTIONS *
 389              **************************/
 390
 391             /* REACTION-FIELD ELECTROSTATICS */
 392             velec            = _mm_mul_ps(qq02,_mm_sub_ps(_mm_macc_ps(krf,rsq02,rinv02),crf));
 393             felec            = _mm_mul_ps(qq02,_mm_msub_ps(rinv02,rinvsq02,krf2));
 394
 395             /* Update potential sum for this i atom from the interaction with this j atom. */
 396             velecsum         = _mm_add_ps(velecsum,velec);
 397
 398             fscal            = felec;
 399
 400              /* Update vectorial force */
 401             fix0             = _mm_macc_ps(dx02,fscal,fix0);
 402             fiy0             = _mm_macc_ps(dy02,fscal,fiy0);
 403             fiz0             = _mm_macc_ps(dz02,fscal,fiz0);
 404
 405             fjx2             = _mm_macc_ps(dx02,fscal,fjx2);
 406             fjy2             = _mm_macc_ps(dy02,fscal,fjy2);
 407             fjz2             = _mm_macc_ps(dz02,fscal,fjz2);
 408
 409             /**************************
 410              * CALCULATE INTERACTIONS *
 411              **************************/
 412
 413             /* REACTION-FIELD ELECTROSTATICS */
 414             velec            = _mm_mul_ps(qq10,_mm_sub_ps(_mm_macc_ps(krf,rsq10,rinv10),crf));
 415             felec            = _mm_mul_ps(qq10,_mm_msub_ps(rinv10,rinvsq10,krf2));
 416
 417             /* Update potential sum for this i atom from the interaction with this j atom. */
 418             velecsum         = _mm_add_ps(velecsum,velec);
 419
 420             fscal            = felec;
 421
 422              /* Update vectorial force */
 423             fix1             = _mm_macc_ps(dx10,fscal,fix1);
 424             fiy1             = _mm_macc_ps(dy10,fscal,fiy1);
 425             fiz1             = _mm_macc_ps(dz10,fscal,fiz1);
 426
 427             fjx0             = _mm_macc_ps(dx10,fscal,fjx0);
 428             fjy0             = _mm_macc_ps(dy10,fscal,fjy0);
 429             fjz0             = _mm_macc_ps(dz10,fscal,fjz0);
 430
 431             /**************************
 432              * CALCULATE INTERACTIONS *
 433              **************************/
 434
 435             /* REACTION-FIELD ELECTROSTATICS */
 436             velec            = _mm_mul_ps(qq11,_mm_sub_ps(_mm_macc_ps(krf,rsq11,rinv11),crf));
 437             felec            = _mm_mul_ps(qq11,_mm_msub_ps(rinv11,rinvsq11,krf2));
 438
 439             /* Update potential sum for this i atom from the interaction with this j atom. */
 440             velecsum         = _mm_add_ps(velecsum,velec);
 441
 442             fscal            = felec;
 443
 444              /* Update vectorial force */
 445             fix1             = _mm_macc_ps(dx11,fscal,fix1);
 446             fiy1             = _mm_macc_ps(dy11,fscal,fiy1);
 447             fiz1             = _mm_macc_ps(dz11,fscal,fiz1);
 448
 449             fjx1             = _mm_macc_ps(dx11,fscal,fjx1);
 450             fjy1             = _mm_macc_ps(dy11,fscal,fjy1);
 451             fjz1             = _mm_macc_ps(dz11,fscal,fjz1);
 452
 453             /**************************
 454              * CALCULATE INTERACTIONS *
 455              **************************/
 456
 457             /* REACTION-FIELD ELECTROSTATICS */
 458             velec            = _mm_mul_ps(qq12,_mm_sub_ps(_mm_macc_ps(krf,rsq12,rinv12),crf));
 459             felec            = _mm_mul_ps(qq12,_mm_msub_ps(rinv12,rinvsq12,krf2));
 460
 461             /* Update potential sum for this i atom from the interaction with this j atom. */
 462             velecsum         = _mm_add_ps(velecsum,velec);
 463
 464             fscal            = felec;
 465
 466              /* Update vectorial force */
 467             fix1             = _mm_macc_ps(dx12,fscal,fix1);
 468             fiy1             = _mm_macc_ps(dy12,fscal,fiy1);
 469             fiz1             = _mm_macc_ps(dz12,fscal,fiz1);
 470
 471             fjx2             = _mm_macc_ps(dx12,fscal,fjx2);
 472             fjy2             = _mm_macc_ps(dy12,fscal,fjy2);
 473             fjz2             = _mm_macc_ps(dz12,fscal,fjz2);
 474
 475             /**************************
 476              * CALCULATE INTERACTIONS *
 477              **************************/
 478
 479             /* REACTION-FIELD ELECTROSTATICS */
 480             velec            = _mm_mul_ps(qq20,_mm_sub_ps(_mm_macc_ps(krf,rsq20,rinv20),crf));
 481             felec            = _mm_mul_ps(qq20,_mm_msub_ps(rinv20,rinvsq20,krf2));
 482
 483             /* Update potential sum for this i atom from the interaction with this j atom. */
 484             velecsum         = _mm_add_ps(velecsum,velec);
 485
 486             fscal            = felec;
 487
 488              /* Update vectorial force */
 489             fix2             = _mm_macc_ps(dx20,fscal,fix2);
 490             fiy2             = _mm_macc_ps(dy20,fscal,fiy2);
 491             fiz2             = _mm_macc_ps(dz20,fscal,fiz2);
 492
 493             fjx0             = _mm_macc_ps(dx20,fscal,fjx0);
 494             fjy0             = _mm_macc_ps(dy20,fscal,fjy0);
 495             fjz0             = _mm_macc_ps(dz20,fscal,fjz0);
 496
 497             /**************************
 498              * CALCULATE INTERACTIONS *
 499              **************************/
 500
 501             /* REACTION-FIELD ELECTROSTATICS */
 502             velec            = _mm_mul_ps(qq21,_mm_sub_ps(_mm_macc_ps(krf,rsq21,rinv21),crf));
 503             felec            = _mm_mul_ps(qq21,_mm_msub_ps(rinv21,rinvsq21,krf2));
 504
 505             /* Update potential sum for this i atom from the interaction with this j atom. */
 506             velecsum         = _mm_add_ps(velecsum,velec);
 507
 508             fscal            = felec;
 509
 510              /* Update vectorial force */
 511             fix2             = _mm_macc_ps(dx21,fscal,fix2);
 512             fiy2             = _mm_macc_ps(dy21,fscal,fiy2);
 513             fiz2             = _mm_macc_ps(dz21,fscal,fiz2);
 514
 515             fjx1             = _mm_macc_ps(dx21,fscal,fjx1);
 516             fjy1             = _mm_macc_ps(dy21,fscal,fjy1);
 517             fjz1             = _mm_macc_ps(dz21,fscal,fjz1);
 518
 519             /**************************
 520              * CALCULATE INTERACTIONS *
 521              **************************/
 522
 523             /* REACTION-FIELD ELECTROSTATICS */
 524             velec            = _mm_mul_ps(qq22,_mm_sub_ps(_mm_macc_ps(krf,rsq22,rinv22),crf));
 525             felec            = _mm_mul_ps(qq22,_mm_msub_ps(rinv22,rinvsq22,krf2));
 526
 527             /* Update potential sum for this i atom from the interaction with this j atom. */
 528             velecsum         = _mm_add_ps(velecsum,velec);
 529
 530             fscal            = felec;
 531
 532              /* Update vectorial force */
 533             fix2             = _mm_macc_ps(dx22,fscal,fix2);
 534             fiy2             = _mm_macc_ps(dy22,fscal,fiy2);
 535             fiz2             = _mm_macc_ps(dz22,fscal,fiz2);
 536
 537             fjx2             = _mm_macc_ps(dx22,fscal,fjx2);
 538             fjy2             = _mm_macc_ps(dy22,fscal,fjy2);
 539             fjz2             = _mm_macc_ps(dz22,fscal,fjz2);
 540
 541             fjptrA             = f+j_coord_offsetA;
 542             fjptrB             = f+j_coord_offsetB;
 543             fjptrC             = f+j_coord_offsetC;
 544             fjptrD             = f+j_coord_offsetD;
 545
 546             gmx_mm_decrement_3rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,
 547                                                    fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
 548
 549             /* Inner loop uses 350 flops */
 550         }
 551
 552         if(jidx<j_index_end)
 553         {
 554
 555             /* Get j neighbor index, and coordinate index */
 556             jnrlistA         = jjnr[jidx];
 557             jnrlistB         = jjnr[jidx+1];
 558             jnrlistC         = jjnr[jidx+2];
 559             jnrlistD         = jjnr[jidx+3];
 560             /* Sign of each element will be negative for non-real atoms.
 561              * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
 562              * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
 563              */
 564             dummy_mask = gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
 565             jnrA       = (jnrlistA>=0) ? jnrlistA : 0;
 566             jnrB       = (jnrlistB>=0) ? jnrlistB : 0;
 567             jnrC       = (jnrlistC>=0) ? jnrlistC : 0;
 568             jnrD       = (jnrlistD>=0) ? jnrlistD : 0;
 569             j_coord_offsetA  = DIM*jnrA;
 570             j_coord_offsetB  = DIM*jnrB;
 571             j_coord_offsetC  = DIM*jnrC;
 572             j_coord_offsetD  = DIM*jnrD;
 573
 574             /* load j atom coordinates */
 575             gmx_mm_load_3rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
 576                                               x+j_coord_offsetC,x+j_coord_offsetD,
 577                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
 578
 579             /* Calculate displacement vector */
 580             dx00             = _mm_sub_ps(ix0,jx0);
 581             dy00             = _mm_sub_ps(iy0,jy0);
 582             dz00             = _mm_sub_ps(iz0,jz0);
 583             dx01             = _mm_sub_ps(ix0,jx1);
 584             dy01             = _mm_sub_ps(iy0,jy1);
 585             dz01             = _mm_sub_ps(iz0,jz1);
 586             dx02             = _mm_sub_ps(ix0,jx2);
 587             dy02             = _mm_sub_ps(iy0,jy2);
 588             dz02             = _mm_sub_ps(iz0,jz2);
 589             dx10             = _mm_sub_ps(ix1,jx0);
 590             dy10             = _mm_sub_ps(iy1,jy0);
 591             dz10             = _mm_sub_ps(iz1,jz0);
 592             dx11             = _mm_sub_ps(ix1,jx1);
 593             dy11             = _mm_sub_ps(iy1,jy1);
 594             dz11             = _mm_sub_ps(iz1,jz1);
 595             dx12             = _mm_sub_ps(ix1,jx2);
 596             dy12             = _mm_sub_ps(iy1,jy2);
 597             dz12             = _mm_sub_ps(iz1,jz2);
 598             dx20             = _mm_sub_ps(ix2,jx0);
 599             dy20             = _mm_sub_ps(iy2,jy0);
 600             dz20             = _mm_sub_ps(iz2,jz0);
 601             dx21             = _mm_sub_ps(ix2,jx1);
 602             dy21             = _mm_sub_ps(iy2,jy1);
 603             dz21             = _mm_sub_ps(iz2,jz1);
 604             dx22             = _mm_sub_ps(ix2,jx2);
 605             dy22             = _mm_sub_ps(iy2,jy2);
 606             dz22             = _mm_sub_ps(iz2,jz2);
 607
 608             /* Calculate squared distance and things based on it */
 609             rsq00            = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
 610             rsq01            = gmx_mm_calc_rsq_ps(dx01,dy01,dz01);
 611             rsq02            = gmx_mm_calc_rsq_ps(dx02,dy02,dz02);
 612             rsq10            = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
 613             rsq11            = gmx_mm_calc_rsq_ps(dx11,dy11,dz11);
 614             rsq12            = gmx_mm_calc_rsq_ps(dx12,dy12,dz12);
 615             rsq20            = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
 616             rsq21            = gmx_mm_calc_rsq_ps(dx21,dy21,dz21);
 617             rsq22            = gmx_mm_calc_rsq_ps(dx22,dy22,dz22);
 618
 619             rinv00           = gmx_mm_invsqrt_ps(rsq00);
 620             rinv01           = gmx_mm_invsqrt_ps(rsq01);
 621             rinv02           = gmx_mm_invsqrt_ps(rsq02);
 622             rinv10           = gmx_mm_invsqrt_ps(rsq10);
 623             rinv11           = gmx_mm_invsqrt_ps(rsq11);
 624             rinv12           = gmx_mm_invsqrt_ps(rsq12);
 625             rinv20           = gmx_mm_invsqrt_ps(rsq20);
 626             rinv21           = gmx_mm_invsqrt_ps(rsq21);
 627             rinv22           = gmx_mm_invsqrt_ps(rsq22);
 628
 629             rinvsq00         = _mm_mul_ps(rinv00,rinv00);
 630             rinvsq01         = _mm_mul_ps(rinv01,rinv01);
 631             rinvsq02         = _mm_mul_ps(rinv02,rinv02);
 632             rinvsq10         = _mm_mul_ps(rinv10,rinv10);
 633             rinvsq11         = _mm_mul_ps(rinv11,rinv11);
 634             rinvsq12         = _mm_mul_ps(rinv12,rinv12);
 635             rinvsq20         = _mm_mul_ps(rinv20,rinv20);
 636             rinvsq21         = _mm_mul_ps(rinv21,rinv21);
 637             rinvsq22         = _mm_mul_ps(rinv22,rinv22);
 638
 639             fjx0             = _mm_setzero_ps();
 640             fjy0             = _mm_setzero_ps();
 641             fjz0             = _mm_setzero_ps();
 642             fjx1             = _mm_setzero_ps();
 643             fjy1             = _mm_setzero_ps();
 644             fjz1             = _mm_setzero_ps();
 645             fjx2             = _mm_setzero_ps();
 646             fjy2             = _mm_setzero_ps();
 647             fjz2             = _mm_setzero_ps();
 648
 649             /**************************
 650              * CALCULATE INTERACTIONS *
 651              **************************/
 652
 653             r00              = _mm_mul_ps(rsq00,rinv00);
 654             r00              = _mm_andnot_ps(dummy_mask,r00);
 655
 656             /* Calculate table index by multiplying r with table scale and truncate to integer */
 657             rt               = _mm_mul_ps(r00,vftabscale);
 658             vfitab           = _mm_cvttps_epi32(rt);
 659 #ifdef __XOP__
 660             vfeps            = _mm_frcz_ps(rt);
 661 #else
 662             vfeps            = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
 663 #endif
 664             twovfeps         = _mm_add_ps(vfeps,vfeps);
 665             vfitab           = _mm_slli_epi32(vfitab,3);
 666
 667             /* REACTION-FIELD ELECTROSTATICS */
 668             velec            = _mm_mul_ps(qq00,_mm_sub_ps(_mm_macc_ps(krf,rsq00,rinv00),crf));
 669             felec            = _mm_mul_ps(qq00,_mm_msub_ps(rinv00,rinvsq00,krf2));
 670
 671             /* CUBIC SPLINE TABLE DISPERSION */
 672             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
 673             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
 674             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
 675             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
 676             _MM_TRANSPOSE4_PS(Y,F,G,H);
 677             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
 678             VV               = _mm_macc_ps(vfeps,Fp,Y);
 679             vvdw6            = _mm_mul_ps(c6_00,VV);
 680             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
 681             fvdw6            = _mm_mul_ps(c6_00,FF);
 682
 683             /* CUBIC SPLINE TABLE REPULSION */
 684             vfitab           = _mm_add_epi32(vfitab,ifour);
 685             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
 686             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
 687             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
 688             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
 689             _MM_TRANSPOSE4_PS(Y,F,G,H);
 690             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
 691             VV               = _mm_macc_ps(vfeps,Fp,Y);
 692             vvdw12           = _mm_mul_ps(c12_00,VV);
 693             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
 694             fvdw12           = _mm_mul_ps(c12_00,FF);
 695             vvdw             = _mm_add_ps(vvdw12,vvdw6);
 696             fvdw             = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
 697
 698             /* Update potential sum for this i atom from the interaction with this j atom. */
 699             velec            = _mm_andnot_ps(dummy_mask,velec);
 700             velecsum         = _mm_add_ps(velecsum,velec);
 701             vvdw             = _mm_andnot_ps(dummy_mask,vvdw);
 702             vvdwsum          = _mm_add_ps(vvdwsum,vvdw);
 703
 704             fscal            = _mm_add_ps(felec,fvdw);
 705
 706             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 707
 708              /* Update vectorial force */
 709             fix0             = _mm_macc_ps(dx00,fscal,fix0);
 710             fiy0             = _mm_macc_ps(dy00,fscal,fiy0);
 711             fiz0             = _mm_macc_ps(dz00,fscal,fiz0);
 712
 713             fjx0             = _mm_macc_ps(dx00,fscal,fjx0);
 714             fjy0             = _mm_macc_ps(dy00,fscal,fjy0);
 715             fjz0             = _mm_macc_ps(dz00,fscal,fjz0);
 716
 717             /**************************
 718              * CALCULATE INTERACTIONS *
 719              **************************/
 720
 721             /* REACTION-FIELD ELECTROSTATICS */
 722             velec            = _mm_mul_ps(qq01,_mm_sub_ps(_mm_macc_ps(krf,rsq01,rinv01),crf));
 723             felec            = _mm_mul_ps(qq01,_mm_msub_ps(rinv01,rinvsq01,krf2));
 724
 725             /* Update potential sum for this i atom from the interaction with this j atom. */
 726             velec            = _mm_andnot_ps(dummy_mask,velec);
 727             velecsum         = _mm_add_ps(velecsum,velec);
 728
 729             fscal            = felec;
 730
 731             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 732
 733              /* Update vectorial force */
 734             fix0             = _mm_macc_ps(dx01,fscal,fix0);
 735             fiy0             = _mm_macc_ps(dy01,fscal,fiy0);
 736             fiz0             = _mm_macc_ps(dz01,fscal,fiz0);
 737
 738             fjx1             = _mm_macc_ps(dx01,fscal,fjx1);
 739             fjy1             = _mm_macc_ps(dy01,fscal,fjy1);
 740             fjz1             = _mm_macc_ps(dz01,fscal,fjz1);
 741
 742             /**************************
 743              * CALCULATE INTERACTIONS *
 744              **************************/
 745
 746             /* REACTION-FIELD ELECTROSTATICS */
 747             velec            = _mm_mul_ps(qq02,_mm_sub_ps(_mm_macc_ps(krf,rsq02,rinv02),crf));
 748             felec            = _mm_mul_ps(qq02,_mm_msub_ps(rinv02,rinvsq02,krf2));
 749
 750             /* Update potential sum for this i atom from the interaction with this j atom. */
 751             velec            = _mm_andnot_ps(dummy_mask,velec);
 752             velecsum         = _mm_add_ps(velecsum,velec);
 753
 754             fscal            = felec;
 755
 756             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 757
 758              /* Update vectorial force */
 759             fix0             = _mm_macc_ps(dx02,fscal,fix0);
 760             fiy0             = _mm_macc_ps(dy02,fscal,fiy0);
 761             fiz0             = _mm_macc_ps(dz02,fscal,fiz0);
 762
 763             fjx2             = _mm_macc_ps(dx02,fscal,fjx2);
 764             fjy2             = _mm_macc_ps(dy02,fscal,fjy2);
 765             fjz2             = _mm_macc_ps(dz02,fscal,fjz2);
 766
 767             /**************************
 768              * CALCULATE INTERACTIONS *
 769              **************************/
 770
 771             /* REACTION-FIELD ELECTROSTATICS */
 772             velec            = _mm_mul_ps(qq10,_mm_sub_ps(_mm_macc_ps(krf,rsq10,rinv10),crf));
 773             felec            = _mm_mul_ps(qq10,_mm_msub_ps(rinv10,rinvsq10,krf2));
 774
 775             /* Update potential sum for this i atom from the interaction with this j atom. */
 776             velec            = _mm_andnot_ps(dummy_mask,velec);
 777             velecsum         = _mm_add_ps(velecsum,velec);
 778
 779             fscal            = felec;
 780
 781             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 782
 783              /* Update vectorial force */
 784             fix1             = _mm_macc_ps(dx10,fscal,fix1);
 785             fiy1             = _mm_macc_ps(dy10,fscal,fiy1);
 786             fiz1             = _mm_macc_ps(dz10,fscal,fiz1);
 787
 788             fjx0             = _mm_macc_ps(dx10,fscal,fjx0);
 789             fjy0             = _mm_macc_ps(dy10,fscal,fjy0);
 790             fjz0             = _mm_macc_ps(dz10,fscal,fjz0);
 791
 792             /**************************
 793              * CALCULATE INTERACTIONS *
 794              **************************/
 795
 796             /* REACTION-FIELD ELECTROSTATICS */
 797             velec            = _mm_mul_ps(qq11,_mm_sub_ps(_mm_macc_ps(krf,rsq11,rinv11),crf));
 798             felec            = _mm_mul_ps(qq11,_mm_msub_ps(rinv11,rinvsq11,krf2));
 799
 800             /* Update potential sum for this i atom from the interaction with this j atom. */
 801             velec            = _mm_andnot_ps(dummy_mask,velec);
 802             velecsum         = _mm_add_ps(velecsum,velec);
 803
 804             fscal            = felec;
 805
 806             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 807
 808              /* Update vectorial force */
 809             fix1             = _mm_macc_ps(dx11,fscal,fix1);
 810             fiy1             = _mm_macc_ps(dy11,fscal,fiy1);
 811             fiz1             = _mm_macc_ps(dz11,fscal,fiz1);
 812
 813             fjx1             = _mm_macc_ps(dx11,fscal,fjx1);
 814             fjy1             = _mm_macc_ps(dy11,fscal,fjy1);
 815             fjz1             = _mm_macc_ps(dz11,fscal,fjz1);
 816
 817             /**************************
 818              * CALCULATE INTERACTIONS *
 819              **************************/
 820
 821             /* REACTION-FIELD ELECTROSTATICS */
 822             velec            = _mm_mul_ps(qq12,_mm_sub_ps(_mm_macc_ps(krf,rsq12,rinv12),crf));
 823             felec            = _mm_mul_ps(qq12,_mm_msub_ps(rinv12,rinvsq12,krf2));
 824
 825             /* Update potential sum for this i atom from the interaction with this j atom. */
 826             velec            = _mm_andnot_ps(dummy_mask,velec);
 827             velecsum         = _mm_add_ps(velecsum,velec);
 828
 829             fscal            = felec;
 830
 831             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 832
 833              /* Update vectorial force */
 834             fix1             = _mm_macc_ps(dx12,fscal,fix1);
 835             fiy1             = _mm_macc_ps(dy12,fscal,fiy1);
 836             fiz1             = _mm_macc_ps(dz12,fscal,fiz1);
 837
 838             fjx2             = _mm_macc_ps(dx12,fscal,fjx2);
 839             fjy2             = _mm_macc_ps(dy12,fscal,fjy2);
 840             fjz2             = _mm_macc_ps(dz12,fscal,fjz2);
 841
 842             /**************************
 843              * CALCULATE INTERACTIONS *
 844              **************************/
 845
 846             /* REACTION-FIELD ELECTROSTATICS */
 847             velec            = _mm_mul_ps(qq20,_mm_sub_ps(_mm_macc_ps(krf,rsq20,rinv20),crf));
 848             felec            = _mm_mul_ps(qq20,_mm_msub_ps(rinv20,rinvsq20,krf2));
 849
 850             /* Update potential sum for this i atom from the interaction with this j atom. */
 851             velec            = _mm_andnot_ps(dummy_mask,velec);
 852             velecsum         = _mm_add_ps(velecsum,velec);
 853
 854             fscal            = felec;
 855
 856             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 857
 858              /* Update vectorial force */
 859             fix2             = _mm_macc_ps(dx20,fscal,fix2);
 860             fiy2             = _mm_macc_ps(dy20,fscal,fiy2);
 861             fiz2             = _mm_macc_ps(dz20,fscal,fiz2);
 862
 863             fjx0             = _mm_macc_ps(dx20,fscal,fjx0);
 864             fjy0             = _mm_macc_ps(dy20,fscal,fjy0);
 865             fjz0             = _mm_macc_ps(dz20,fscal,fjz0);
 866
 867             /**************************
 868              * CALCULATE INTERACTIONS *
 869              **************************/
 870
 871             /* REACTION-FIELD ELECTROSTATICS */
 872             velec            = _mm_mul_ps(qq21,_mm_sub_ps(_mm_macc_ps(krf,rsq21,rinv21),crf));
 873             felec            = _mm_mul_ps(qq21,_mm_msub_ps(rinv21,rinvsq21,krf2));
 874
 875             /* Update potential sum for this i atom from the interaction with this j atom. */
 876             velec            = _mm_andnot_ps(dummy_mask,velec);
 877             velecsum         = _mm_add_ps(velecsum,velec);
 878
 879             fscal            = felec;
 880
 881             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 882
 883              /* Update vectorial force */
 884             fix2             = _mm_macc_ps(dx21,fscal,fix2);
 885             fiy2             = _mm_macc_ps(dy21,fscal,fiy2);
 886             fiz2             = _mm_macc_ps(dz21,fscal,fiz2);
 887
 888             fjx1             = _mm_macc_ps(dx21,fscal,fjx1);
 889             fjy1             = _mm_macc_ps(dy21,fscal,fjy1);
 890             fjz1             = _mm_macc_ps(dz21,fscal,fjz1);
 891
 892             /**************************
 893              * CALCULATE INTERACTIONS *
 894              **************************/
 895
 896             /* REACTION-FIELD ELECTROSTATICS */
 897             velec            = _mm_mul_ps(qq22,_mm_sub_ps(_mm_macc_ps(krf,rsq22,rinv22),crf));
 898             felec            = _mm_mul_ps(qq22,_mm_msub_ps(rinv22,rinvsq22,krf2));
 899
 900             /* Update potential sum for this i atom from the interaction with this j atom. */
 901             velec            = _mm_andnot_ps(dummy_mask,velec);
 902             velecsum         = _mm_add_ps(velecsum,velec);
 903
 904             fscal            = felec;
 905
 906             fscal            = _mm_andnot_ps(dummy_mask,fscal);
 907
 908              /* Update vectorial force */
 909             fix2             = _mm_macc_ps(dx22,fscal,fix2);
 910             fiy2             = _mm_macc_ps(dy22,fscal,fiy2);
 911             fiz2             = _mm_macc_ps(dz22,fscal,fiz2);
 912
 913             fjx2             = _mm_macc_ps(dx22,fscal,fjx2);
 914             fjy2             = _mm_macc_ps(dy22,fscal,fjy2);
 915             fjz2             = _mm_macc_ps(dz22,fscal,fjz2);
 916
 917             fjptrA             = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
 918             fjptrB             = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
 919             fjptrC             = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
 920             fjptrD             = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
 921
 922             gmx_mm_decrement_3rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,
 923                                                    fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
 924
 925             /* Inner loop uses 351 flops */
 926         }
 927
 928         /* End of innermost loop */
 929
 930         gmx_mm_update_iforce_3atom_swizzle_ps(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,
 931                                               f+i_coord_offset,fshift+i_shift_offset);
 932
 933         ggid                        = gid[iidx];
 934         /* Update potential energies */
 935         gmx_mm_update_1pot_ps(velecsum,kernel_data->energygrp_elec+ggid);
 936         gmx_mm_update_1pot_ps(vvdwsum,kernel_data->energygrp_vdw+ggid);
 937
 938         /* Increment number of inner iterations */
 939         inneriter                  += j_index_end - j_index_start;
 940
 941         /* Outer loop uses 20 flops */
 942     }
 943
 944     /* Increment number of outer iterations */
 945     outeriter        += nri;
 946
 947     /* Update outer/inner flops */
 948
 949     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W3W3_VF,outeriter*20 + inneriter*351);
 950 }
 951 /*
 952  * Gromacs nonbonded kernel:   nb_kernel_ElecRF_VdwCSTab_GeomW3W3_F_avx_128_fma_single
 953  * Electrostatics interaction: ReactionField
 954  * VdW interaction:            CubicSplineTable
 955  * Geometry:                   Water3-Water3
 956  * Calculate force/pot:        Force
 957  */
 958 void
 959 nb_kernel_ElecRF_VdwCSTab_GeomW3W3_F_avx_128_fma_single
 960                     (t_nblist                    * gmx_restrict       nlist,
 961                      rvec                        * gmx_restrict          xx,
 962                      rvec                        * gmx_restrict          ff,
 963                      t_forcerec                  * gmx_restrict          fr,
 964                      t_mdatoms                   * gmx_restrict     mdatoms,
 965                      nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
 966                      t_nrnb                      * gmx_restrict        nrnb)
 967 {
 968     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
 969      * just 0 for non-waters.
 970      * Suffixes A,B,C,D refer to j loop unrolling done with AVX_128, e.g. for the four different
 971      * jnr indices corresponding to data put in the four positions in the SIMD register.
 972      */
 973     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
 974     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
 975     int              jnrA,jnrB,jnrC,jnrD;
 976     int              jnrlistA,jnrlistB,jnrlistC,jnrlistD;
 977     int              j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
 978     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
 979     real             rcutoff_scalar;
 980     real             *shiftvec,*fshift,*x,*f;
 981     real             *fjptrA,*fjptrB,*fjptrC,*fjptrD;
 982     real             scratch[4*DIM];
 983     __m128           fscal,rcutoff,rcutoff2,jidxall;
 984     int              vdwioffset0;
 985     __m128           ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
 986     int              vdwioffset1;
 987     __m128           ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
 988     int              vdwioffset2;
 989     __m128           ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
 990     int              vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
 991     __m128           jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
 992     int              vdwjidx1A,vdwjidx1B,vdwjidx1C,vdwjidx1D;
 993     __m128           jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
 994     int              vdwjidx2A,vdwjidx2B,vdwjidx2C,vdwjidx2D;
 995     __m128           jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
 996     __m128           dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
 997     __m128           dx01,dy01,dz01,rsq01,rinv01,rinvsq01,r01,qq01,c6_01,c12_01;
 998     __m128           dx02,dy02,dz02,rsq02,rinv02,rinvsq02,r02,qq02,c6_02,c12_02;
 999     __m128           dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
1000     __m128           dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
1001     __m128           dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
1002     __m128           dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
1003     __m128           dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
1004     __m128           dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
1005     __m128           velec,felec,velecsum,facel,crf,krf,krf2;
1006     real             *charge;
1007     int              nvdwtype;
1008     __m128           rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
1009     int              *vdwtype;
1010     real             *vdwparam;
1011     __m128           one_sixth   = _mm_set1_ps(1.0/6.0);
1012     __m128           one_twelfth = _mm_set1_ps(1.0/12.0);
1013     __m128i          vfitab;
1014     __m128i          ifour       = _mm_set1_epi32(4);
1015     __m128           rt,vfeps,twovfeps,vftabscale,Y,F,G,H,Fp,VV,FF;
1016     real             *vftab;
1017     __m128           dummy_mask,cutoff_mask;
1018     __m128           signbit = _mm_castsi128_ps( _mm_set1_epi32(0x80000000) );
1019     __m128           one     = _mm_set1_ps(1.0);
1020     __m128           two     = _mm_set1_ps(2.0);
1021     x                = xx[0];
1022     f                = ff[0];
1023
1024     nri              = nlist->nri;
1025     iinr             = nlist->iinr;
1026     jindex           = nlist->jindex;
1027     jjnr             = nlist->jjnr;
1028     shiftidx         = nlist->shift;
1029     gid              = nlist->gid;
1030     shiftvec         = fr->shift_vec[0];
1031     fshift           = fr->fshift[0];
1032     facel            = _mm_set1_ps(fr->epsfac);
1033     charge           = mdatoms->chargeA;
1034     krf              = _mm_set1_ps(fr->ic->k_rf);
1035     krf2             = _mm_set1_ps(fr->ic->k_rf*2.0);
1036     crf              = _mm_set1_ps(fr->ic->c_rf);
1037     nvdwtype         = fr->ntype;
1038     vdwparam         = fr->nbfp;
1039     vdwtype          = mdatoms->typeA;
1040
1041     vftab            = kernel_data->table_vdw->data;
1042     vftabscale       = _mm_set1_ps(kernel_data->table_vdw->scale);
1043
1044     /* Setup water-specific parameters */
1045     inr              = nlist->iinr[0];
1046     iq0              = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+0]));
1047     iq1              = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+1]));
1048     iq2              = _mm_mul_ps(facel,_mm_set1_ps(charge[inr+2]));
1049     vdwioffset0      = 2*nvdwtype*vdwtype[inr+0];
1050
1051     jq0              = _mm_set1_ps(charge[inr+0]);
1052     jq1              = _mm_set1_ps(charge[inr+1]);
1053     jq2              = _mm_set1_ps(charge[inr+2]);
1054     vdwjidx0A        = 2*vdwtype[inr+0];
1055     qq00             = _mm_mul_ps(iq0,jq0);
1056     c6_00            = _mm_set1_ps(vdwparam[vdwioffset0+vdwjidx0A]);
1057     c12_00           = _mm_set1_ps(vdwparam[vdwioffset0+vdwjidx0A+1]);
1058     qq01             = _mm_mul_ps(iq0,jq1);
1059     qq02             = _mm_mul_ps(iq0,jq2);
1060     qq10             = _mm_mul_ps(iq1,jq0);
1061     qq11             = _mm_mul_ps(iq1,jq1);
1062     qq12             = _mm_mul_ps(iq1,jq2);
1063     qq20             = _mm_mul_ps(iq2,jq0);
1064     qq21             = _mm_mul_ps(iq2,jq1);
1065     qq22             = _mm_mul_ps(iq2,jq2);
1066
1067     /* Avoid stupid compiler warnings */
1068     jnrA = jnrB = jnrC = jnrD = 0;
1069     j_coord_offsetA = 0;
1070     j_coord_offsetB = 0;
1071     j_coord_offsetC = 0;
1072     j_coord_offsetD = 0;
1073
1074     outeriter        = 0;
1075     inneriter        = 0;
1076
1077     for(iidx=0;iidx<4*DIM;iidx++)
1078     {
1079         scratch[iidx] = 0.0;
1080     }
1081
1082     /* Start outer loop over neighborlists */
1083     for(iidx=0; iidx<nri; iidx++)
1084     {
1085         /* Load shift vector for this list */
1086         i_shift_offset   = DIM*shiftidx[iidx];
1087
1088         /* Load limits for loop over neighbors */
1089         j_index_start    = jindex[iidx];
1090         j_index_end      = jindex[iidx+1];
1091
1092         /* Get outer coordinate index */
1093         inr              = iinr[iidx];
1094         i_coord_offset   = DIM*inr;
1095
1096         /* Load i particle coords and add shift vector */
1097         gmx_mm_load_shift_and_3rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,
1098                                                  &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2);
1099
1100         fix0             = _mm_setzero_ps();
1101         fiy0             = _mm_setzero_ps();
1102         fiz0             = _mm_setzero_ps();
1103         fix1             = _mm_setzero_ps();
1104         fiy1             = _mm_setzero_ps();
1105         fiz1             = _mm_setzero_ps();
1106         fix2             = _mm_setzero_ps();
1107         fiy2             = _mm_setzero_ps();
1108         fiz2             = _mm_setzero_ps();
1109
1110         /* Start inner kernel loop */
1111         for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
1112         {
1113
1114             /* Get j neighbor index, and coordinate index */
1115             jnrA             = jjnr[jidx];
1116             jnrB             = jjnr[jidx+1];
1117             jnrC             = jjnr[jidx+2];
1118             jnrD             = jjnr[jidx+3];
1119             j_coord_offsetA  = DIM*jnrA;
1120             j_coord_offsetB  = DIM*jnrB;
1121             j_coord_offsetC  = DIM*jnrC;
1122             j_coord_offsetD  = DIM*jnrD;
1123
1124             /* load j atom coordinates */
1125             gmx_mm_load_3rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
1126                                               x+j_coord_offsetC,x+j_coord_offsetD,
1127                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
1128
1129             /* Calculate displacement vector */
1130             dx00             = _mm_sub_ps(ix0,jx0);
1131             dy00             = _mm_sub_ps(iy0,jy0);
1132             dz00             = _mm_sub_ps(iz0,jz0);
1133             dx01             = _mm_sub_ps(ix0,jx1);
1134             dy01             = _mm_sub_ps(iy0,jy1);
1135             dz01             = _mm_sub_ps(iz0,jz1);
1136             dx02             = _mm_sub_ps(ix0,jx2);
1137             dy02             = _mm_sub_ps(iy0,jy2);
1138             dz02             = _mm_sub_ps(iz0,jz2);
1139             dx10             = _mm_sub_ps(ix1,jx0);
1140             dy10             = _mm_sub_ps(iy1,jy0);
1141             dz10             = _mm_sub_ps(iz1,jz0);
1142             dx11             = _mm_sub_ps(ix1,jx1);
1143             dy11             = _mm_sub_ps(iy1,jy1);
1144             dz11             = _mm_sub_ps(iz1,jz1);
1145             dx12             = _mm_sub_ps(ix1,jx2);
1146             dy12             = _mm_sub_ps(iy1,jy2);
1147             dz12             = _mm_sub_ps(iz1,jz2);
1148             dx20             = _mm_sub_ps(ix2,jx0);
1149             dy20             = _mm_sub_ps(iy2,jy0);
1150             dz20             = _mm_sub_ps(iz2,jz0);
1151             dx21             = _mm_sub_ps(ix2,jx1);
1152             dy21             = _mm_sub_ps(iy2,jy1);
1153             dz21             = _mm_sub_ps(iz2,jz1);
1154             dx22             = _mm_sub_ps(ix2,jx2);
1155             dy22             = _mm_sub_ps(iy2,jy2);
1156             dz22             = _mm_sub_ps(iz2,jz2);
1157
1158             /* Calculate squared distance and things based on it */
1159             rsq00            = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
1160             rsq01            = gmx_mm_calc_rsq_ps(dx01,dy01,dz01);
1161             rsq02            = gmx_mm_calc_rsq_ps(dx02,dy02,dz02);
1162             rsq10            = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
1163             rsq11            = gmx_mm_calc_rsq_ps(dx11,dy11,dz11);
1164             rsq12            = gmx_mm_calc_rsq_ps(dx12,dy12,dz12);
1165             rsq20            = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
1166             rsq21            = gmx_mm_calc_rsq_ps(dx21,dy21,dz21);
1167             rsq22            = gmx_mm_calc_rsq_ps(dx22,dy22,dz22);
1168
1169             rinv00           = gmx_mm_invsqrt_ps(rsq00);
1170             rinv01           = gmx_mm_invsqrt_ps(rsq01);
1171             rinv02           = gmx_mm_invsqrt_ps(rsq02);
1172             rinv10           = gmx_mm_invsqrt_ps(rsq10);
1173             rinv11           = gmx_mm_invsqrt_ps(rsq11);
1174             rinv12           = gmx_mm_invsqrt_ps(rsq12);
1175             rinv20           = gmx_mm_invsqrt_ps(rsq20);
1176             rinv21           = gmx_mm_invsqrt_ps(rsq21);
1177             rinv22           = gmx_mm_invsqrt_ps(rsq22);
1178
1179             rinvsq00         = _mm_mul_ps(rinv00,rinv00);
1180             rinvsq01         = _mm_mul_ps(rinv01,rinv01);
1181             rinvsq02         = _mm_mul_ps(rinv02,rinv02);
1182             rinvsq10         = _mm_mul_ps(rinv10,rinv10);
1183             rinvsq11         = _mm_mul_ps(rinv11,rinv11);
1184             rinvsq12         = _mm_mul_ps(rinv12,rinv12);
1185             rinvsq20         = _mm_mul_ps(rinv20,rinv20);
1186             rinvsq21         = _mm_mul_ps(rinv21,rinv21);
1187             rinvsq22         = _mm_mul_ps(rinv22,rinv22);
1188
1189             fjx0             = _mm_setzero_ps();
1190             fjy0             = _mm_setzero_ps();
1191             fjz0             = _mm_setzero_ps();
1192             fjx1             = _mm_setzero_ps();
1193             fjy1             = _mm_setzero_ps();
1194             fjz1             = _mm_setzero_ps();
1195             fjx2             = _mm_setzero_ps();
1196             fjy2             = _mm_setzero_ps();
1197             fjz2             = _mm_setzero_ps();
1198
1199             /**************************
1200              * CALCULATE INTERACTIONS *
1201              **************************/
1202
1203             r00              = _mm_mul_ps(rsq00,rinv00);
1204
1205             /* Calculate table index by multiplying r with table scale and truncate to integer */
1206             rt               = _mm_mul_ps(r00,vftabscale);
1207             vfitab           = _mm_cvttps_epi32(rt);
1208 #ifdef __XOP__
1209             vfeps            = _mm_frcz_ps(rt);
1210 #else
1211             vfeps            = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
1212 #endif
1213             twovfeps         = _mm_add_ps(vfeps,vfeps);
1214             vfitab           = _mm_slli_epi32(vfitab,3);
1215
1216             /* REACTION-FIELD ELECTROSTATICS */
1217             felec            = _mm_mul_ps(qq00,_mm_msub_ps(rinv00,rinvsq00,krf2));
1218
1219             /* CUBIC SPLINE TABLE DISPERSION */
1220             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
1221             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
1222             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
1223             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
1224             _MM_TRANSPOSE4_PS(Y,F,G,H);
1225             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
1226             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
1227             fvdw6            = _mm_mul_ps(c6_00,FF);
1228
1229             /* CUBIC SPLINE TABLE REPULSION */
1230             vfitab           = _mm_add_epi32(vfitab,ifour);
1231             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
1232             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
1233             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
1234             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
1235             _MM_TRANSPOSE4_PS(Y,F,G,H);
1236             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
1237             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
1238             fvdw12           = _mm_mul_ps(c12_00,FF);
1239             fvdw             = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
1240
1241             fscal            = _mm_add_ps(felec,fvdw);
1242
1243              /* Update vectorial force */
1244             fix0             = _mm_macc_ps(dx00,fscal,fix0);
1245             fiy0             = _mm_macc_ps(dy00,fscal,fiy0);
1246             fiz0             = _mm_macc_ps(dz00,fscal,fiz0);
1247
1248             fjx0             = _mm_macc_ps(dx00,fscal,fjx0);
1249             fjy0             = _mm_macc_ps(dy00,fscal,fjy0);
1250             fjz0             = _mm_macc_ps(dz00,fscal,fjz0);
1251
1252             /**************************
1253              * CALCULATE INTERACTIONS *
1254              **************************/
1255
1256             /* REACTION-FIELD ELECTROSTATICS */
1257             felec            = _mm_mul_ps(qq01,_mm_msub_ps(rinv01,rinvsq01,krf2));
1258
1259             fscal            = felec;
1260
1261              /* Update vectorial force */
1262             fix0             = _mm_macc_ps(dx01,fscal,fix0);
1263             fiy0             = _mm_macc_ps(dy01,fscal,fiy0);
1264             fiz0             = _mm_macc_ps(dz01,fscal,fiz0);
1265
1266             fjx1             = _mm_macc_ps(dx01,fscal,fjx1);
1267             fjy1             = _mm_macc_ps(dy01,fscal,fjy1);
1268             fjz1             = _mm_macc_ps(dz01,fscal,fjz1);
1269
1270             /**************************
1271              * CALCULATE INTERACTIONS *
1272              **************************/
1273
1274             /* REACTION-FIELD ELECTROSTATICS */
1275             felec            = _mm_mul_ps(qq02,_mm_msub_ps(rinv02,rinvsq02,krf2));
1276
1277             fscal            = felec;
1278
1279              /* Update vectorial force */
1280             fix0             = _mm_macc_ps(dx02,fscal,fix0);
1281             fiy0             = _mm_macc_ps(dy02,fscal,fiy0);
1282             fiz0             = _mm_macc_ps(dz02,fscal,fiz0);
1283
1284             fjx2             = _mm_macc_ps(dx02,fscal,fjx2);
1285             fjy2             = _mm_macc_ps(dy02,fscal,fjy2);
1286             fjz2             = _mm_macc_ps(dz02,fscal,fjz2);
1287
1288             /**************************
1289              * CALCULATE INTERACTIONS *
1290              **************************/
1291
1292             /* REACTION-FIELD ELECTROSTATICS */
1293             felec            = _mm_mul_ps(qq10,_mm_msub_ps(rinv10,rinvsq10,krf2));
1294
1295             fscal            = felec;
1296
1297              /* Update vectorial force */
1298             fix1             = _mm_macc_ps(dx10,fscal,fix1);
1299             fiy1             = _mm_macc_ps(dy10,fscal,fiy1);
1300             fiz1             = _mm_macc_ps(dz10,fscal,fiz1);
1301
1302             fjx0             = _mm_macc_ps(dx10,fscal,fjx0);
1303             fjy0             = _mm_macc_ps(dy10,fscal,fjy0);
1304             fjz0             = _mm_macc_ps(dz10,fscal,fjz0);
1305
1306             /**************************
1307              * CALCULATE INTERACTIONS *
1308              **************************/
1309
1310             /* REACTION-FIELD ELECTROSTATICS */
1311             felec            = _mm_mul_ps(qq11,_mm_msub_ps(rinv11,rinvsq11,krf2));
1312
1313             fscal            = felec;
1314
1315              /* Update vectorial force */
1316             fix1             = _mm_macc_ps(dx11,fscal,fix1);
1317             fiy1             = _mm_macc_ps(dy11,fscal,fiy1);
1318             fiz1             = _mm_macc_ps(dz11,fscal,fiz1);
1319
1320             fjx1             = _mm_macc_ps(dx11,fscal,fjx1);
1321             fjy1             = _mm_macc_ps(dy11,fscal,fjy1);
1322             fjz1             = _mm_macc_ps(dz11,fscal,fjz1);
1323
1324             /**************************
1325              * CALCULATE INTERACTIONS *
1326              **************************/
1327
1328             /* REACTION-FIELD ELECTROSTATICS */
1329             felec            = _mm_mul_ps(qq12,_mm_msub_ps(rinv12,rinvsq12,krf2));
1330
1331             fscal            = felec;
1332
1333              /* Update vectorial force */
1334             fix1             = _mm_macc_ps(dx12,fscal,fix1);
1335             fiy1             = _mm_macc_ps(dy12,fscal,fiy1);
1336             fiz1             = _mm_macc_ps(dz12,fscal,fiz1);
1337
1338             fjx2             = _mm_macc_ps(dx12,fscal,fjx2);
1339             fjy2             = _mm_macc_ps(dy12,fscal,fjy2);
1340             fjz2             = _mm_macc_ps(dz12,fscal,fjz2);
1341
1342             /**************************
1343              * CALCULATE INTERACTIONS *
1344              **************************/
1345
1346             /* REACTION-FIELD ELECTROSTATICS */
1347             felec            = _mm_mul_ps(qq20,_mm_msub_ps(rinv20,rinvsq20,krf2));
1348
1349             fscal            = felec;
1350
1351              /* Update vectorial force */
1352             fix2             = _mm_macc_ps(dx20,fscal,fix2);
1353             fiy2             = _mm_macc_ps(dy20,fscal,fiy2);
1354             fiz2             = _mm_macc_ps(dz20,fscal,fiz2);
1355
1356             fjx0             = _mm_macc_ps(dx20,fscal,fjx0);
1357             fjy0             = _mm_macc_ps(dy20,fscal,fjy0);
1358             fjz0             = _mm_macc_ps(dz20,fscal,fjz0);
1359
1360             /**************************
1361              * CALCULATE INTERACTIONS *
1362              **************************/
1363
1364             /* REACTION-FIELD ELECTROSTATICS */
1365             felec            = _mm_mul_ps(qq21,_mm_msub_ps(rinv21,rinvsq21,krf2));
1366
1367             fscal            = felec;
1368
1369              /* Update vectorial force */
1370             fix2             = _mm_macc_ps(dx21,fscal,fix2);
1371             fiy2             = _mm_macc_ps(dy21,fscal,fiy2);
1372             fiz2             = _mm_macc_ps(dz21,fscal,fiz2);
1373
1374             fjx1             = _mm_macc_ps(dx21,fscal,fjx1);
1375             fjy1             = _mm_macc_ps(dy21,fscal,fjy1);
1376             fjz1             = _mm_macc_ps(dz21,fscal,fjz1);
1377
1378             /**************************
1379              * CALCULATE INTERACTIONS *
1380              **************************/
1381
1382             /* REACTION-FIELD ELECTROSTATICS */
1383             felec            = _mm_mul_ps(qq22,_mm_msub_ps(rinv22,rinvsq22,krf2));
1384
1385             fscal            = felec;
1386
1387              /* Update vectorial force */
1388             fix2             = _mm_macc_ps(dx22,fscal,fix2);
1389             fiy2             = _mm_macc_ps(dy22,fscal,fiy2);
1390             fiz2             = _mm_macc_ps(dz22,fscal,fiz2);
1391
1392             fjx2             = _mm_macc_ps(dx22,fscal,fjx2);
1393             fjy2             = _mm_macc_ps(dy22,fscal,fjy2);
1394             fjz2             = _mm_macc_ps(dz22,fscal,fjz2);
1395
1396             fjptrA             = f+j_coord_offsetA;
1397             fjptrB             = f+j_coord_offsetB;
1398             fjptrC             = f+j_coord_offsetC;
1399             fjptrD             = f+j_coord_offsetD;
1400
1401             gmx_mm_decrement_3rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,
1402                                                    fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
1403
1404             /* Inner loop uses 297 flops */
1405         }
1406
1407         if(jidx<j_index_end)
1408         {
1409
1410             /* Get j neighbor index, and coordinate index */
1411             jnrlistA         = jjnr[jidx];
1412             jnrlistB         = jjnr[jidx+1];
1413             jnrlistC         = jjnr[jidx+2];
1414             jnrlistD         = jjnr[jidx+3];
1415             /* Sign of each element will be negative for non-real atoms.
1416              * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
1417              * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
1418              */
1419             dummy_mask = gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
1420             jnrA       = (jnrlistA>=0) ? jnrlistA : 0;
1421             jnrB       = (jnrlistB>=0) ? jnrlistB : 0;
1422             jnrC       = (jnrlistC>=0) ? jnrlistC : 0;
1423             jnrD       = (jnrlistD>=0) ? jnrlistD : 0;
1424             j_coord_offsetA  = DIM*jnrA;
1425             j_coord_offsetB  = DIM*jnrB;
1426             j_coord_offsetC  = DIM*jnrC;
1427             j_coord_offsetD  = DIM*jnrD;
1428
1429             /* load j atom coordinates */
1430             gmx_mm_load_3rvec_4ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
1431                                               x+j_coord_offsetC,x+j_coord_offsetD,
1432                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
1433
1434             /* Calculate displacement vector */
1435             dx00             = _mm_sub_ps(ix0,jx0);
1436             dy00             = _mm_sub_ps(iy0,jy0);
1437             dz00             = _mm_sub_ps(iz0,jz0);
1438             dx01             = _mm_sub_ps(ix0,jx1);
1439             dy01             = _mm_sub_ps(iy0,jy1);
1440             dz01             = _mm_sub_ps(iz0,jz1);
1441             dx02             = _mm_sub_ps(ix0,jx2);
1442             dy02             = _mm_sub_ps(iy0,jy2);
1443             dz02             = _mm_sub_ps(iz0,jz2);
1444             dx10             = _mm_sub_ps(ix1,jx0);
1445             dy10             = _mm_sub_ps(iy1,jy0);
1446             dz10             = _mm_sub_ps(iz1,jz0);
1447             dx11             = _mm_sub_ps(ix1,jx1);
1448             dy11             = _mm_sub_ps(iy1,jy1);
1449             dz11             = _mm_sub_ps(iz1,jz1);
1450             dx12             = _mm_sub_ps(ix1,jx2);
1451             dy12             = _mm_sub_ps(iy1,jy2);
1452             dz12             = _mm_sub_ps(iz1,jz2);
1453             dx20             = _mm_sub_ps(ix2,jx0);
1454             dy20             = _mm_sub_ps(iy2,jy0);
1455             dz20             = _mm_sub_ps(iz2,jz0);
1456             dx21             = _mm_sub_ps(ix2,jx1);
1457             dy21             = _mm_sub_ps(iy2,jy1);
1458             dz21             = _mm_sub_ps(iz2,jz1);
1459             dx22             = _mm_sub_ps(ix2,jx2);
1460             dy22             = _mm_sub_ps(iy2,jy2);
1461             dz22             = _mm_sub_ps(iz2,jz2);
1462
1463             /* Calculate squared distance and things based on it */
1464             rsq00            = gmx_mm_calc_rsq_ps(dx00,dy00,dz00);
1465             rsq01            = gmx_mm_calc_rsq_ps(dx01,dy01,dz01);
1466             rsq02            = gmx_mm_calc_rsq_ps(dx02,dy02,dz02);
1467             rsq10            = gmx_mm_calc_rsq_ps(dx10,dy10,dz10);
1468             rsq11            = gmx_mm_calc_rsq_ps(dx11,dy11,dz11);
1469             rsq12            = gmx_mm_calc_rsq_ps(dx12,dy12,dz12);
1470             rsq20            = gmx_mm_calc_rsq_ps(dx20,dy20,dz20);
1471             rsq21            = gmx_mm_calc_rsq_ps(dx21,dy21,dz21);
1472             rsq22            = gmx_mm_calc_rsq_ps(dx22,dy22,dz22);
1473
1474             rinv00           = gmx_mm_invsqrt_ps(rsq00);
1475             rinv01           = gmx_mm_invsqrt_ps(rsq01);
1476             rinv02           = gmx_mm_invsqrt_ps(rsq02);
1477             rinv10           = gmx_mm_invsqrt_ps(rsq10);
1478             rinv11           = gmx_mm_invsqrt_ps(rsq11);
1479             rinv12           = gmx_mm_invsqrt_ps(rsq12);
1480             rinv20           = gmx_mm_invsqrt_ps(rsq20);
1481             rinv21           = gmx_mm_invsqrt_ps(rsq21);
1482             rinv22           = gmx_mm_invsqrt_ps(rsq22);
1483
1484             rinvsq00         = _mm_mul_ps(rinv00,rinv00);
1485             rinvsq01         = _mm_mul_ps(rinv01,rinv01);
1486             rinvsq02         = _mm_mul_ps(rinv02,rinv02);
1487             rinvsq10         = _mm_mul_ps(rinv10,rinv10);
1488             rinvsq11         = _mm_mul_ps(rinv11,rinv11);
1489             rinvsq12         = _mm_mul_ps(rinv12,rinv12);
1490             rinvsq20         = _mm_mul_ps(rinv20,rinv20);
1491             rinvsq21         = _mm_mul_ps(rinv21,rinv21);
1492             rinvsq22         = _mm_mul_ps(rinv22,rinv22);
1493
1494             fjx0             = _mm_setzero_ps();
1495             fjy0             = _mm_setzero_ps();
1496             fjz0             = _mm_setzero_ps();
1497             fjx1             = _mm_setzero_ps();
1498             fjy1             = _mm_setzero_ps();
1499             fjz1             = _mm_setzero_ps();
1500             fjx2             = _mm_setzero_ps();
1501             fjy2             = _mm_setzero_ps();
1502             fjz2             = _mm_setzero_ps();
1503
1504             /**************************
1505              * CALCULATE INTERACTIONS *
1506              **************************/
1507
1508             r00              = _mm_mul_ps(rsq00,rinv00);
1509             r00              = _mm_andnot_ps(dummy_mask,r00);
1510
1511             /* Calculate table index by multiplying r with table scale and truncate to integer */
1512             rt               = _mm_mul_ps(r00,vftabscale);
1513             vfitab           = _mm_cvttps_epi32(rt);
1514 #ifdef __XOP__
1515             vfeps            = _mm_frcz_ps(rt);
1516 #else
1517             vfeps            = _mm_sub_ps(rt,_mm_round_ps(rt, _MM_FROUND_FLOOR));
1518 #endif
1519             twovfeps         = _mm_add_ps(vfeps,vfeps);
1520             vfitab           = _mm_slli_epi32(vfitab,3);
1521
1522             /* REACTION-FIELD ELECTROSTATICS */
1523             felec            = _mm_mul_ps(qq00,_mm_msub_ps(rinv00,rinvsq00,krf2));
1524
1525             /* CUBIC SPLINE TABLE DISPERSION */
1526             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
1527             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
1528             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
1529             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
1530             _MM_TRANSPOSE4_PS(Y,F,G,H);
1531             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
1532             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
1533             fvdw6            = _mm_mul_ps(c6_00,FF);
1534
1535             /* CUBIC SPLINE TABLE REPULSION */
1536             vfitab           = _mm_add_epi32(vfitab,ifour);
1537             Y                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,0) );
1538             F                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,1) );
1539             G                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,2) );
1540             H                = _mm_load_ps( vftab + _mm_extract_epi32(vfitab,3) );
1541             _MM_TRANSPOSE4_PS(Y,F,G,H);
1542             Fp               = _mm_macc_ps(vfeps,_mm_macc_ps(H,vfeps,G),F);
1543             FF               = _mm_macc_ps(vfeps,_mm_macc_ps(twovfeps,H,G),Fp);
1544             fvdw12           = _mm_mul_ps(c12_00,FF);
1545             fvdw             = _mm_xor_ps(signbit,_mm_mul_ps(_mm_add_ps(fvdw6,fvdw12),_mm_mul_ps(vftabscale,rinv00)));
1546
1547             fscal            = _mm_add_ps(felec,fvdw);
1548
1549             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1550
1551              /* Update vectorial force */
1552             fix0             = _mm_macc_ps(dx00,fscal,fix0);
1553             fiy0             = _mm_macc_ps(dy00,fscal,fiy0);
1554             fiz0             = _mm_macc_ps(dz00,fscal,fiz0);
1555
1556             fjx0             = _mm_macc_ps(dx00,fscal,fjx0);
1557             fjy0             = _mm_macc_ps(dy00,fscal,fjy0);
1558             fjz0             = _mm_macc_ps(dz00,fscal,fjz0);
1559
1560             /**************************
1561              * CALCULATE INTERACTIONS *
1562              **************************/
1563
1564             /* REACTION-FIELD ELECTROSTATICS */
1565             felec            = _mm_mul_ps(qq01,_mm_msub_ps(rinv01,rinvsq01,krf2));
1566
1567             fscal            = felec;
1568
1569             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1570
1571              /* Update vectorial force */
1572             fix0             = _mm_macc_ps(dx01,fscal,fix0);
1573             fiy0             = _mm_macc_ps(dy01,fscal,fiy0);
1574             fiz0             = _mm_macc_ps(dz01,fscal,fiz0);
1575
1576             fjx1             = _mm_macc_ps(dx01,fscal,fjx1);
1577             fjy1             = _mm_macc_ps(dy01,fscal,fjy1);
1578             fjz1             = _mm_macc_ps(dz01,fscal,fjz1);
1579
1580             /**************************
1581              * CALCULATE INTERACTIONS *
1582              **************************/
1583
1584             /* REACTION-FIELD ELECTROSTATICS */
1585             felec            = _mm_mul_ps(qq02,_mm_msub_ps(rinv02,rinvsq02,krf2));
1586
1587             fscal            = felec;
1588
1589             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1590
1591              /* Update vectorial force */
1592             fix0             = _mm_macc_ps(dx02,fscal,fix0);
1593             fiy0             = _mm_macc_ps(dy02,fscal,fiy0);
1594             fiz0             = _mm_macc_ps(dz02,fscal,fiz0);
1595
1596             fjx2             = _mm_macc_ps(dx02,fscal,fjx2);
1597             fjy2             = _mm_macc_ps(dy02,fscal,fjy2);
1598             fjz2             = _mm_macc_ps(dz02,fscal,fjz2);
1599
1600             /**************************
1601              * CALCULATE INTERACTIONS *
1602              **************************/
1603
1604             /* REACTION-FIELD ELECTROSTATICS */
1605             felec            = _mm_mul_ps(qq10,_mm_msub_ps(rinv10,rinvsq10,krf2));
1606
1607             fscal            = felec;
1608
1609             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1610
1611              /* Update vectorial force */
1612             fix1             = _mm_macc_ps(dx10,fscal,fix1);
1613             fiy1             = _mm_macc_ps(dy10,fscal,fiy1);
1614             fiz1             = _mm_macc_ps(dz10,fscal,fiz1);
1615
1616             fjx0             = _mm_macc_ps(dx10,fscal,fjx0);
1617             fjy0             = _mm_macc_ps(dy10,fscal,fjy0);
1618             fjz0             = _mm_macc_ps(dz10,fscal,fjz0);
1619
1620             /**************************
1621              * CALCULATE INTERACTIONS *
1622              **************************/
1623
1624             /* REACTION-FIELD ELECTROSTATICS */
1625             felec            = _mm_mul_ps(qq11,_mm_msub_ps(rinv11,rinvsq11,krf2));
1626
1627             fscal            = felec;
1628
1629             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1630
1631              /* Update vectorial force */
1632             fix1             = _mm_macc_ps(dx11,fscal,fix1);
1633             fiy1             = _mm_macc_ps(dy11,fscal,fiy1);
1634             fiz1             = _mm_macc_ps(dz11,fscal,fiz1);
1635
1636             fjx1             = _mm_macc_ps(dx11,fscal,fjx1);
1637             fjy1             = _mm_macc_ps(dy11,fscal,fjy1);
1638             fjz1             = _mm_macc_ps(dz11,fscal,fjz1);
1639
1640             /**************************
1641              * CALCULATE INTERACTIONS *
1642              **************************/
1643
1644             /* REACTION-FIELD ELECTROSTATICS */
1645             felec            = _mm_mul_ps(qq12,_mm_msub_ps(rinv12,rinvsq12,krf2));
1646
1647             fscal            = felec;
1648
1649             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1650
1651              /* Update vectorial force */
1652             fix1             = _mm_macc_ps(dx12,fscal,fix1);
1653             fiy1             = _mm_macc_ps(dy12,fscal,fiy1);
1654             fiz1             = _mm_macc_ps(dz12,fscal,fiz1);
1655
1656             fjx2             = _mm_macc_ps(dx12,fscal,fjx2);
1657             fjy2             = _mm_macc_ps(dy12,fscal,fjy2);
1658             fjz2             = _mm_macc_ps(dz12,fscal,fjz2);
1659
1660             /**************************
1661              * CALCULATE INTERACTIONS *
1662              **************************/
1663
1664             /* REACTION-FIELD ELECTROSTATICS */
1665             felec            = _mm_mul_ps(qq20,_mm_msub_ps(rinv20,rinvsq20,krf2));
1666
1667             fscal            = felec;
1668
1669             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1670
1671              /* Update vectorial force */
1672             fix2             = _mm_macc_ps(dx20,fscal,fix2);
1673             fiy2             = _mm_macc_ps(dy20,fscal,fiy2);
1674             fiz2             = _mm_macc_ps(dz20,fscal,fiz2);
1675
1676             fjx0             = _mm_macc_ps(dx20,fscal,fjx0);
1677             fjy0             = _mm_macc_ps(dy20,fscal,fjy0);
1678             fjz0             = _mm_macc_ps(dz20,fscal,fjz0);
1679
1680             /**************************
1681              * CALCULATE INTERACTIONS *
1682              **************************/
1683
1684             /* REACTION-FIELD ELECTROSTATICS */
1685             felec            = _mm_mul_ps(qq21,_mm_msub_ps(rinv21,rinvsq21,krf2));
1686
1687             fscal            = felec;
1688
1689             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1690
1691              /* Update vectorial force */
1692             fix2             = _mm_macc_ps(dx21,fscal,fix2);
1693             fiy2             = _mm_macc_ps(dy21,fscal,fiy2);
1694             fiz2             = _mm_macc_ps(dz21,fscal,fiz2);
1695
1696             fjx1             = _mm_macc_ps(dx21,fscal,fjx1);
1697             fjy1             = _mm_macc_ps(dy21,fscal,fjy1);
1698             fjz1             = _mm_macc_ps(dz21,fscal,fjz1);
1699
1700             /**************************
1701              * CALCULATE INTERACTIONS *
1702              **************************/
1703
1704             /* REACTION-FIELD ELECTROSTATICS */
1705             felec            = _mm_mul_ps(qq22,_mm_msub_ps(rinv22,rinvsq22,krf2));
1706
1707             fscal            = felec;
1708
1709             fscal            = _mm_andnot_ps(dummy_mask,fscal);
1710
1711              /* Update vectorial force */
1712             fix2             = _mm_macc_ps(dx22,fscal,fix2);
1713             fiy2             = _mm_macc_ps(dy22,fscal,fiy2);
1714             fiz2             = _mm_macc_ps(dz22,fscal,fiz2);
1715
1716             fjx2             = _mm_macc_ps(dx22,fscal,fjx2);
1717             fjy2             = _mm_macc_ps(dy22,fscal,fjy2);
1718             fjz2             = _mm_macc_ps(dz22,fscal,fjz2);
1719
1720             fjptrA             = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
1721             fjptrB             = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
1722             fjptrC             = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
1723             fjptrD             = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
1724
1725             gmx_mm_decrement_3rvec_4ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,
1726                                                    fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
1727
1728             /* Inner loop uses 298 flops */
1729         }
1730
1731         /* End of innermost loop */
1732
1733         gmx_mm_update_iforce_3atom_swizzle_ps(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,
1734                                               f+i_coord_offset,fshift+i_shift_offset);
1735
1736         /* Increment number of inner iterations */
1737         inneriter                  += j_index_end - j_index_start;
1738
1739         /* Outer loop uses 18 flops */
1740     }
1741
1742     /* Increment number of outer iterations */
1743     outeriter        += nri;
1744
1745     /* Update outer/inner flops */
1746
1747     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W3W3_F,outeriter*18 + inneriter*298);
1748 }