src/gmxlib/nonbonded/nb_kernel_avx_128_fma_single/nb_kernel_ElecRFCut_VdwCSTab_GeomW3W3_avx_128_fma_single.c

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