src/gromacs/gmxlib/nonbonded/nb_kernel_avx_256_single/nb_kernel_ElecRF_VdwCSTab_GeomW4W4_avx_256_single.c

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