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

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2012,2013,2014, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 /*
  36  * Note: this file was generated by the GROMACS avx_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 "gromacs/math/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_ElecCSTab_VdwLJ_GeomW4W4_VF_avx_256_single
  54  * Electrostatics interaction: CubicSplineTable
  55  * VdW interaction:            LennardJones
  56  * Geometry:                   Water4-Water4
  57  * Calculate force/pot:        PotentialAndForce
  58  */
  59 void
  60 nb_kernel_ElecCSTab_VdwLJ_GeomW4W4_VF_avx_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     nvdwtype         = fr->ntype;
 145     vdwparam         = fr->nbfp;
 146     vdwtype          = mdatoms->typeA;
 147
 148     vftab            = kernel_data->table_elec->data;
 149     vftabscale       = _mm256_set1_ps(kernel_data->table_elec->scale);
 150
 151     /* Setup water-specific parameters */
 152     inr              = nlist->iinr[0];
 153     iq1              = _mm256_mul_ps(facel,_mm256_set1_ps(charge[inr+1]));
 154     iq2              = _mm256_mul_ps(facel,_mm256_set1_ps(charge[inr+2]));
 155     iq3              = _mm256_mul_ps(facel,_mm256_set1_ps(charge[inr+3]));
 156     vdwioffsetptr0   = vdwparam+2*nvdwtype*vdwtype[inr+0];
 157
 158     jq1              = _mm256_set1_ps(charge[inr+1]);
 159     jq2              = _mm256_set1_ps(charge[inr+2]);
 160     jq3              = _mm256_set1_ps(charge[inr+3]);
 161     vdwjidx0A        = 2*vdwtype[inr+0];
 162     c6_00            = _mm256_set1_ps(vdwioffsetptr0[vdwjidx0A]);
 163     c12_00           = _mm256_set1_ps(vdwioffsetptr0[vdwjidx0A+1]);
 164     qq11             = _mm256_mul_ps(iq1,jq1);
 165     qq12             = _mm256_mul_ps(iq1,jq2);
 166     qq13             = _mm256_mul_ps(iq1,jq3);
 167     qq21             = _mm256_mul_ps(iq2,jq1);
 168     qq22             = _mm256_mul_ps(iq2,jq2);
 169     qq23             = _mm256_mul_ps(iq2,jq3);
 170     qq31             = _mm256_mul_ps(iq3,jq1);
 171     qq32             = _mm256_mul_ps(iq3,jq2);
 172     qq33             = _mm256_mul_ps(iq3,jq3);
 173
 174     /* Avoid stupid compiler warnings */
 175     jnrA = jnrB = jnrC = jnrD = jnrE = jnrF = jnrG = jnrH = 0;
 176     j_coord_offsetA = 0;
 177     j_coord_offsetB = 0;
 178     j_coord_offsetC = 0;
 179     j_coord_offsetD = 0;
 180     j_coord_offsetE = 0;
 181     j_coord_offsetF = 0;
 182     j_coord_offsetG = 0;
 183     j_coord_offsetH = 0;
 184
 185     outeriter        = 0;
 186     inneriter        = 0;
 187
 188     for(iidx=0;iidx<4*DIM;iidx++)
 189     {
 190         scratch[iidx] = 0.0;
 191     }
 192
 193     /* Start outer loop over neighborlists */
 194     for(iidx=0; iidx<nri; iidx++)
 195     {
 196         /* Load shift vector for this list */
 197         i_shift_offset   = DIM*shiftidx[iidx];
 198
 199         /* Load limits for loop over neighbors */
 200         j_index_start    = jindex[iidx];
 201         j_index_end      = jindex[iidx+1];
 202
 203         /* Get outer coordinate index */
 204         inr              = iinr[iidx];
 205         i_coord_offset   = DIM*inr;
 206
 207         /* Load i particle coords and add shift vector */
 208         gmx_mm256_load_shift_and_4rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,
 209                                                     &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
 210
 211         fix0             = _mm256_setzero_ps();
 212         fiy0             = _mm256_setzero_ps();
 213         fiz0             = _mm256_setzero_ps();
 214         fix1             = _mm256_setzero_ps();
 215         fiy1             = _mm256_setzero_ps();
 216         fiz1             = _mm256_setzero_ps();
 217         fix2             = _mm256_setzero_ps();
 218         fiy2             = _mm256_setzero_ps();
 219         fiz2             = _mm256_setzero_ps();
 220         fix3             = _mm256_setzero_ps();
 221         fiy3             = _mm256_setzero_ps();
 222         fiz3             = _mm256_setzero_ps();
 223
 224         /* Reset potential sums */
 225         velecsum         = _mm256_setzero_ps();
 226         vvdwsum          = _mm256_setzero_ps();
 227
 228         /* Start inner kernel loop */
 229         for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+7]>=0; jidx+=8)
 230         {
 231
 232             /* Get j neighbor index, and coordinate index */
 233             jnrA             = jjnr[jidx];
 234             jnrB             = jjnr[jidx+1];
 235             jnrC             = jjnr[jidx+2];
 236             jnrD             = jjnr[jidx+3];
 237             jnrE             = jjnr[jidx+4];
 238             jnrF             = jjnr[jidx+5];
 239             jnrG             = jjnr[jidx+6];
 240             jnrH             = jjnr[jidx+7];
 241             j_coord_offsetA  = DIM*jnrA;
 242             j_coord_offsetB  = DIM*jnrB;
 243             j_coord_offsetC  = DIM*jnrC;
 244             j_coord_offsetD  = DIM*jnrD;
 245             j_coord_offsetE  = DIM*jnrE;
 246             j_coord_offsetF  = DIM*jnrF;
 247             j_coord_offsetG  = DIM*jnrG;
 248             j_coord_offsetH  = DIM*jnrH;
 249
 250             /* load j atom coordinates */
 251             gmx_mm256_load_4rvec_8ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
 252                                                  x+j_coord_offsetC,x+j_coord_offsetD,
 253                                                  x+j_coord_offsetE,x+j_coord_offsetF,
 254                                                  x+j_coord_offsetG,x+j_coord_offsetH,
 255                                                  &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
 256                                                  &jy2,&jz2,&jx3,&jy3,&jz3);
 257
 258             /* Calculate displacement vector */
 259             dx00             = _mm256_sub_ps(ix0,jx0);
 260             dy00             = _mm256_sub_ps(iy0,jy0);
 261             dz00             = _mm256_sub_ps(iz0,jz0);
 262             dx11             = _mm256_sub_ps(ix1,jx1);
 263             dy11             = _mm256_sub_ps(iy1,jy1);
 264             dz11             = _mm256_sub_ps(iz1,jz1);
 265             dx12             = _mm256_sub_ps(ix1,jx2);
 266             dy12             = _mm256_sub_ps(iy1,jy2);
 267             dz12             = _mm256_sub_ps(iz1,jz2);
 268             dx13             = _mm256_sub_ps(ix1,jx3);
 269             dy13             = _mm256_sub_ps(iy1,jy3);
 270             dz13             = _mm256_sub_ps(iz1,jz3);
 271             dx21             = _mm256_sub_ps(ix2,jx1);
 272             dy21             = _mm256_sub_ps(iy2,jy1);
 273             dz21             = _mm256_sub_ps(iz2,jz1);
 274             dx22             = _mm256_sub_ps(ix2,jx2);
 275             dy22             = _mm256_sub_ps(iy2,jy2);
 276             dz22             = _mm256_sub_ps(iz2,jz2);
 277             dx23             = _mm256_sub_ps(ix2,jx3);
 278             dy23             = _mm256_sub_ps(iy2,jy3);
 279             dz23             = _mm256_sub_ps(iz2,jz3);
 280             dx31             = _mm256_sub_ps(ix3,jx1);
 281             dy31             = _mm256_sub_ps(iy3,jy1);
 282             dz31             = _mm256_sub_ps(iz3,jz1);
 283             dx32             = _mm256_sub_ps(ix3,jx2);
 284             dy32             = _mm256_sub_ps(iy3,jy2);
 285             dz32             = _mm256_sub_ps(iz3,jz2);
 286             dx33             = _mm256_sub_ps(ix3,jx3);
 287             dy33             = _mm256_sub_ps(iy3,jy3);
 288             dz33             = _mm256_sub_ps(iz3,jz3);
 289
 290             /* Calculate squared distance and things based on it */
 291             rsq00            = gmx_mm256_calc_rsq_ps(dx00,dy00,dz00);
 292             rsq11            = gmx_mm256_calc_rsq_ps(dx11,dy11,dz11);
 293             rsq12            = gmx_mm256_calc_rsq_ps(dx12,dy12,dz12);
 294             rsq13            = gmx_mm256_calc_rsq_ps(dx13,dy13,dz13);
 295             rsq21            = gmx_mm256_calc_rsq_ps(dx21,dy21,dz21);
 296             rsq22            = gmx_mm256_calc_rsq_ps(dx22,dy22,dz22);
 297             rsq23            = gmx_mm256_calc_rsq_ps(dx23,dy23,dz23);
 298             rsq31            = gmx_mm256_calc_rsq_ps(dx31,dy31,dz31);
 299             rsq32            = gmx_mm256_calc_rsq_ps(dx32,dy32,dz32);
 300             rsq33            = gmx_mm256_calc_rsq_ps(dx33,dy33,dz33);
 301
 302             rinv11           = gmx_mm256_invsqrt_ps(rsq11);
 303             rinv12           = gmx_mm256_invsqrt_ps(rsq12);
 304             rinv13           = gmx_mm256_invsqrt_ps(rsq13);
 305             rinv21           = gmx_mm256_invsqrt_ps(rsq21);
 306             rinv22           = gmx_mm256_invsqrt_ps(rsq22);
 307             rinv23           = gmx_mm256_invsqrt_ps(rsq23);
 308             rinv31           = gmx_mm256_invsqrt_ps(rsq31);
 309             rinv32           = gmx_mm256_invsqrt_ps(rsq32);
 310             rinv33           = gmx_mm256_invsqrt_ps(rsq33);
 311
 312             rinvsq00         = gmx_mm256_inv_ps(rsq00);
 313
 314             fjx0             = _mm256_setzero_ps();
 315             fjy0             = _mm256_setzero_ps();
 316             fjz0             = _mm256_setzero_ps();
 317             fjx1             = _mm256_setzero_ps();
 318             fjy1             = _mm256_setzero_ps();
 319             fjz1             = _mm256_setzero_ps();
 320             fjx2             = _mm256_setzero_ps();
 321             fjy2             = _mm256_setzero_ps();
 322             fjz2             = _mm256_setzero_ps();
 323             fjx3             = _mm256_setzero_ps();
 324             fjy3             = _mm256_setzero_ps();
 325             fjz3             = _mm256_setzero_ps();
 326
 327             /**************************
 328              * CALCULATE INTERACTIONS *
 329              **************************/
 330
 331             /* LENNARD-JONES DISPERSION/REPULSION */
 332
 333             rinvsix          = _mm256_mul_ps(_mm256_mul_ps(rinvsq00,rinvsq00),rinvsq00);
 334             vvdw6            = _mm256_mul_ps(c6_00,rinvsix);
 335             vvdw12           = _mm256_mul_ps(c12_00,_mm256_mul_ps(rinvsix,rinvsix));
 336             vvdw             = _mm256_sub_ps( _mm256_mul_ps(vvdw12,one_twelfth) , _mm256_mul_ps(vvdw6,one_sixth) );
 337             fvdw             = _mm256_mul_ps(_mm256_sub_ps(vvdw12,vvdw6),rinvsq00);
 338
 339             /* Update potential sum for this i atom from the interaction with this j atom. */
 340             vvdwsum          = _mm256_add_ps(vvdwsum,vvdw);
 341
 342             fscal            = fvdw;
 343
 344             /* Calculate temporary vectorial force */
 345             tx               = _mm256_mul_ps(fscal,dx00);
 346             ty               = _mm256_mul_ps(fscal,dy00);
 347             tz               = _mm256_mul_ps(fscal,dz00);
 348
 349             /* Update vectorial force */
 350             fix0             = _mm256_add_ps(fix0,tx);
 351             fiy0             = _mm256_add_ps(fiy0,ty);
 352             fiz0             = _mm256_add_ps(fiz0,tz);
 353
 354             fjx0             = _mm256_add_ps(fjx0,tx);
 355             fjy0             = _mm256_add_ps(fjy0,ty);
 356             fjz0             = _mm256_add_ps(fjz0,tz);
 357
 358             /**************************
 359              * CALCULATE INTERACTIONS *
 360              **************************/
 361
 362             r11              = _mm256_mul_ps(rsq11,rinv11);
 363
 364             /* Calculate table index by multiplying r with table scale and truncate to integer */
 365             rt               = _mm256_mul_ps(r11,vftabscale);
 366             vfitab           = _mm256_cvttps_epi32(rt);
 367             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 368             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 369             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 370             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 371             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 372             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 373
 374             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 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             velec            = _mm256_mul_ps(qq11,VV);
 388             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 389             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq11,FF),_mm256_mul_ps(vftabscale,rinv11)));
 390
 391             /* Update potential sum for this i atom from the interaction with this j atom. */
 392             velecsum         = _mm256_add_ps(velecsum,velec);
 393
 394             fscal            = felec;
 395
 396             /* Calculate temporary vectorial force */
 397             tx               = _mm256_mul_ps(fscal,dx11);
 398             ty               = _mm256_mul_ps(fscal,dy11);
 399             tz               = _mm256_mul_ps(fscal,dz11);
 400
 401             /* Update vectorial force */
 402             fix1             = _mm256_add_ps(fix1,tx);
 403             fiy1             = _mm256_add_ps(fiy1,ty);
 404             fiz1             = _mm256_add_ps(fiz1,tz);
 405
 406             fjx1             = _mm256_add_ps(fjx1,tx);
 407             fjy1             = _mm256_add_ps(fjy1,ty);
 408             fjz1             = _mm256_add_ps(fjz1,tz);
 409
 410             /**************************
 411              * CALCULATE INTERACTIONS *
 412              **************************/
 413
 414             r12              = _mm256_mul_ps(rsq12,rinv12);
 415
 416             /* Calculate table index by multiplying r with table scale and truncate to integer */
 417             rt               = _mm256_mul_ps(r12,vftabscale);
 418             vfitab           = _mm256_cvttps_epi32(rt);
 419             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 420             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 421             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 422             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 423             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 424             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 425
 426             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 427             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 428                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 429             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 430                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 431             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 432                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 433             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 434                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 435             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 436             Heps             = _mm256_mul_ps(vfeps,H);
 437             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 438             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 439             velec            = _mm256_mul_ps(qq12,VV);
 440             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 441             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq12,FF),_mm256_mul_ps(vftabscale,rinv12)));
 442
 443             /* Update potential sum for this i atom from the interaction with this j atom. */
 444             velecsum         = _mm256_add_ps(velecsum,velec);
 445
 446             fscal            = felec;
 447
 448             /* Calculate temporary vectorial force */
 449             tx               = _mm256_mul_ps(fscal,dx12);
 450             ty               = _mm256_mul_ps(fscal,dy12);
 451             tz               = _mm256_mul_ps(fscal,dz12);
 452
 453             /* Update vectorial force */
 454             fix1             = _mm256_add_ps(fix1,tx);
 455             fiy1             = _mm256_add_ps(fiy1,ty);
 456             fiz1             = _mm256_add_ps(fiz1,tz);
 457
 458             fjx2             = _mm256_add_ps(fjx2,tx);
 459             fjy2             = _mm256_add_ps(fjy2,ty);
 460             fjz2             = _mm256_add_ps(fjz2,tz);
 461
 462             /**************************
 463              * CALCULATE INTERACTIONS *
 464              **************************/
 465
 466             r13              = _mm256_mul_ps(rsq13,rinv13);
 467
 468             /* Calculate table index by multiplying r with table scale and truncate to integer */
 469             rt               = _mm256_mul_ps(r13,vftabscale);
 470             vfitab           = _mm256_cvttps_epi32(rt);
 471             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 472             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 473             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 474             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 475             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 476             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 477
 478             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 479             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 480                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 481             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 482                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 483             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 484                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 485             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 486                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 487             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 488             Heps             = _mm256_mul_ps(vfeps,H);
 489             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 490             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 491             velec            = _mm256_mul_ps(qq13,VV);
 492             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 493             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq13,FF),_mm256_mul_ps(vftabscale,rinv13)));
 494
 495             /* Update potential sum for this i atom from the interaction with this j atom. */
 496             velecsum         = _mm256_add_ps(velecsum,velec);
 497
 498             fscal            = felec;
 499
 500             /* Calculate temporary vectorial force */
 501             tx               = _mm256_mul_ps(fscal,dx13);
 502             ty               = _mm256_mul_ps(fscal,dy13);
 503             tz               = _mm256_mul_ps(fscal,dz13);
 504
 505             /* Update vectorial force */
 506             fix1             = _mm256_add_ps(fix1,tx);
 507             fiy1             = _mm256_add_ps(fiy1,ty);
 508             fiz1             = _mm256_add_ps(fiz1,tz);
 509
 510             fjx3             = _mm256_add_ps(fjx3,tx);
 511             fjy3             = _mm256_add_ps(fjy3,ty);
 512             fjz3             = _mm256_add_ps(fjz3,tz);
 513
 514             /**************************
 515              * CALCULATE INTERACTIONS *
 516              **************************/
 517
 518             r21              = _mm256_mul_ps(rsq21,rinv21);
 519
 520             /* Calculate table index by multiplying r with table scale and truncate to integer */
 521             rt               = _mm256_mul_ps(r21,vftabscale);
 522             vfitab           = _mm256_cvttps_epi32(rt);
 523             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 524             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 525             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 526             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 527             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 528             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 529
 530             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 531             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 532                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 533             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 534                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 535             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 536                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 537             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 538                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 539             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 540             Heps             = _mm256_mul_ps(vfeps,H);
 541             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 542             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 543             velec            = _mm256_mul_ps(qq21,VV);
 544             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 545             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq21,FF),_mm256_mul_ps(vftabscale,rinv21)));
 546
 547             /* Update potential sum for this i atom from the interaction with this j atom. */
 548             velecsum         = _mm256_add_ps(velecsum,velec);
 549
 550             fscal            = felec;
 551
 552             /* Calculate temporary vectorial force */
 553             tx               = _mm256_mul_ps(fscal,dx21);
 554             ty               = _mm256_mul_ps(fscal,dy21);
 555             tz               = _mm256_mul_ps(fscal,dz21);
 556
 557             /* Update vectorial force */
 558             fix2             = _mm256_add_ps(fix2,tx);
 559             fiy2             = _mm256_add_ps(fiy2,ty);
 560             fiz2             = _mm256_add_ps(fiz2,tz);
 561
 562             fjx1             = _mm256_add_ps(fjx1,tx);
 563             fjy1             = _mm256_add_ps(fjy1,ty);
 564             fjz1             = _mm256_add_ps(fjz1,tz);
 565
 566             /**************************
 567              * CALCULATE INTERACTIONS *
 568              **************************/
 569
 570             r22              = _mm256_mul_ps(rsq22,rinv22);
 571
 572             /* Calculate table index by multiplying r with table scale and truncate to integer */
 573             rt               = _mm256_mul_ps(r22,vftabscale);
 574             vfitab           = _mm256_cvttps_epi32(rt);
 575             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 576             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 577             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 578             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 579             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 580             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 581
 582             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 583             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 584                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 585             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 586                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 587             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 588                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 589             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 590                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 591             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 592             Heps             = _mm256_mul_ps(vfeps,H);
 593             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 594             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 595             velec            = _mm256_mul_ps(qq22,VV);
 596             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 597             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq22,FF),_mm256_mul_ps(vftabscale,rinv22)));
 598
 599             /* Update potential sum for this i atom from the interaction with this j atom. */
 600             velecsum         = _mm256_add_ps(velecsum,velec);
 601
 602             fscal            = felec;
 603
 604             /* Calculate temporary vectorial force */
 605             tx               = _mm256_mul_ps(fscal,dx22);
 606             ty               = _mm256_mul_ps(fscal,dy22);
 607             tz               = _mm256_mul_ps(fscal,dz22);
 608
 609             /* Update vectorial force */
 610             fix2             = _mm256_add_ps(fix2,tx);
 611             fiy2             = _mm256_add_ps(fiy2,ty);
 612             fiz2             = _mm256_add_ps(fiz2,tz);
 613
 614             fjx2             = _mm256_add_ps(fjx2,tx);
 615             fjy2             = _mm256_add_ps(fjy2,ty);
 616             fjz2             = _mm256_add_ps(fjz2,tz);
 617
 618             /**************************
 619              * CALCULATE INTERACTIONS *
 620              **************************/
 621
 622             r23              = _mm256_mul_ps(rsq23,rinv23);
 623
 624             /* Calculate table index by multiplying r with table scale and truncate to integer */
 625             rt               = _mm256_mul_ps(r23,vftabscale);
 626             vfitab           = _mm256_cvttps_epi32(rt);
 627             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 628             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 629             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 630             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 631             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 632             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 633
 634             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 635             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 636                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 637             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 638                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 639             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 640                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 641             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 642                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 643             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 644             Heps             = _mm256_mul_ps(vfeps,H);
 645             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 646             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 647             velec            = _mm256_mul_ps(qq23,VV);
 648             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 649             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq23,FF),_mm256_mul_ps(vftabscale,rinv23)));
 650
 651             /* Update potential sum for this i atom from the interaction with this j atom. */
 652             velecsum         = _mm256_add_ps(velecsum,velec);
 653
 654             fscal            = felec;
 655
 656             /* Calculate temporary vectorial force */
 657             tx               = _mm256_mul_ps(fscal,dx23);
 658             ty               = _mm256_mul_ps(fscal,dy23);
 659             tz               = _mm256_mul_ps(fscal,dz23);
 660
 661             /* Update vectorial force */
 662             fix2             = _mm256_add_ps(fix2,tx);
 663             fiy2             = _mm256_add_ps(fiy2,ty);
 664             fiz2             = _mm256_add_ps(fiz2,tz);
 665
 666             fjx3             = _mm256_add_ps(fjx3,tx);
 667             fjy3             = _mm256_add_ps(fjy3,ty);
 668             fjz3             = _mm256_add_ps(fjz3,tz);
 669
 670             /**************************
 671              * CALCULATE INTERACTIONS *
 672              **************************/
 673
 674             r31              = _mm256_mul_ps(rsq31,rinv31);
 675
 676             /* Calculate table index by multiplying r with table scale and truncate to integer */
 677             rt               = _mm256_mul_ps(r31,vftabscale);
 678             vfitab           = _mm256_cvttps_epi32(rt);
 679             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 680             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 681             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 682             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 683             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 684             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 685
 686             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 687             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 688                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 689             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 690                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 691             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 692                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 693             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 694                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 695             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 696             Heps             = _mm256_mul_ps(vfeps,H);
 697             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 698             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 699             velec            = _mm256_mul_ps(qq31,VV);
 700             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 701             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq31,FF),_mm256_mul_ps(vftabscale,rinv31)));
 702
 703             /* Update potential sum for this i atom from the interaction with this j atom. */
 704             velecsum         = _mm256_add_ps(velecsum,velec);
 705
 706             fscal            = felec;
 707
 708             /* Calculate temporary vectorial force */
 709             tx               = _mm256_mul_ps(fscal,dx31);
 710             ty               = _mm256_mul_ps(fscal,dy31);
 711             tz               = _mm256_mul_ps(fscal,dz31);
 712
 713             /* Update vectorial force */
 714             fix3             = _mm256_add_ps(fix3,tx);
 715             fiy3             = _mm256_add_ps(fiy3,ty);
 716             fiz3             = _mm256_add_ps(fiz3,tz);
 717
 718             fjx1             = _mm256_add_ps(fjx1,tx);
 719             fjy1             = _mm256_add_ps(fjy1,ty);
 720             fjz1             = _mm256_add_ps(fjz1,tz);
 721
 722             /**************************
 723              * CALCULATE INTERACTIONS *
 724              **************************/
 725
 726             r32              = _mm256_mul_ps(rsq32,rinv32);
 727
 728             /* Calculate table index by multiplying r with table scale and truncate to integer */
 729             rt               = _mm256_mul_ps(r32,vftabscale);
 730             vfitab           = _mm256_cvttps_epi32(rt);
 731             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 732             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 733             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 734             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 735             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 736             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 737
 738             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 739             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 740                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 741             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 742                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 743             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 744                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 745             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 746                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 747             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 748             Heps             = _mm256_mul_ps(vfeps,H);
 749             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 750             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 751             velec            = _mm256_mul_ps(qq32,VV);
 752             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 753             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq32,FF),_mm256_mul_ps(vftabscale,rinv32)));
 754
 755             /* Update potential sum for this i atom from the interaction with this j atom. */
 756             velecsum         = _mm256_add_ps(velecsum,velec);
 757
 758             fscal            = felec;
 759
 760             /* Calculate temporary vectorial force */
 761             tx               = _mm256_mul_ps(fscal,dx32);
 762             ty               = _mm256_mul_ps(fscal,dy32);
 763             tz               = _mm256_mul_ps(fscal,dz32);
 764
 765             /* Update vectorial force */
 766             fix3             = _mm256_add_ps(fix3,tx);
 767             fiy3             = _mm256_add_ps(fiy3,ty);
 768             fiz3             = _mm256_add_ps(fiz3,tz);
 769
 770             fjx2             = _mm256_add_ps(fjx2,tx);
 771             fjy2             = _mm256_add_ps(fjy2,ty);
 772             fjz2             = _mm256_add_ps(fjz2,tz);
 773
 774             /**************************
 775              * CALCULATE INTERACTIONS *
 776              **************************/
 777
 778             r33              = _mm256_mul_ps(rsq33,rinv33);
 779
 780             /* Calculate table index by multiplying r with table scale and truncate to integer */
 781             rt               = _mm256_mul_ps(r33,vftabscale);
 782             vfitab           = _mm256_cvttps_epi32(rt);
 783             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
 784             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
 785             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
 786             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
 787             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
 788             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
 789
 790             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 791             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
 792                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
 793             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
 794                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
 795             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
 796                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
 797             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
 798                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
 799             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
 800             Heps             = _mm256_mul_ps(vfeps,H);
 801             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
 802             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
 803             velec            = _mm256_mul_ps(qq33,VV);
 804             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
 805             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq33,FF),_mm256_mul_ps(vftabscale,rinv33)));
 806
 807             /* Update potential sum for this i atom from the interaction with this j atom. */
 808             velecsum         = _mm256_add_ps(velecsum,velec);
 809
 810             fscal            = felec;
 811
 812             /* Calculate temporary vectorial force */
 813             tx               = _mm256_mul_ps(fscal,dx33);
 814             ty               = _mm256_mul_ps(fscal,dy33);
 815             tz               = _mm256_mul_ps(fscal,dz33);
 816
 817             /* Update vectorial force */
 818             fix3             = _mm256_add_ps(fix3,tx);
 819             fiy3             = _mm256_add_ps(fiy3,ty);
 820             fiz3             = _mm256_add_ps(fiz3,tz);
 821
 822             fjx3             = _mm256_add_ps(fjx3,tx);
 823             fjy3             = _mm256_add_ps(fjy3,ty);
 824             fjz3             = _mm256_add_ps(fjz3,tz);
 825
 826             fjptrA             = f+j_coord_offsetA;
 827             fjptrB             = f+j_coord_offsetB;
 828             fjptrC             = f+j_coord_offsetC;
 829             fjptrD             = f+j_coord_offsetD;
 830             fjptrE             = f+j_coord_offsetE;
 831             fjptrF             = f+j_coord_offsetF;
 832             fjptrG             = f+j_coord_offsetG;
 833             fjptrH             = f+j_coord_offsetH;
 834
 835             gmx_mm256_decrement_4rvec_8ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjptrE,fjptrF,fjptrG,fjptrH,
 836                                                       fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,
 837                                                       fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
 838
 839             /* Inner loop uses 422 flops */
 840         }
 841
 842         if(jidx<j_index_end)
 843         {
 844
 845             /* Get j neighbor index, and coordinate index */
 846             jnrlistA         = jjnr[jidx];
 847             jnrlistB         = jjnr[jidx+1];
 848             jnrlistC         = jjnr[jidx+2];
 849             jnrlistD         = jjnr[jidx+3];
 850             jnrlistE         = jjnr[jidx+4];
 851             jnrlistF         = jjnr[jidx+5];
 852             jnrlistG         = jjnr[jidx+6];
 853             jnrlistH         = jjnr[jidx+7];
 854             /* Sign of each element will be negative for non-real atoms.
 855              * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
 856              * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
 857              */
 858             dummy_mask = gmx_mm256_set_m128(gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx+4)),_mm_setzero_si128())),
 859                                             gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128())));
 860
 861             jnrA       = (jnrlistA>=0) ? jnrlistA : 0;
 862             jnrB       = (jnrlistB>=0) ? jnrlistB : 0;
 863             jnrC       = (jnrlistC>=0) ? jnrlistC : 0;
 864             jnrD       = (jnrlistD>=0) ? jnrlistD : 0;
 865             jnrE       = (jnrlistE>=0) ? jnrlistE : 0;
 866             jnrF       = (jnrlistF>=0) ? jnrlistF : 0;
 867             jnrG       = (jnrlistG>=0) ? jnrlistG : 0;
 868             jnrH       = (jnrlistH>=0) ? jnrlistH : 0;
 869             j_coord_offsetA  = DIM*jnrA;
 870             j_coord_offsetB  = DIM*jnrB;
 871             j_coord_offsetC  = DIM*jnrC;
 872             j_coord_offsetD  = DIM*jnrD;
 873             j_coord_offsetE  = DIM*jnrE;
 874             j_coord_offsetF  = DIM*jnrF;
 875             j_coord_offsetG  = DIM*jnrG;
 876             j_coord_offsetH  = DIM*jnrH;
 877
 878             /* load j atom coordinates */
 879             gmx_mm256_load_4rvec_8ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
 880                                                  x+j_coord_offsetC,x+j_coord_offsetD,
 881                                                  x+j_coord_offsetE,x+j_coord_offsetF,
 882                                                  x+j_coord_offsetG,x+j_coord_offsetH,
 883                                                  &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
 884                                                  &jy2,&jz2,&jx3,&jy3,&jz3);
 885
 886             /* Calculate displacement vector */
 887             dx00             = _mm256_sub_ps(ix0,jx0);
 888             dy00             = _mm256_sub_ps(iy0,jy0);
 889             dz00             = _mm256_sub_ps(iz0,jz0);
 890             dx11             = _mm256_sub_ps(ix1,jx1);
 891             dy11             = _mm256_sub_ps(iy1,jy1);
 892             dz11             = _mm256_sub_ps(iz1,jz1);
 893             dx12             = _mm256_sub_ps(ix1,jx2);
 894             dy12             = _mm256_sub_ps(iy1,jy2);
 895             dz12             = _mm256_sub_ps(iz1,jz2);
 896             dx13             = _mm256_sub_ps(ix1,jx3);
 897             dy13             = _mm256_sub_ps(iy1,jy3);
 898             dz13             = _mm256_sub_ps(iz1,jz3);
 899             dx21             = _mm256_sub_ps(ix2,jx1);
 900             dy21             = _mm256_sub_ps(iy2,jy1);
 901             dz21             = _mm256_sub_ps(iz2,jz1);
 902             dx22             = _mm256_sub_ps(ix2,jx2);
 903             dy22             = _mm256_sub_ps(iy2,jy2);
 904             dz22             = _mm256_sub_ps(iz2,jz2);
 905             dx23             = _mm256_sub_ps(ix2,jx3);
 906             dy23             = _mm256_sub_ps(iy2,jy3);
 907             dz23             = _mm256_sub_ps(iz2,jz3);
 908             dx31             = _mm256_sub_ps(ix3,jx1);
 909             dy31             = _mm256_sub_ps(iy3,jy1);
 910             dz31             = _mm256_sub_ps(iz3,jz1);
 911             dx32             = _mm256_sub_ps(ix3,jx2);
 912             dy32             = _mm256_sub_ps(iy3,jy2);
 913             dz32             = _mm256_sub_ps(iz3,jz2);
 914             dx33             = _mm256_sub_ps(ix3,jx3);
 915             dy33             = _mm256_sub_ps(iy3,jy3);
 916             dz33             = _mm256_sub_ps(iz3,jz3);
 917
 918             /* Calculate squared distance and things based on it */
 919             rsq00            = gmx_mm256_calc_rsq_ps(dx00,dy00,dz00);
 920             rsq11            = gmx_mm256_calc_rsq_ps(dx11,dy11,dz11);
 921             rsq12            = gmx_mm256_calc_rsq_ps(dx12,dy12,dz12);
 922             rsq13            = gmx_mm256_calc_rsq_ps(dx13,dy13,dz13);
 923             rsq21            = gmx_mm256_calc_rsq_ps(dx21,dy21,dz21);
 924             rsq22            = gmx_mm256_calc_rsq_ps(dx22,dy22,dz22);
 925             rsq23            = gmx_mm256_calc_rsq_ps(dx23,dy23,dz23);
 926             rsq31            = gmx_mm256_calc_rsq_ps(dx31,dy31,dz31);
 927             rsq32            = gmx_mm256_calc_rsq_ps(dx32,dy32,dz32);
 928             rsq33            = gmx_mm256_calc_rsq_ps(dx33,dy33,dz33);
 929
 930             rinv11           = gmx_mm256_invsqrt_ps(rsq11);
 931             rinv12           = gmx_mm256_invsqrt_ps(rsq12);
 932             rinv13           = gmx_mm256_invsqrt_ps(rsq13);
 933             rinv21           = gmx_mm256_invsqrt_ps(rsq21);
 934             rinv22           = gmx_mm256_invsqrt_ps(rsq22);
 935             rinv23           = gmx_mm256_invsqrt_ps(rsq23);
 936             rinv31           = gmx_mm256_invsqrt_ps(rsq31);
 937             rinv32           = gmx_mm256_invsqrt_ps(rsq32);
 938             rinv33           = gmx_mm256_invsqrt_ps(rsq33);
 939
 940             rinvsq00         = gmx_mm256_inv_ps(rsq00);
 941
 942             fjx0             = _mm256_setzero_ps();
 943             fjy0             = _mm256_setzero_ps();
 944             fjz0             = _mm256_setzero_ps();
 945             fjx1             = _mm256_setzero_ps();
 946             fjy1             = _mm256_setzero_ps();
 947             fjz1             = _mm256_setzero_ps();
 948             fjx2             = _mm256_setzero_ps();
 949             fjy2             = _mm256_setzero_ps();
 950             fjz2             = _mm256_setzero_ps();
 951             fjx3             = _mm256_setzero_ps();
 952             fjy3             = _mm256_setzero_ps();
 953             fjz3             = _mm256_setzero_ps();
 954
 955             /**************************
 956              * CALCULATE INTERACTIONS *
 957              **************************/
 958
 959             /* LENNARD-JONES DISPERSION/REPULSION */
 960
 961             rinvsix          = _mm256_mul_ps(_mm256_mul_ps(rinvsq00,rinvsq00),rinvsq00);
 962             vvdw6            = _mm256_mul_ps(c6_00,rinvsix);
 963             vvdw12           = _mm256_mul_ps(c12_00,_mm256_mul_ps(rinvsix,rinvsix));
 964             vvdw             = _mm256_sub_ps( _mm256_mul_ps(vvdw12,one_twelfth) , _mm256_mul_ps(vvdw6,one_sixth) );
 965             fvdw             = _mm256_mul_ps(_mm256_sub_ps(vvdw12,vvdw6),rinvsq00);
 966
 967             /* Update potential sum for this i atom from the interaction with this j atom. */
 968             vvdw             = _mm256_andnot_ps(dummy_mask,vvdw);
 969             vvdwsum          = _mm256_add_ps(vvdwsum,vvdw);
 970
 971             fscal            = fvdw;
 972
 973             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
 974
 975             /* Calculate temporary vectorial force */
 976             tx               = _mm256_mul_ps(fscal,dx00);
 977             ty               = _mm256_mul_ps(fscal,dy00);
 978             tz               = _mm256_mul_ps(fscal,dz00);
 979
 980             /* Update vectorial force */
 981             fix0             = _mm256_add_ps(fix0,tx);
 982             fiy0             = _mm256_add_ps(fiy0,ty);
 983             fiz0             = _mm256_add_ps(fiz0,tz);
 984
 985             fjx0             = _mm256_add_ps(fjx0,tx);
 986             fjy0             = _mm256_add_ps(fjy0,ty);
 987             fjz0             = _mm256_add_ps(fjz0,tz);
 988
 989             /**************************
 990              * CALCULATE INTERACTIONS *
 991              **************************/
 992
 993             r11              = _mm256_mul_ps(rsq11,rinv11);
 994             r11              = _mm256_andnot_ps(dummy_mask,r11);
 995
 996             /* Calculate table index by multiplying r with table scale and truncate to integer */
 997             rt               = _mm256_mul_ps(r11,vftabscale);
 998             vfitab           = _mm256_cvttps_epi32(rt);
 999             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1000             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1001             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1002             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1003             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1004             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1005
1006             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1007             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1008                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1009             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1010                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1011             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1012                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1013             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1014                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1015             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1016             Heps             = _mm256_mul_ps(vfeps,H);
1017             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1018             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1019             velec            = _mm256_mul_ps(qq11,VV);
1020             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1021             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq11,FF),_mm256_mul_ps(vftabscale,rinv11)));
1022
1023             /* Update potential sum for this i atom from the interaction with this j atom. */
1024             velec            = _mm256_andnot_ps(dummy_mask,velec);
1025             velecsum         = _mm256_add_ps(velecsum,velec);
1026
1027             fscal            = felec;
1028
1029             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1030
1031             /* Calculate temporary vectorial force */
1032             tx               = _mm256_mul_ps(fscal,dx11);
1033             ty               = _mm256_mul_ps(fscal,dy11);
1034             tz               = _mm256_mul_ps(fscal,dz11);
1035
1036             /* Update vectorial force */
1037             fix1             = _mm256_add_ps(fix1,tx);
1038             fiy1             = _mm256_add_ps(fiy1,ty);
1039             fiz1             = _mm256_add_ps(fiz1,tz);
1040
1041             fjx1             = _mm256_add_ps(fjx1,tx);
1042             fjy1             = _mm256_add_ps(fjy1,ty);
1043             fjz1             = _mm256_add_ps(fjz1,tz);
1044
1045             /**************************
1046              * CALCULATE INTERACTIONS *
1047              **************************/
1048
1049             r12              = _mm256_mul_ps(rsq12,rinv12);
1050             r12              = _mm256_andnot_ps(dummy_mask,r12);
1051
1052             /* Calculate table index by multiplying r with table scale and truncate to integer */
1053             rt               = _mm256_mul_ps(r12,vftabscale);
1054             vfitab           = _mm256_cvttps_epi32(rt);
1055             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1056             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1057             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1058             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1059             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1060             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1061
1062             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1063             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1064                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1065             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1066                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1067             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1068                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1069             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1070                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1071             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1072             Heps             = _mm256_mul_ps(vfeps,H);
1073             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1074             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1075             velec            = _mm256_mul_ps(qq12,VV);
1076             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1077             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq12,FF),_mm256_mul_ps(vftabscale,rinv12)));
1078
1079             /* Update potential sum for this i atom from the interaction with this j atom. */
1080             velec            = _mm256_andnot_ps(dummy_mask,velec);
1081             velecsum         = _mm256_add_ps(velecsum,velec);
1082
1083             fscal            = felec;
1084
1085             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1086
1087             /* Calculate temporary vectorial force */
1088             tx               = _mm256_mul_ps(fscal,dx12);
1089             ty               = _mm256_mul_ps(fscal,dy12);
1090             tz               = _mm256_mul_ps(fscal,dz12);
1091
1092             /* Update vectorial force */
1093             fix1             = _mm256_add_ps(fix1,tx);
1094             fiy1             = _mm256_add_ps(fiy1,ty);
1095             fiz1             = _mm256_add_ps(fiz1,tz);
1096
1097             fjx2             = _mm256_add_ps(fjx2,tx);
1098             fjy2             = _mm256_add_ps(fjy2,ty);
1099             fjz2             = _mm256_add_ps(fjz2,tz);
1100
1101             /**************************
1102              * CALCULATE INTERACTIONS *
1103              **************************/
1104
1105             r13              = _mm256_mul_ps(rsq13,rinv13);
1106             r13              = _mm256_andnot_ps(dummy_mask,r13);
1107
1108             /* Calculate table index by multiplying r with table scale and truncate to integer */
1109             rt               = _mm256_mul_ps(r13,vftabscale);
1110             vfitab           = _mm256_cvttps_epi32(rt);
1111             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1112             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1113             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1114             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1115             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1116             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1117
1118             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1119             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1120                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1121             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1122                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1123             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1124                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1125             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1126                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1127             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1128             Heps             = _mm256_mul_ps(vfeps,H);
1129             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1130             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1131             velec            = _mm256_mul_ps(qq13,VV);
1132             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1133             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq13,FF),_mm256_mul_ps(vftabscale,rinv13)));
1134
1135             /* Update potential sum for this i atom from the interaction with this j atom. */
1136             velec            = _mm256_andnot_ps(dummy_mask,velec);
1137             velecsum         = _mm256_add_ps(velecsum,velec);
1138
1139             fscal            = felec;
1140
1141             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1142
1143             /* Calculate temporary vectorial force */
1144             tx               = _mm256_mul_ps(fscal,dx13);
1145             ty               = _mm256_mul_ps(fscal,dy13);
1146             tz               = _mm256_mul_ps(fscal,dz13);
1147
1148             /* Update vectorial force */
1149             fix1             = _mm256_add_ps(fix1,tx);
1150             fiy1             = _mm256_add_ps(fiy1,ty);
1151             fiz1             = _mm256_add_ps(fiz1,tz);
1152
1153             fjx3             = _mm256_add_ps(fjx3,tx);
1154             fjy3             = _mm256_add_ps(fjy3,ty);
1155             fjz3             = _mm256_add_ps(fjz3,tz);
1156
1157             /**************************
1158              * CALCULATE INTERACTIONS *
1159              **************************/
1160
1161             r21              = _mm256_mul_ps(rsq21,rinv21);
1162             r21              = _mm256_andnot_ps(dummy_mask,r21);
1163
1164             /* Calculate table index by multiplying r with table scale and truncate to integer */
1165             rt               = _mm256_mul_ps(r21,vftabscale);
1166             vfitab           = _mm256_cvttps_epi32(rt);
1167             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1168             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1169             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1170             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1171             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1172             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1173
1174             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1175             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1176                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1177             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1178                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1179             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1180                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1181             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1182                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1183             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1184             Heps             = _mm256_mul_ps(vfeps,H);
1185             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1186             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1187             velec            = _mm256_mul_ps(qq21,VV);
1188             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1189             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq21,FF),_mm256_mul_ps(vftabscale,rinv21)));
1190
1191             /* Update potential sum for this i atom from the interaction with this j atom. */
1192             velec            = _mm256_andnot_ps(dummy_mask,velec);
1193             velecsum         = _mm256_add_ps(velecsum,velec);
1194
1195             fscal            = felec;
1196
1197             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1198
1199             /* Calculate temporary vectorial force */
1200             tx               = _mm256_mul_ps(fscal,dx21);
1201             ty               = _mm256_mul_ps(fscal,dy21);
1202             tz               = _mm256_mul_ps(fscal,dz21);
1203
1204             /* Update vectorial force */
1205             fix2             = _mm256_add_ps(fix2,tx);
1206             fiy2             = _mm256_add_ps(fiy2,ty);
1207             fiz2             = _mm256_add_ps(fiz2,tz);
1208
1209             fjx1             = _mm256_add_ps(fjx1,tx);
1210             fjy1             = _mm256_add_ps(fjy1,ty);
1211             fjz1             = _mm256_add_ps(fjz1,tz);
1212
1213             /**************************
1214              * CALCULATE INTERACTIONS *
1215              **************************/
1216
1217             r22              = _mm256_mul_ps(rsq22,rinv22);
1218             r22              = _mm256_andnot_ps(dummy_mask,r22);
1219
1220             /* Calculate table index by multiplying r with table scale and truncate to integer */
1221             rt               = _mm256_mul_ps(r22,vftabscale);
1222             vfitab           = _mm256_cvttps_epi32(rt);
1223             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1224             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1225             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1226             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1227             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1228             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1229
1230             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1231             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1232                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1233             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1234                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1235             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1236                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1237             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1238                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1239             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1240             Heps             = _mm256_mul_ps(vfeps,H);
1241             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1242             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1243             velec            = _mm256_mul_ps(qq22,VV);
1244             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1245             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq22,FF),_mm256_mul_ps(vftabscale,rinv22)));
1246
1247             /* Update potential sum for this i atom from the interaction with this j atom. */
1248             velec            = _mm256_andnot_ps(dummy_mask,velec);
1249             velecsum         = _mm256_add_ps(velecsum,velec);
1250
1251             fscal            = felec;
1252
1253             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1254
1255             /* Calculate temporary vectorial force */
1256             tx               = _mm256_mul_ps(fscal,dx22);
1257             ty               = _mm256_mul_ps(fscal,dy22);
1258             tz               = _mm256_mul_ps(fscal,dz22);
1259
1260             /* Update vectorial force */
1261             fix2             = _mm256_add_ps(fix2,tx);
1262             fiy2             = _mm256_add_ps(fiy2,ty);
1263             fiz2             = _mm256_add_ps(fiz2,tz);
1264
1265             fjx2             = _mm256_add_ps(fjx2,tx);
1266             fjy2             = _mm256_add_ps(fjy2,ty);
1267             fjz2             = _mm256_add_ps(fjz2,tz);
1268
1269             /**************************
1270              * CALCULATE INTERACTIONS *
1271              **************************/
1272
1273             r23              = _mm256_mul_ps(rsq23,rinv23);
1274             r23              = _mm256_andnot_ps(dummy_mask,r23);
1275
1276             /* Calculate table index by multiplying r with table scale and truncate to integer */
1277             rt               = _mm256_mul_ps(r23,vftabscale);
1278             vfitab           = _mm256_cvttps_epi32(rt);
1279             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1280             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1281             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1282             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1283             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1284             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1285
1286             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1287             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1288                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1289             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1290                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1291             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1292                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1293             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1294                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1295             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1296             Heps             = _mm256_mul_ps(vfeps,H);
1297             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1298             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1299             velec            = _mm256_mul_ps(qq23,VV);
1300             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1301             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq23,FF),_mm256_mul_ps(vftabscale,rinv23)));
1302
1303             /* Update potential sum for this i atom from the interaction with this j atom. */
1304             velec            = _mm256_andnot_ps(dummy_mask,velec);
1305             velecsum         = _mm256_add_ps(velecsum,velec);
1306
1307             fscal            = felec;
1308
1309             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1310
1311             /* Calculate temporary vectorial force */
1312             tx               = _mm256_mul_ps(fscal,dx23);
1313             ty               = _mm256_mul_ps(fscal,dy23);
1314             tz               = _mm256_mul_ps(fscal,dz23);
1315
1316             /* Update vectorial force */
1317             fix2             = _mm256_add_ps(fix2,tx);
1318             fiy2             = _mm256_add_ps(fiy2,ty);
1319             fiz2             = _mm256_add_ps(fiz2,tz);
1320
1321             fjx3             = _mm256_add_ps(fjx3,tx);
1322             fjy3             = _mm256_add_ps(fjy3,ty);
1323             fjz3             = _mm256_add_ps(fjz3,tz);
1324
1325             /**************************
1326              * CALCULATE INTERACTIONS *
1327              **************************/
1328
1329             r31              = _mm256_mul_ps(rsq31,rinv31);
1330             r31              = _mm256_andnot_ps(dummy_mask,r31);
1331
1332             /* Calculate table index by multiplying r with table scale and truncate to integer */
1333             rt               = _mm256_mul_ps(r31,vftabscale);
1334             vfitab           = _mm256_cvttps_epi32(rt);
1335             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1336             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1337             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1338             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1339             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1340             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1341
1342             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1343             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1344                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1345             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1346                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1347             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1348                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1349             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1350                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1351             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1352             Heps             = _mm256_mul_ps(vfeps,H);
1353             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1354             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1355             velec            = _mm256_mul_ps(qq31,VV);
1356             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1357             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq31,FF),_mm256_mul_ps(vftabscale,rinv31)));
1358
1359             /* Update potential sum for this i atom from the interaction with this j atom. */
1360             velec            = _mm256_andnot_ps(dummy_mask,velec);
1361             velecsum         = _mm256_add_ps(velecsum,velec);
1362
1363             fscal            = felec;
1364
1365             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1366
1367             /* Calculate temporary vectorial force */
1368             tx               = _mm256_mul_ps(fscal,dx31);
1369             ty               = _mm256_mul_ps(fscal,dy31);
1370             tz               = _mm256_mul_ps(fscal,dz31);
1371
1372             /* Update vectorial force */
1373             fix3             = _mm256_add_ps(fix3,tx);
1374             fiy3             = _mm256_add_ps(fiy3,ty);
1375             fiz3             = _mm256_add_ps(fiz3,tz);
1376
1377             fjx1             = _mm256_add_ps(fjx1,tx);
1378             fjy1             = _mm256_add_ps(fjy1,ty);
1379             fjz1             = _mm256_add_ps(fjz1,tz);
1380
1381             /**************************
1382              * CALCULATE INTERACTIONS *
1383              **************************/
1384
1385             r32              = _mm256_mul_ps(rsq32,rinv32);
1386             r32              = _mm256_andnot_ps(dummy_mask,r32);
1387
1388             /* Calculate table index by multiplying r with table scale and truncate to integer */
1389             rt               = _mm256_mul_ps(r32,vftabscale);
1390             vfitab           = _mm256_cvttps_epi32(rt);
1391             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1392             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1393             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1394             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1395             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1396             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1397
1398             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1399             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1400                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1401             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1402                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1403             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1404                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1405             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1406                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1407             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1408             Heps             = _mm256_mul_ps(vfeps,H);
1409             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1410             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1411             velec            = _mm256_mul_ps(qq32,VV);
1412             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1413             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq32,FF),_mm256_mul_ps(vftabscale,rinv32)));
1414
1415             /* Update potential sum for this i atom from the interaction with this j atom. */
1416             velec            = _mm256_andnot_ps(dummy_mask,velec);
1417             velecsum         = _mm256_add_ps(velecsum,velec);
1418
1419             fscal            = felec;
1420
1421             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1422
1423             /* Calculate temporary vectorial force */
1424             tx               = _mm256_mul_ps(fscal,dx32);
1425             ty               = _mm256_mul_ps(fscal,dy32);
1426             tz               = _mm256_mul_ps(fscal,dz32);
1427
1428             /* Update vectorial force */
1429             fix3             = _mm256_add_ps(fix3,tx);
1430             fiy3             = _mm256_add_ps(fiy3,ty);
1431             fiz3             = _mm256_add_ps(fiz3,tz);
1432
1433             fjx2             = _mm256_add_ps(fjx2,tx);
1434             fjy2             = _mm256_add_ps(fjy2,ty);
1435             fjz2             = _mm256_add_ps(fjz2,tz);
1436
1437             /**************************
1438              * CALCULATE INTERACTIONS *
1439              **************************/
1440
1441             r33              = _mm256_mul_ps(rsq33,rinv33);
1442             r33              = _mm256_andnot_ps(dummy_mask,r33);
1443
1444             /* Calculate table index by multiplying r with table scale and truncate to integer */
1445             rt               = _mm256_mul_ps(r33,vftabscale);
1446             vfitab           = _mm256_cvttps_epi32(rt);
1447             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1448             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1449             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1450             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1451             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1452             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1453
1454             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1455             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1456                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1457             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1458                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1459             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1460                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1461             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1462                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1463             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1464             Heps             = _mm256_mul_ps(vfeps,H);
1465             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1466             VV               = _mm256_add_ps(Y,_mm256_mul_ps(vfeps,Fp));
1467             velec            = _mm256_mul_ps(qq33,VV);
1468             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1469             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq33,FF),_mm256_mul_ps(vftabscale,rinv33)));
1470
1471             /* Update potential sum for this i atom from the interaction with this j atom. */
1472             velec            = _mm256_andnot_ps(dummy_mask,velec);
1473             velecsum         = _mm256_add_ps(velecsum,velec);
1474
1475             fscal            = felec;
1476
1477             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
1478
1479             /* Calculate temporary vectorial force */
1480             tx               = _mm256_mul_ps(fscal,dx33);
1481             ty               = _mm256_mul_ps(fscal,dy33);
1482             tz               = _mm256_mul_ps(fscal,dz33);
1483
1484             /* Update vectorial force */
1485             fix3             = _mm256_add_ps(fix3,tx);
1486             fiy3             = _mm256_add_ps(fiy3,ty);
1487             fiz3             = _mm256_add_ps(fiz3,tz);
1488
1489             fjx3             = _mm256_add_ps(fjx3,tx);
1490             fjy3             = _mm256_add_ps(fjy3,ty);
1491             fjz3             = _mm256_add_ps(fjz3,tz);
1492
1493             fjptrA             = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
1494             fjptrB             = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
1495             fjptrC             = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
1496             fjptrD             = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
1497             fjptrE             = (jnrlistE>=0) ? f+j_coord_offsetE : scratch;
1498             fjptrF             = (jnrlistF>=0) ? f+j_coord_offsetF : scratch;
1499             fjptrG             = (jnrlistG>=0) ? f+j_coord_offsetG : scratch;
1500             fjptrH             = (jnrlistH>=0) ? f+j_coord_offsetH : scratch;
1501
1502             gmx_mm256_decrement_4rvec_8ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjptrE,fjptrF,fjptrG,fjptrH,
1503                                                       fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,
1504                                                       fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
1505
1506             /* Inner loop uses 431 flops */
1507         }
1508
1509         /* End of innermost loop */
1510
1511         gmx_mm256_update_iforce_4atom_swizzle_ps(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
1512                                                  f+i_coord_offset,fshift+i_shift_offset);
1513
1514         ggid                        = gid[iidx];
1515         /* Update potential energies */
1516         gmx_mm256_update_1pot_ps(velecsum,kernel_data->energygrp_elec+ggid);
1517         gmx_mm256_update_1pot_ps(vvdwsum,kernel_data->energygrp_vdw+ggid);
1518
1519         /* Increment number of inner iterations */
1520         inneriter                  += j_index_end - j_index_start;
1521
1522         /* Outer loop uses 26 flops */
1523     }
1524
1525     /* Increment number of outer iterations */
1526     outeriter        += nri;
1527
1528     /* Update outer/inner flops */
1529
1530     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4W4_VF,outeriter*26 + inneriter*431);
1531 }
1532 /*
1533  * Gromacs nonbonded kernel:   nb_kernel_ElecCSTab_VdwLJ_GeomW4W4_F_avx_256_single
1534  * Electrostatics interaction: CubicSplineTable
1535  * VdW interaction:            LennardJones
1536  * Geometry:                   Water4-Water4
1537  * Calculate force/pot:        Force
1538  */
1539 void
1540 nb_kernel_ElecCSTab_VdwLJ_GeomW4W4_F_avx_256_single
1541                     (t_nblist                    * gmx_restrict       nlist,
1542                      rvec                        * gmx_restrict          xx,
1543                      rvec                        * gmx_restrict          ff,
1544                      t_forcerec                  * gmx_restrict          fr,
1545                      t_mdatoms                   * gmx_restrict     mdatoms,
1546                      nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
1547                      t_nrnb                      * gmx_restrict        nrnb)
1548 {
1549     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
1550      * just 0 for non-waters.
1551      * Suffixes A,B,C,D,E,F,G,H refer to j loop unrolling done with AVX, e.g. for the eight different
1552      * jnr indices corresponding to data put in the four positions in the SIMD register.
1553      */
1554     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
1555     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
1556     int              jnrA,jnrB,jnrC,jnrD;
1557     int              jnrE,jnrF,jnrG,jnrH;
1558     int              jnrlistA,jnrlistB,jnrlistC,jnrlistD;
1559     int              jnrlistE,jnrlistF,jnrlistG,jnrlistH;
1560     int              j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
1561     int              j_coord_offsetE,j_coord_offsetF,j_coord_offsetG,j_coord_offsetH;
1562     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
1563     real             rcutoff_scalar;
1564     real             *shiftvec,*fshift,*x,*f;
1565     real             *fjptrA,*fjptrB,*fjptrC,*fjptrD,*fjptrE,*fjptrF,*fjptrG,*fjptrH;
1566     real             scratch[4*DIM];
1567     __m256           tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
1568     real *           vdwioffsetptr0;
1569     __m256           ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
1570     real *           vdwioffsetptr1;
1571     __m256           ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
1572     real *           vdwioffsetptr2;
1573     __m256           ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
1574     real *           vdwioffsetptr3;
1575     __m256           ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
1576     int              vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D,vdwjidx0E,vdwjidx0F,vdwjidx0G,vdwjidx0H;
1577     __m256           jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
1578     int              vdwjidx1A,vdwjidx1B,vdwjidx1C,vdwjidx1D,vdwjidx1E,vdwjidx1F,vdwjidx1G,vdwjidx1H;
1579     __m256           jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
1580     int              vdwjidx2A,vdwjidx2B,vdwjidx2C,vdwjidx2D,vdwjidx2E,vdwjidx2F,vdwjidx2G,vdwjidx2H;
1581     __m256           jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
1582     int              vdwjidx3A,vdwjidx3B,vdwjidx3C,vdwjidx3D,vdwjidx3E,vdwjidx3F,vdwjidx3G,vdwjidx3H;
1583     __m256           jx3,jy3,jz3,fjx3,fjy3,fjz3,jq3,isaj3;
1584     __m256           dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
1585     __m256           dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
1586     __m256           dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
1587     __m256           dx13,dy13,dz13,rsq13,rinv13,rinvsq13,r13,qq13,c6_13,c12_13;
1588     __m256           dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
1589     __m256           dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
1590     __m256           dx23,dy23,dz23,rsq23,rinv23,rinvsq23,r23,qq23,c6_23,c12_23;
1591     __m256           dx31,dy31,dz31,rsq31,rinv31,rinvsq31,r31,qq31,c6_31,c12_31;
1592     __m256           dx32,dy32,dz32,rsq32,rinv32,rinvsq32,r32,qq32,c6_32,c12_32;
1593     __m256           dx33,dy33,dz33,rsq33,rinv33,rinvsq33,r33,qq33,c6_33,c12_33;
1594     __m256           velec,felec,velecsum,facel,crf,krf,krf2;
1595     real             *charge;
1596     int              nvdwtype;
1597     __m256           rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
1598     int              *vdwtype;
1599     real             *vdwparam;
1600     __m256           one_sixth   = _mm256_set1_ps(1.0/6.0);
1601     __m256           one_twelfth = _mm256_set1_ps(1.0/12.0);
1602     __m256i          vfitab;
1603     __m128i          vfitab_lo,vfitab_hi;
1604     __m128i          ifour       = _mm_set1_epi32(4);
1605     __m256           rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF;
1606     real             *vftab;
1607     __m256           dummy_mask,cutoff_mask;
1608     __m256           signbit = _mm256_castsi256_ps( _mm256_set1_epi32(0x80000000) );
1609     __m256           one     = _mm256_set1_ps(1.0);
1610     __m256           two     = _mm256_set1_ps(2.0);
1611     x                = xx[0];
1612     f                = ff[0];
1613
1614     nri              = nlist->nri;
1615     iinr             = nlist->iinr;
1616     jindex           = nlist->jindex;
1617     jjnr             = nlist->jjnr;
1618     shiftidx         = nlist->shift;
1619     gid              = nlist->gid;
1620     shiftvec         = fr->shift_vec[0];
1621     fshift           = fr->fshift[0];
1622     facel            = _mm256_set1_ps(fr->epsfac);
1623     charge           = mdatoms->chargeA;
1624     nvdwtype         = fr->ntype;
1625     vdwparam         = fr->nbfp;
1626     vdwtype          = mdatoms->typeA;
1627
1628     vftab            = kernel_data->table_elec->data;
1629     vftabscale       = _mm256_set1_ps(kernel_data->table_elec->scale);
1630
1631     /* Setup water-specific parameters */
1632     inr              = nlist->iinr[0];
1633     iq1              = _mm256_mul_ps(facel,_mm256_set1_ps(charge[inr+1]));
1634     iq2              = _mm256_mul_ps(facel,_mm256_set1_ps(charge[inr+2]));
1635     iq3              = _mm256_mul_ps(facel,_mm256_set1_ps(charge[inr+3]));
1636     vdwioffsetptr0   = vdwparam+2*nvdwtype*vdwtype[inr+0];
1637
1638     jq1              = _mm256_set1_ps(charge[inr+1]);
1639     jq2              = _mm256_set1_ps(charge[inr+2]);
1640     jq3              = _mm256_set1_ps(charge[inr+3]);
1641     vdwjidx0A        = 2*vdwtype[inr+0];
1642     c6_00            = _mm256_set1_ps(vdwioffsetptr0[vdwjidx0A]);
1643     c12_00           = _mm256_set1_ps(vdwioffsetptr0[vdwjidx0A+1]);
1644     qq11             = _mm256_mul_ps(iq1,jq1);
1645     qq12             = _mm256_mul_ps(iq1,jq2);
1646     qq13             = _mm256_mul_ps(iq1,jq3);
1647     qq21             = _mm256_mul_ps(iq2,jq1);
1648     qq22             = _mm256_mul_ps(iq2,jq2);
1649     qq23             = _mm256_mul_ps(iq2,jq3);
1650     qq31             = _mm256_mul_ps(iq3,jq1);
1651     qq32             = _mm256_mul_ps(iq3,jq2);
1652     qq33             = _mm256_mul_ps(iq3,jq3);
1653
1654     /* Avoid stupid compiler warnings */
1655     jnrA = jnrB = jnrC = jnrD = jnrE = jnrF = jnrG = jnrH = 0;
1656     j_coord_offsetA = 0;
1657     j_coord_offsetB = 0;
1658     j_coord_offsetC = 0;
1659     j_coord_offsetD = 0;
1660     j_coord_offsetE = 0;
1661     j_coord_offsetF = 0;
1662     j_coord_offsetG = 0;
1663     j_coord_offsetH = 0;
1664
1665     outeriter        = 0;
1666     inneriter        = 0;
1667
1668     for(iidx=0;iidx<4*DIM;iidx++)
1669     {
1670         scratch[iidx] = 0.0;
1671     }
1672
1673     /* Start outer loop over neighborlists */
1674     for(iidx=0; iidx<nri; iidx++)
1675     {
1676         /* Load shift vector for this list */
1677         i_shift_offset   = DIM*shiftidx[iidx];
1678
1679         /* Load limits for loop over neighbors */
1680         j_index_start    = jindex[iidx];
1681         j_index_end      = jindex[iidx+1];
1682
1683         /* Get outer coordinate index */
1684         inr              = iinr[iidx];
1685         i_coord_offset   = DIM*inr;
1686
1687         /* Load i particle coords and add shift vector */
1688         gmx_mm256_load_shift_and_4rvec_broadcast_ps(shiftvec+i_shift_offset,x+i_coord_offset,
1689                                                     &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
1690
1691         fix0             = _mm256_setzero_ps();
1692         fiy0             = _mm256_setzero_ps();
1693         fiz0             = _mm256_setzero_ps();
1694         fix1             = _mm256_setzero_ps();
1695         fiy1             = _mm256_setzero_ps();
1696         fiz1             = _mm256_setzero_ps();
1697         fix2             = _mm256_setzero_ps();
1698         fiy2             = _mm256_setzero_ps();
1699         fiz2             = _mm256_setzero_ps();
1700         fix3             = _mm256_setzero_ps();
1701         fiy3             = _mm256_setzero_ps();
1702         fiz3             = _mm256_setzero_ps();
1703
1704         /* Start inner kernel loop */
1705         for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+7]>=0; jidx+=8)
1706         {
1707
1708             /* Get j neighbor index, and coordinate index */
1709             jnrA             = jjnr[jidx];
1710             jnrB             = jjnr[jidx+1];
1711             jnrC             = jjnr[jidx+2];
1712             jnrD             = jjnr[jidx+3];
1713             jnrE             = jjnr[jidx+4];
1714             jnrF             = jjnr[jidx+5];
1715             jnrG             = jjnr[jidx+6];
1716             jnrH             = jjnr[jidx+7];
1717             j_coord_offsetA  = DIM*jnrA;
1718             j_coord_offsetB  = DIM*jnrB;
1719             j_coord_offsetC  = DIM*jnrC;
1720             j_coord_offsetD  = DIM*jnrD;
1721             j_coord_offsetE  = DIM*jnrE;
1722             j_coord_offsetF  = DIM*jnrF;
1723             j_coord_offsetG  = DIM*jnrG;
1724             j_coord_offsetH  = DIM*jnrH;
1725
1726             /* load j atom coordinates */
1727             gmx_mm256_load_4rvec_8ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
1728                                                  x+j_coord_offsetC,x+j_coord_offsetD,
1729                                                  x+j_coord_offsetE,x+j_coord_offsetF,
1730                                                  x+j_coord_offsetG,x+j_coord_offsetH,
1731                                                  &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
1732                                                  &jy2,&jz2,&jx3,&jy3,&jz3);
1733
1734             /* Calculate displacement vector */
1735             dx00             = _mm256_sub_ps(ix0,jx0);
1736             dy00             = _mm256_sub_ps(iy0,jy0);
1737             dz00             = _mm256_sub_ps(iz0,jz0);
1738             dx11             = _mm256_sub_ps(ix1,jx1);
1739             dy11             = _mm256_sub_ps(iy1,jy1);
1740             dz11             = _mm256_sub_ps(iz1,jz1);
1741             dx12             = _mm256_sub_ps(ix1,jx2);
1742             dy12             = _mm256_sub_ps(iy1,jy2);
1743             dz12             = _mm256_sub_ps(iz1,jz2);
1744             dx13             = _mm256_sub_ps(ix1,jx3);
1745             dy13             = _mm256_sub_ps(iy1,jy3);
1746             dz13             = _mm256_sub_ps(iz1,jz3);
1747             dx21             = _mm256_sub_ps(ix2,jx1);
1748             dy21             = _mm256_sub_ps(iy2,jy1);
1749             dz21             = _mm256_sub_ps(iz2,jz1);
1750             dx22             = _mm256_sub_ps(ix2,jx2);
1751             dy22             = _mm256_sub_ps(iy2,jy2);
1752             dz22             = _mm256_sub_ps(iz2,jz2);
1753             dx23             = _mm256_sub_ps(ix2,jx3);
1754             dy23             = _mm256_sub_ps(iy2,jy3);
1755             dz23             = _mm256_sub_ps(iz2,jz3);
1756             dx31             = _mm256_sub_ps(ix3,jx1);
1757             dy31             = _mm256_sub_ps(iy3,jy1);
1758             dz31             = _mm256_sub_ps(iz3,jz1);
1759             dx32             = _mm256_sub_ps(ix3,jx2);
1760             dy32             = _mm256_sub_ps(iy3,jy2);
1761             dz32             = _mm256_sub_ps(iz3,jz2);
1762             dx33             = _mm256_sub_ps(ix3,jx3);
1763             dy33             = _mm256_sub_ps(iy3,jy3);
1764             dz33             = _mm256_sub_ps(iz3,jz3);
1765
1766             /* Calculate squared distance and things based on it */
1767             rsq00            = gmx_mm256_calc_rsq_ps(dx00,dy00,dz00);
1768             rsq11            = gmx_mm256_calc_rsq_ps(dx11,dy11,dz11);
1769             rsq12            = gmx_mm256_calc_rsq_ps(dx12,dy12,dz12);
1770             rsq13            = gmx_mm256_calc_rsq_ps(dx13,dy13,dz13);
1771             rsq21            = gmx_mm256_calc_rsq_ps(dx21,dy21,dz21);
1772             rsq22            = gmx_mm256_calc_rsq_ps(dx22,dy22,dz22);
1773             rsq23            = gmx_mm256_calc_rsq_ps(dx23,dy23,dz23);
1774             rsq31            = gmx_mm256_calc_rsq_ps(dx31,dy31,dz31);
1775             rsq32            = gmx_mm256_calc_rsq_ps(dx32,dy32,dz32);
1776             rsq33            = gmx_mm256_calc_rsq_ps(dx33,dy33,dz33);
1777
1778             rinv11           = gmx_mm256_invsqrt_ps(rsq11);
1779             rinv12           = gmx_mm256_invsqrt_ps(rsq12);
1780             rinv13           = gmx_mm256_invsqrt_ps(rsq13);
1781             rinv21           = gmx_mm256_invsqrt_ps(rsq21);
1782             rinv22           = gmx_mm256_invsqrt_ps(rsq22);
1783             rinv23           = gmx_mm256_invsqrt_ps(rsq23);
1784             rinv31           = gmx_mm256_invsqrt_ps(rsq31);
1785             rinv32           = gmx_mm256_invsqrt_ps(rsq32);
1786             rinv33           = gmx_mm256_invsqrt_ps(rsq33);
1787
1788             rinvsq00         = gmx_mm256_inv_ps(rsq00);
1789
1790             fjx0             = _mm256_setzero_ps();
1791             fjy0             = _mm256_setzero_ps();
1792             fjz0             = _mm256_setzero_ps();
1793             fjx1             = _mm256_setzero_ps();
1794             fjy1             = _mm256_setzero_ps();
1795             fjz1             = _mm256_setzero_ps();
1796             fjx2             = _mm256_setzero_ps();
1797             fjy2             = _mm256_setzero_ps();
1798             fjz2             = _mm256_setzero_ps();
1799             fjx3             = _mm256_setzero_ps();
1800             fjy3             = _mm256_setzero_ps();
1801             fjz3             = _mm256_setzero_ps();
1802
1803             /**************************
1804              * CALCULATE INTERACTIONS *
1805              **************************/
1806
1807             /* LENNARD-JONES DISPERSION/REPULSION */
1808
1809             rinvsix          = _mm256_mul_ps(_mm256_mul_ps(rinvsq00,rinvsq00),rinvsq00);
1810             fvdw             = _mm256_mul_ps(_mm256_sub_ps(_mm256_mul_ps(c12_00,rinvsix),c6_00),_mm256_mul_ps(rinvsix,rinvsq00));
1811
1812             fscal            = fvdw;
1813
1814             /* Calculate temporary vectorial force */
1815             tx               = _mm256_mul_ps(fscal,dx00);
1816             ty               = _mm256_mul_ps(fscal,dy00);
1817             tz               = _mm256_mul_ps(fscal,dz00);
1818
1819             /* Update vectorial force */
1820             fix0             = _mm256_add_ps(fix0,tx);
1821             fiy0             = _mm256_add_ps(fiy0,ty);
1822             fiz0             = _mm256_add_ps(fiz0,tz);
1823
1824             fjx0             = _mm256_add_ps(fjx0,tx);
1825             fjy0             = _mm256_add_ps(fjy0,ty);
1826             fjz0             = _mm256_add_ps(fjz0,tz);
1827
1828             /**************************
1829              * CALCULATE INTERACTIONS *
1830              **************************/
1831
1832             r11              = _mm256_mul_ps(rsq11,rinv11);
1833
1834             /* Calculate table index by multiplying r with table scale and truncate to integer */
1835             rt               = _mm256_mul_ps(r11,vftabscale);
1836             vfitab           = _mm256_cvttps_epi32(rt);
1837             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1838             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1839             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1840             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1841             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1842             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1843
1844             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1845             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1846                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1847             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1848                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1849             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1850                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1851             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1852                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1853             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1854             Heps             = _mm256_mul_ps(vfeps,H);
1855             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1856             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1857             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq11,FF),_mm256_mul_ps(vftabscale,rinv11)));
1858
1859             fscal            = felec;
1860
1861             /* Calculate temporary vectorial force */
1862             tx               = _mm256_mul_ps(fscal,dx11);
1863             ty               = _mm256_mul_ps(fscal,dy11);
1864             tz               = _mm256_mul_ps(fscal,dz11);
1865
1866             /* Update vectorial force */
1867             fix1             = _mm256_add_ps(fix1,tx);
1868             fiy1             = _mm256_add_ps(fiy1,ty);
1869             fiz1             = _mm256_add_ps(fiz1,tz);
1870
1871             fjx1             = _mm256_add_ps(fjx1,tx);
1872             fjy1             = _mm256_add_ps(fjy1,ty);
1873             fjz1             = _mm256_add_ps(fjz1,tz);
1874
1875             /**************************
1876              * CALCULATE INTERACTIONS *
1877              **************************/
1878
1879             r12              = _mm256_mul_ps(rsq12,rinv12);
1880
1881             /* Calculate table index by multiplying r with table scale and truncate to integer */
1882             rt               = _mm256_mul_ps(r12,vftabscale);
1883             vfitab           = _mm256_cvttps_epi32(rt);
1884             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1885             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1886             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1887             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1888             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1889             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1890
1891             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1892             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1893                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1894             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1895                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1896             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1897                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1898             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1899                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1900             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1901             Heps             = _mm256_mul_ps(vfeps,H);
1902             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1903             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1904             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq12,FF),_mm256_mul_ps(vftabscale,rinv12)));
1905
1906             fscal            = felec;
1907
1908             /* Calculate temporary vectorial force */
1909             tx               = _mm256_mul_ps(fscal,dx12);
1910             ty               = _mm256_mul_ps(fscal,dy12);
1911             tz               = _mm256_mul_ps(fscal,dz12);
1912
1913             /* Update vectorial force */
1914             fix1             = _mm256_add_ps(fix1,tx);
1915             fiy1             = _mm256_add_ps(fiy1,ty);
1916             fiz1             = _mm256_add_ps(fiz1,tz);
1917
1918             fjx2             = _mm256_add_ps(fjx2,tx);
1919             fjy2             = _mm256_add_ps(fjy2,ty);
1920             fjz2             = _mm256_add_ps(fjz2,tz);
1921
1922             /**************************
1923              * CALCULATE INTERACTIONS *
1924              **************************/
1925
1926             r13              = _mm256_mul_ps(rsq13,rinv13);
1927
1928             /* Calculate table index by multiplying r with table scale and truncate to integer */
1929             rt               = _mm256_mul_ps(r13,vftabscale);
1930             vfitab           = _mm256_cvttps_epi32(rt);
1931             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1932             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1933             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1934             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1935             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1936             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1937
1938             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1939             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1940                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1941             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1942                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1943             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1944                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1945             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1946                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1947             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1948             Heps             = _mm256_mul_ps(vfeps,H);
1949             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1950             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1951             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq13,FF),_mm256_mul_ps(vftabscale,rinv13)));
1952
1953             fscal            = felec;
1954
1955             /* Calculate temporary vectorial force */
1956             tx               = _mm256_mul_ps(fscal,dx13);
1957             ty               = _mm256_mul_ps(fscal,dy13);
1958             tz               = _mm256_mul_ps(fscal,dz13);
1959
1960             /* Update vectorial force */
1961             fix1             = _mm256_add_ps(fix1,tx);
1962             fiy1             = _mm256_add_ps(fiy1,ty);
1963             fiz1             = _mm256_add_ps(fiz1,tz);
1964
1965             fjx3             = _mm256_add_ps(fjx3,tx);
1966             fjy3             = _mm256_add_ps(fjy3,ty);
1967             fjz3             = _mm256_add_ps(fjz3,tz);
1968
1969             /**************************
1970              * CALCULATE INTERACTIONS *
1971              **************************/
1972
1973             r21              = _mm256_mul_ps(rsq21,rinv21);
1974
1975             /* Calculate table index by multiplying r with table scale and truncate to integer */
1976             rt               = _mm256_mul_ps(r21,vftabscale);
1977             vfitab           = _mm256_cvttps_epi32(rt);
1978             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
1979             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
1980             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
1981             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
1982             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
1983             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
1984
1985             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1986             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
1987                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
1988             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
1989                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
1990             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
1991                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
1992             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
1993                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
1994             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
1995             Heps             = _mm256_mul_ps(vfeps,H);
1996             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
1997             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
1998             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq21,FF),_mm256_mul_ps(vftabscale,rinv21)));
1999
2000             fscal            = felec;
2001
2002             /* Calculate temporary vectorial force */
2003             tx               = _mm256_mul_ps(fscal,dx21);
2004             ty               = _mm256_mul_ps(fscal,dy21);
2005             tz               = _mm256_mul_ps(fscal,dz21);
2006
2007             /* Update vectorial force */
2008             fix2             = _mm256_add_ps(fix2,tx);
2009             fiy2             = _mm256_add_ps(fiy2,ty);
2010             fiz2             = _mm256_add_ps(fiz2,tz);
2011
2012             fjx1             = _mm256_add_ps(fjx1,tx);
2013             fjy1             = _mm256_add_ps(fjy1,ty);
2014             fjz1             = _mm256_add_ps(fjz1,tz);
2015
2016             /**************************
2017              * CALCULATE INTERACTIONS *
2018              **************************/
2019
2020             r22              = _mm256_mul_ps(rsq22,rinv22);
2021
2022             /* Calculate table index by multiplying r with table scale and truncate to integer */
2023             rt               = _mm256_mul_ps(r22,vftabscale);
2024             vfitab           = _mm256_cvttps_epi32(rt);
2025             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2026             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2027             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2028             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2029             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2030             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2031
2032             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2033             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2034                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2035             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2036                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2037             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2038                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2039             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2040                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2041             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2042             Heps             = _mm256_mul_ps(vfeps,H);
2043             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2044             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2045             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq22,FF),_mm256_mul_ps(vftabscale,rinv22)));
2046
2047             fscal            = felec;
2048
2049             /* Calculate temporary vectorial force */
2050             tx               = _mm256_mul_ps(fscal,dx22);
2051             ty               = _mm256_mul_ps(fscal,dy22);
2052             tz               = _mm256_mul_ps(fscal,dz22);
2053
2054             /* Update vectorial force */
2055             fix2             = _mm256_add_ps(fix2,tx);
2056             fiy2             = _mm256_add_ps(fiy2,ty);
2057             fiz2             = _mm256_add_ps(fiz2,tz);
2058
2059             fjx2             = _mm256_add_ps(fjx2,tx);
2060             fjy2             = _mm256_add_ps(fjy2,ty);
2061             fjz2             = _mm256_add_ps(fjz2,tz);
2062
2063             /**************************
2064              * CALCULATE INTERACTIONS *
2065              **************************/
2066
2067             r23              = _mm256_mul_ps(rsq23,rinv23);
2068
2069             /* Calculate table index by multiplying r with table scale and truncate to integer */
2070             rt               = _mm256_mul_ps(r23,vftabscale);
2071             vfitab           = _mm256_cvttps_epi32(rt);
2072             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2073             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2074             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2075             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2076             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2077             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2078
2079             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2080             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2081                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2082             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2083                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2084             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2085                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2086             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2087                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2088             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2089             Heps             = _mm256_mul_ps(vfeps,H);
2090             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2091             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2092             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq23,FF),_mm256_mul_ps(vftabscale,rinv23)));
2093
2094             fscal            = felec;
2095
2096             /* Calculate temporary vectorial force */
2097             tx               = _mm256_mul_ps(fscal,dx23);
2098             ty               = _mm256_mul_ps(fscal,dy23);
2099             tz               = _mm256_mul_ps(fscal,dz23);
2100
2101             /* Update vectorial force */
2102             fix2             = _mm256_add_ps(fix2,tx);
2103             fiy2             = _mm256_add_ps(fiy2,ty);
2104             fiz2             = _mm256_add_ps(fiz2,tz);
2105
2106             fjx3             = _mm256_add_ps(fjx3,tx);
2107             fjy3             = _mm256_add_ps(fjy3,ty);
2108             fjz3             = _mm256_add_ps(fjz3,tz);
2109
2110             /**************************
2111              * CALCULATE INTERACTIONS *
2112              **************************/
2113
2114             r31              = _mm256_mul_ps(rsq31,rinv31);
2115
2116             /* Calculate table index by multiplying r with table scale and truncate to integer */
2117             rt               = _mm256_mul_ps(r31,vftabscale);
2118             vfitab           = _mm256_cvttps_epi32(rt);
2119             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2120             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2121             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2122             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2123             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2124             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2125
2126             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2127             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2128                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2129             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2130                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2131             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2132                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2133             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2134                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2135             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2136             Heps             = _mm256_mul_ps(vfeps,H);
2137             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2138             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2139             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq31,FF),_mm256_mul_ps(vftabscale,rinv31)));
2140
2141             fscal            = felec;
2142
2143             /* Calculate temporary vectorial force */
2144             tx               = _mm256_mul_ps(fscal,dx31);
2145             ty               = _mm256_mul_ps(fscal,dy31);
2146             tz               = _mm256_mul_ps(fscal,dz31);
2147
2148             /* Update vectorial force */
2149             fix3             = _mm256_add_ps(fix3,tx);
2150             fiy3             = _mm256_add_ps(fiy3,ty);
2151             fiz3             = _mm256_add_ps(fiz3,tz);
2152
2153             fjx1             = _mm256_add_ps(fjx1,tx);
2154             fjy1             = _mm256_add_ps(fjy1,ty);
2155             fjz1             = _mm256_add_ps(fjz1,tz);
2156
2157             /**************************
2158              * CALCULATE INTERACTIONS *
2159              **************************/
2160
2161             r32              = _mm256_mul_ps(rsq32,rinv32);
2162
2163             /* Calculate table index by multiplying r with table scale and truncate to integer */
2164             rt               = _mm256_mul_ps(r32,vftabscale);
2165             vfitab           = _mm256_cvttps_epi32(rt);
2166             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2167             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2168             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2169             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2170             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2171             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2172
2173             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2174             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2175                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2176             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2177                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2178             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2179                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2180             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2181                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2182             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2183             Heps             = _mm256_mul_ps(vfeps,H);
2184             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2185             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2186             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq32,FF),_mm256_mul_ps(vftabscale,rinv32)));
2187
2188             fscal            = felec;
2189
2190             /* Calculate temporary vectorial force */
2191             tx               = _mm256_mul_ps(fscal,dx32);
2192             ty               = _mm256_mul_ps(fscal,dy32);
2193             tz               = _mm256_mul_ps(fscal,dz32);
2194
2195             /* Update vectorial force */
2196             fix3             = _mm256_add_ps(fix3,tx);
2197             fiy3             = _mm256_add_ps(fiy3,ty);
2198             fiz3             = _mm256_add_ps(fiz3,tz);
2199
2200             fjx2             = _mm256_add_ps(fjx2,tx);
2201             fjy2             = _mm256_add_ps(fjy2,ty);
2202             fjz2             = _mm256_add_ps(fjz2,tz);
2203
2204             /**************************
2205              * CALCULATE INTERACTIONS *
2206              **************************/
2207
2208             r33              = _mm256_mul_ps(rsq33,rinv33);
2209
2210             /* Calculate table index by multiplying r with table scale and truncate to integer */
2211             rt               = _mm256_mul_ps(r33,vftabscale);
2212             vfitab           = _mm256_cvttps_epi32(rt);
2213             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2214             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2215             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2216             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2217             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2218             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2219
2220             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2221             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2222                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2223             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2224                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2225             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2226                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2227             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2228                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2229             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2230             Heps             = _mm256_mul_ps(vfeps,H);
2231             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2232             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2233             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq33,FF),_mm256_mul_ps(vftabscale,rinv33)));
2234
2235             fscal            = felec;
2236
2237             /* Calculate temporary vectorial force */
2238             tx               = _mm256_mul_ps(fscal,dx33);
2239             ty               = _mm256_mul_ps(fscal,dy33);
2240             tz               = _mm256_mul_ps(fscal,dz33);
2241
2242             /* Update vectorial force */
2243             fix3             = _mm256_add_ps(fix3,tx);
2244             fiy3             = _mm256_add_ps(fiy3,ty);
2245             fiz3             = _mm256_add_ps(fiz3,tz);
2246
2247             fjx3             = _mm256_add_ps(fjx3,tx);
2248             fjy3             = _mm256_add_ps(fjy3,ty);
2249             fjz3             = _mm256_add_ps(fjz3,tz);
2250
2251             fjptrA             = f+j_coord_offsetA;
2252             fjptrB             = f+j_coord_offsetB;
2253             fjptrC             = f+j_coord_offsetC;
2254             fjptrD             = f+j_coord_offsetD;
2255             fjptrE             = f+j_coord_offsetE;
2256             fjptrF             = f+j_coord_offsetF;
2257             fjptrG             = f+j_coord_offsetG;
2258             fjptrH             = f+j_coord_offsetH;
2259
2260             gmx_mm256_decrement_4rvec_8ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjptrE,fjptrF,fjptrG,fjptrH,
2261                                                       fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,
2262                                                       fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
2263
2264             /* Inner loop uses 381 flops */
2265         }
2266
2267         if(jidx<j_index_end)
2268         {
2269
2270             /* Get j neighbor index, and coordinate index */
2271             jnrlistA         = jjnr[jidx];
2272             jnrlistB         = jjnr[jidx+1];
2273             jnrlistC         = jjnr[jidx+2];
2274             jnrlistD         = jjnr[jidx+3];
2275             jnrlistE         = jjnr[jidx+4];
2276             jnrlistF         = jjnr[jidx+5];
2277             jnrlistG         = jjnr[jidx+6];
2278             jnrlistH         = jjnr[jidx+7];
2279             /* Sign of each element will be negative for non-real atoms.
2280              * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
2281              * so use it as val = _mm_andnot_ps(mask,val) to clear dummy entries.
2282              */
2283             dummy_mask = gmx_mm256_set_m128(gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx+4)),_mm_setzero_si128())),
2284                                             gmx_mm_castsi128_ps(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128())));
2285
2286             jnrA       = (jnrlistA>=0) ? jnrlistA : 0;
2287             jnrB       = (jnrlistB>=0) ? jnrlistB : 0;
2288             jnrC       = (jnrlistC>=0) ? jnrlistC : 0;
2289             jnrD       = (jnrlistD>=0) ? jnrlistD : 0;
2290             jnrE       = (jnrlistE>=0) ? jnrlistE : 0;
2291             jnrF       = (jnrlistF>=0) ? jnrlistF : 0;
2292             jnrG       = (jnrlistG>=0) ? jnrlistG : 0;
2293             jnrH       = (jnrlistH>=0) ? jnrlistH : 0;
2294             j_coord_offsetA  = DIM*jnrA;
2295             j_coord_offsetB  = DIM*jnrB;
2296             j_coord_offsetC  = DIM*jnrC;
2297             j_coord_offsetD  = DIM*jnrD;
2298             j_coord_offsetE  = DIM*jnrE;
2299             j_coord_offsetF  = DIM*jnrF;
2300             j_coord_offsetG  = DIM*jnrG;
2301             j_coord_offsetH  = DIM*jnrH;
2302
2303             /* load j atom coordinates */
2304             gmx_mm256_load_4rvec_8ptr_swizzle_ps(x+j_coord_offsetA,x+j_coord_offsetB,
2305                                                  x+j_coord_offsetC,x+j_coord_offsetD,
2306                                                  x+j_coord_offsetE,x+j_coord_offsetF,
2307                                                  x+j_coord_offsetG,x+j_coord_offsetH,
2308                                                  &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,
2309                                                  &jy2,&jz2,&jx3,&jy3,&jz3);
2310
2311             /* Calculate displacement vector */
2312             dx00             = _mm256_sub_ps(ix0,jx0);
2313             dy00             = _mm256_sub_ps(iy0,jy0);
2314             dz00             = _mm256_sub_ps(iz0,jz0);
2315             dx11             = _mm256_sub_ps(ix1,jx1);
2316             dy11             = _mm256_sub_ps(iy1,jy1);
2317             dz11             = _mm256_sub_ps(iz1,jz1);
2318             dx12             = _mm256_sub_ps(ix1,jx2);
2319             dy12             = _mm256_sub_ps(iy1,jy2);
2320             dz12             = _mm256_sub_ps(iz1,jz2);
2321             dx13             = _mm256_sub_ps(ix1,jx3);
2322             dy13             = _mm256_sub_ps(iy1,jy3);
2323             dz13             = _mm256_sub_ps(iz1,jz3);
2324             dx21             = _mm256_sub_ps(ix2,jx1);
2325             dy21             = _mm256_sub_ps(iy2,jy1);
2326             dz21             = _mm256_sub_ps(iz2,jz1);
2327             dx22             = _mm256_sub_ps(ix2,jx2);
2328             dy22             = _mm256_sub_ps(iy2,jy2);
2329             dz22             = _mm256_sub_ps(iz2,jz2);
2330             dx23             = _mm256_sub_ps(ix2,jx3);
2331             dy23             = _mm256_sub_ps(iy2,jy3);
2332             dz23             = _mm256_sub_ps(iz2,jz3);
2333             dx31             = _mm256_sub_ps(ix3,jx1);
2334             dy31             = _mm256_sub_ps(iy3,jy1);
2335             dz31             = _mm256_sub_ps(iz3,jz1);
2336             dx32             = _mm256_sub_ps(ix3,jx2);
2337             dy32             = _mm256_sub_ps(iy3,jy2);
2338             dz32             = _mm256_sub_ps(iz3,jz2);
2339             dx33             = _mm256_sub_ps(ix3,jx3);
2340             dy33             = _mm256_sub_ps(iy3,jy3);
2341             dz33             = _mm256_sub_ps(iz3,jz3);
2342
2343             /* Calculate squared distance and things based on it */
2344             rsq00            = gmx_mm256_calc_rsq_ps(dx00,dy00,dz00);
2345             rsq11            = gmx_mm256_calc_rsq_ps(dx11,dy11,dz11);
2346             rsq12            = gmx_mm256_calc_rsq_ps(dx12,dy12,dz12);
2347             rsq13            = gmx_mm256_calc_rsq_ps(dx13,dy13,dz13);
2348             rsq21            = gmx_mm256_calc_rsq_ps(dx21,dy21,dz21);
2349             rsq22            = gmx_mm256_calc_rsq_ps(dx22,dy22,dz22);
2350             rsq23            = gmx_mm256_calc_rsq_ps(dx23,dy23,dz23);
2351             rsq31            = gmx_mm256_calc_rsq_ps(dx31,dy31,dz31);
2352             rsq32            = gmx_mm256_calc_rsq_ps(dx32,dy32,dz32);
2353             rsq33            = gmx_mm256_calc_rsq_ps(dx33,dy33,dz33);
2354
2355             rinv11           = gmx_mm256_invsqrt_ps(rsq11);
2356             rinv12           = gmx_mm256_invsqrt_ps(rsq12);
2357             rinv13           = gmx_mm256_invsqrt_ps(rsq13);
2358             rinv21           = gmx_mm256_invsqrt_ps(rsq21);
2359             rinv22           = gmx_mm256_invsqrt_ps(rsq22);
2360             rinv23           = gmx_mm256_invsqrt_ps(rsq23);
2361             rinv31           = gmx_mm256_invsqrt_ps(rsq31);
2362             rinv32           = gmx_mm256_invsqrt_ps(rsq32);
2363             rinv33           = gmx_mm256_invsqrt_ps(rsq33);
2364
2365             rinvsq00         = gmx_mm256_inv_ps(rsq00);
2366
2367             fjx0             = _mm256_setzero_ps();
2368             fjy0             = _mm256_setzero_ps();
2369             fjz0             = _mm256_setzero_ps();
2370             fjx1             = _mm256_setzero_ps();
2371             fjy1             = _mm256_setzero_ps();
2372             fjz1             = _mm256_setzero_ps();
2373             fjx2             = _mm256_setzero_ps();
2374             fjy2             = _mm256_setzero_ps();
2375             fjz2             = _mm256_setzero_ps();
2376             fjx3             = _mm256_setzero_ps();
2377             fjy3             = _mm256_setzero_ps();
2378             fjz3             = _mm256_setzero_ps();
2379
2380             /**************************
2381              * CALCULATE INTERACTIONS *
2382              **************************/
2383
2384             /* LENNARD-JONES DISPERSION/REPULSION */
2385
2386             rinvsix          = _mm256_mul_ps(_mm256_mul_ps(rinvsq00,rinvsq00),rinvsq00);
2387             fvdw             = _mm256_mul_ps(_mm256_sub_ps(_mm256_mul_ps(c12_00,rinvsix),c6_00),_mm256_mul_ps(rinvsix,rinvsq00));
2388
2389             fscal            = fvdw;
2390
2391             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2392
2393             /* Calculate temporary vectorial force */
2394             tx               = _mm256_mul_ps(fscal,dx00);
2395             ty               = _mm256_mul_ps(fscal,dy00);
2396             tz               = _mm256_mul_ps(fscal,dz00);
2397
2398             /* Update vectorial force */
2399             fix0             = _mm256_add_ps(fix0,tx);
2400             fiy0             = _mm256_add_ps(fiy0,ty);
2401             fiz0             = _mm256_add_ps(fiz0,tz);
2402
2403             fjx0             = _mm256_add_ps(fjx0,tx);
2404             fjy0             = _mm256_add_ps(fjy0,ty);
2405             fjz0             = _mm256_add_ps(fjz0,tz);
2406
2407             /**************************
2408              * CALCULATE INTERACTIONS *
2409              **************************/
2410
2411             r11              = _mm256_mul_ps(rsq11,rinv11);
2412             r11              = _mm256_andnot_ps(dummy_mask,r11);
2413
2414             /* Calculate table index by multiplying r with table scale and truncate to integer */
2415             rt               = _mm256_mul_ps(r11,vftabscale);
2416             vfitab           = _mm256_cvttps_epi32(rt);
2417             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2418             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2419             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2420             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2421             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2422             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2423
2424             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2425             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2426                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2427             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2428                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2429             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2430                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2431             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2432                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2433             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2434             Heps             = _mm256_mul_ps(vfeps,H);
2435             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2436             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2437             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq11,FF),_mm256_mul_ps(vftabscale,rinv11)));
2438
2439             fscal            = felec;
2440
2441             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2442
2443             /* Calculate temporary vectorial force */
2444             tx               = _mm256_mul_ps(fscal,dx11);
2445             ty               = _mm256_mul_ps(fscal,dy11);
2446             tz               = _mm256_mul_ps(fscal,dz11);
2447
2448             /* Update vectorial force */
2449             fix1             = _mm256_add_ps(fix1,tx);
2450             fiy1             = _mm256_add_ps(fiy1,ty);
2451             fiz1             = _mm256_add_ps(fiz1,tz);
2452
2453             fjx1             = _mm256_add_ps(fjx1,tx);
2454             fjy1             = _mm256_add_ps(fjy1,ty);
2455             fjz1             = _mm256_add_ps(fjz1,tz);
2456
2457             /**************************
2458              * CALCULATE INTERACTIONS *
2459              **************************/
2460
2461             r12              = _mm256_mul_ps(rsq12,rinv12);
2462             r12              = _mm256_andnot_ps(dummy_mask,r12);
2463
2464             /* Calculate table index by multiplying r with table scale and truncate to integer */
2465             rt               = _mm256_mul_ps(r12,vftabscale);
2466             vfitab           = _mm256_cvttps_epi32(rt);
2467             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2468             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2469             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2470             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2471             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2472             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2473
2474             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2475             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2476                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2477             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2478                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2479             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2480                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2481             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2482                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2483             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2484             Heps             = _mm256_mul_ps(vfeps,H);
2485             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2486             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2487             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq12,FF),_mm256_mul_ps(vftabscale,rinv12)));
2488
2489             fscal            = felec;
2490
2491             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2492
2493             /* Calculate temporary vectorial force */
2494             tx               = _mm256_mul_ps(fscal,dx12);
2495             ty               = _mm256_mul_ps(fscal,dy12);
2496             tz               = _mm256_mul_ps(fscal,dz12);
2497
2498             /* Update vectorial force */
2499             fix1             = _mm256_add_ps(fix1,tx);
2500             fiy1             = _mm256_add_ps(fiy1,ty);
2501             fiz1             = _mm256_add_ps(fiz1,tz);
2502
2503             fjx2             = _mm256_add_ps(fjx2,tx);
2504             fjy2             = _mm256_add_ps(fjy2,ty);
2505             fjz2             = _mm256_add_ps(fjz2,tz);
2506
2507             /**************************
2508              * CALCULATE INTERACTIONS *
2509              **************************/
2510
2511             r13              = _mm256_mul_ps(rsq13,rinv13);
2512             r13              = _mm256_andnot_ps(dummy_mask,r13);
2513
2514             /* Calculate table index by multiplying r with table scale and truncate to integer */
2515             rt               = _mm256_mul_ps(r13,vftabscale);
2516             vfitab           = _mm256_cvttps_epi32(rt);
2517             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2518             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2519             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2520             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2521             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2522             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2523
2524             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2525             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2526                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2527             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2528                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2529             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2530                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2531             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2532                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2533             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2534             Heps             = _mm256_mul_ps(vfeps,H);
2535             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2536             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2537             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq13,FF),_mm256_mul_ps(vftabscale,rinv13)));
2538
2539             fscal            = felec;
2540
2541             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2542
2543             /* Calculate temporary vectorial force */
2544             tx               = _mm256_mul_ps(fscal,dx13);
2545             ty               = _mm256_mul_ps(fscal,dy13);
2546             tz               = _mm256_mul_ps(fscal,dz13);
2547
2548             /* Update vectorial force */
2549             fix1             = _mm256_add_ps(fix1,tx);
2550             fiy1             = _mm256_add_ps(fiy1,ty);
2551             fiz1             = _mm256_add_ps(fiz1,tz);
2552
2553             fjx3             = _mm256_add_ps(fjx3,tx);
2554             fjy3             = _mm256_add_ps(fjy3,ty);
2555             fjz3             = _mm256_add_ps(fjz3,tz);
2556
2557             /**************************
2558              * CALCULATE INTERACTIONS *
2559              **************************/
2560
2561             r21              = _mm256_mul_ps(rsq21,rinv21);
2562             r21              = _mm256_andnot_ps(dummy_mask,r21);
2563
2564             /* Calculate table index by multiplying r with table scale and truncate to integer */
2565             rt               = _mm256_mul_ps(r21,vftabscale);
2566             vfitab           = _mm256_cvttps_epi32(rt);
2567             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2568             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2569             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2570             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2571             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2572             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2573
2574             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2575             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2576                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2577             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2578                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2579             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2580                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2581             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2582                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2583             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2584             Heps             = _mm256_mul_ps(vfeps,H);
2585             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2586             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2587             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq21,FF),_mm256_mul_ps(vftabscale,rinv21)));
2588
2589             fscal            = felec;
2590
2591             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2592
2593             /* Calculate temporary vectorial force */
2594             tx               = _mm256_mul_ps(fscal,dx21);
2595             ty               = _mm256_mul_ps(fscal,dy21);
2596             tz               = _mm256_mul_ps(fscal,dz21);
2597
2598             /* Update vectorial force */
2599             fix2             = _mm256_add_ps(fix2,tx);
2600             fiy2             = _mm256_add_ps(fiy2,ty);
2601             fiz2             = _mm256_add_ps(fiz2,tz);
2602
2603             fjx1             = _mm256_add_ps(fjx1,tx);
2604             fjy1             = _mm256_add_ps(fjy1,ty);
2605             fjz1             = _mm256_add_ps(fjz1,tz);
2606
2607             /**************************
2608              * CALCULATE INTERACTIONS *
2609              **************************/
2610
2611             r22              = _mm256_mul_ps(rsq22,rinv22);
2612             r22              = _mm256_andnot_ps(dummy_mask,r22);
2613
2614             /* Calculate table index by multiplying r with table scale and truncate to integer */
2615             rt               = _mm256_mul_ps(r22,vftabscale);
2616             vfitab           = _mm256_cvttps_epi32(rt);
2617             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2618             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2619             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2620             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2621             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2622             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2623
2624             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2625             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2626                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2627             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2628                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2629             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2630                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2631             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2632                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2633             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2634             Heps             = _mm256_mul_ps(vfeps,H);
2635             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2636             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2637             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq22,FF),_mm256_mul_ps(vftabscale,rinv22)));
2638
2639             fscal            = felec;
2640
2641             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2642
2643             /* Calculate temporary vectorial force */
2644             tx               = _mm256_mul_ps(fscal,dx22);
2645             ty               = _mm256_mul_ps(fscal,dy22);
2646             tz               = _mm256_mul_ps(fscal,dz22);
2647
2648             /* Update vectorial force */
2649             fix2             = _mm256_add_ps(fix2,tx);
2650             fiy2             = _mm256_add_ps(fiy2,ty);
2651             fiz2             = _mm256_add_ps(fiz2,tz);
2652
2653             fjx2             = _mm256_add_ps(fjx2,tx);
2654             fjy2             = _mm256_add_ps(fjy2,ty);
2655             fjz2             = _mm256_add_ps(fjz2,tz);
2656
2657             /**************************
2658              * CALCULATE INTERACTIONS *
2659              **************************/
2660
2661             r23              = _mm256_mul_ps(rsq23,rinv23);
2662             r23              = _mm256_andnot_ps(dummy_mask,r23);
2663
2664             /* Calculate table index by multiplying r with table scale and truncate to integer */
2665             rt               = _mm256_mul_ps(r23,vftabscale);
2666             vfitab           = _mm256_cvttps_epi32(rt);
2667             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2668             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2669             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2670             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2671             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2672             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2673
2674             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2675             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2676                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2677             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2678                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2679             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2680                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2681             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2682                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2683             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2684             Heps             = _mm256_mul_ps(vfeps,H);
2685             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2686             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2687             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq23,FF),_mm256_mul_ps(vftabscale,rinv23)));
2688
2689             fscal            = felec;
2690
2691             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2692
2693             /* Calculate temporary vectorial force */
2694             tx               = _mm256_mul_ps(fscal,dx23);
2695             ty               = _mm256_mul_ps(fscal,dy23);
2696             tz               = _mm256_mul_ps(fscal,dz23);
2697
2698             /* Update vectorial force */
2699             fix2             = _mm256_add_ps(fix2,tx);
2700             fiy2             = _mm256_add_ps(fiy2,ty);
2701             fiz2             = _mm256_add_ps(fiz2,tz);
2702
2703             fjx3             = _mm256_add_ps(fjx3,tx);
2704             fjy3             = _mm256_add_ps(fjy3,ty);
2705             fjz3             = _mm256_add_ps(fjz3,tz);
2706
2707             /**************************
2708              * CALCULATE INTERACTIONS *
2709              **************************/
2710
2711             r31              = _mm256_mul_ps(rsq31,rinv31);
2712             r31              = _mm256_andnot_ps(dummy_mask,r31);
2713
2714             /* Calculate table index by multiplying r with table scale and truncate to integer */
2715             rt               = _mm256_mul_ps(r31,vftabscale);
2716             vfitab           = _mm256_cvttps_epi32(rt);
2717             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2718             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2719             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2720             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2721             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2722             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2723
2724             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2725             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2726                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2727             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2728                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2729             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2730                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2731             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2732                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2733             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2734             Heps             = _mm256_mul_ps(vfeps,H);
2735             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2736             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2737             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq31,FF),_mm256_mul_ps(vftabscale,rinv31)));
2738
2739             fscal            = felec;
2740
2741             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2742
2743             /* Calculate temporary vectorial force */
2744             tx               = _mm256_mul_ps(fscal,dx31);
2745             ty               = _mm256_mul_ps(fscal,dy31);
2746             tz               = _mm256_mul_ps(fscal,dz31);
2747
2748             /* Update vectorial force */
2749             fix3             = _mm256_add_ps(fix3,tx);
2750             fiy3             = _mm256_add_ps(fiy3,ty);
2751             fiz3             = _mm256_add_ps(fiz3,tz);
2752
2753             fjx1             = _mm256_add_ps(fjx1,tx);
2754             fjy1             = _mm256_add_ps(fjy1,ty);
2755             fjz1             = _mm256_add_ps(fjz1,tz);
2756
2757             /**************************
2758              * CALCULATE INTERACTIONS *
2759              **************************/
2760
2761             r32              = _mm256_mul_ps(rsq32,rinv32);
2762             r32              = _mm256_andnot_ps(dummy_mask,r32);
2763
2764             /* Calculate table index by multiplying r with table scale and truncate to integer */
2765             rt               = _mm256_mul_ps(r32,vftabscale);
2766             vfitab           = _mm256_cvttps_epi32(rt);
2767             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2768             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2769             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2770             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2771             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2772             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2773
2774             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2775             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2776                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2777             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2778                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2779             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2780                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2781             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2782                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2783             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2784             Heps             = _mm256_mul_ps(vfeps,H);
2785             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2786             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2787             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq32,FF),_mm256_mul_ps(vftabscale,rinv32)));
2788
2789             fscal            = felec;
2790
2791             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2792
2793             /* Calculate temporary vectorial force */
2794             tx               = _mm256_mul_ps(fscal,dx32);
2795             ty               = _mm256_mul_ps(fscal,dy32);
2796             tz               = _mm256_mul_ps(fscal,dz32);
2797
2798             /* Update vectorial force */
2799             fix3             = _mm256_add_ps(fix3,tx);
2800             fiy3             = _mm256_add_ps(fiy3,ty);
2801             fiz3             = _mm256_add_ps(fiz3,tz);
2802
2803             fjx2             = _mm256_add_ps(fjx2,tx);
2804             fjy2             = _mm256_add_ps(fjy2,ty);
2805             fjz2             = _mm256_add_ps(fjz2,tz);
2806
2807             /**************************
2808              * CALCULATE INTERACTIONS *
2809              **************************/
2810
2811             r33              = _mm256_mul_ps(rsq33,rinv33);
2812             r33              = _mm256_andnot_ps(dummy_mask,r33);
2813
2814             /* Calculate table index by multiplying r with table scale and truncate to integer */
2815             rt               = _mm256_mul_ps(r33,vftabscale);
2816             vfitab           = _mm256_cvttps_epi32(rt);
2817             vfeps            = _mm256_sub_ps(rt,_mm256_round_ps(rt, _MM_FROUND_FLOOR));
2818             /*         AVX1 does not support 256-bit integer operations, so now we go to 128-bit mode... */
2819             vfitab_lo        = _mm256_extractf128_si256(vfitab,0x0);
2820             vfitab_hi        = _mm256_extractf128_si256(vfitab,0x1);
2821             vfitab_lo        = _mm_slli_epi32(vfitab_lo,2);
2822             vfitab_hi        = _mm_slli_epi32(vfitab_hi,2);
2823
2824             /* CUBIC SPLINE TABLE ELECTROSTATICS */
2825             Y                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,0)),
2826                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,0)));
2827             F                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,1)),
2828                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,1)));
2829             G                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,2)),
2830                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,2)));
2831             H                = gmx_mm256_set_m128(_mm_load_ps(vftab + _mm_extract_epi32(vfitab_hi,3)),
2832                                                   _mm_load_ps(vftab + _mm_extract_epi32(vfitab_lo,3)));
2833             GMX_MM256_HALFTRANSPOSE4_PS(Y,F,G,H);
2834             Heps             = _mm256_mul_ps(vfeps,H);
2835             Fp               = _mm256_add_ps(F,_mm256_mul_ps(vfeps,_mm256_add_ps(G,Heps)));
2836             FF               = _mm256_add_ps(Fp,_mm256_mul_ps(vfeps,_mm256_add_ps(G,_mm256_add_ps(Heps,Heps))));
2837             felec            = _mm256_xor_ps(signbit,_mm256_mul_ps(_mm256_mul_ps(qq33,FF),_mm256_mul_ps(vftabscale,rinv33)));
2838
2839             fscal            = felec;
2840
2841             fscal            = _mm256_andnot_ps(dummy_mask,fscal);
2842
2843             /* Calculate temporary vectorial force */
2844             tx               = _mm256_mul_ps(fscal,dx33);
2845             ty               = _mm256_mul_ps(fscal,dy33);
2846             tz               = _mm256_mul_ps(fscal,dz33);
2847
2848             /* Update vectorial force */
2849             fix3             = _mm256_add_ps(fix3,tx);
2850             fiy3             = _mm256_add_ps(fiy3,ty);
2851             fiz3             = _mm256_add_ps(fiz3,tz);
2852
2853             fjx3             = _mm256_add_ps(fjx3,tx);
2854             fjy3             = _mm256_add_ps(fjy3,ty);
2855             fjz3             = _mm256_add_ps(fjz3,tz);
2856
2857             fjptrA             = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
2858             fjptrB             = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
2859             fjptrC             = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
2860             fjptrD             = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
2861             fjptrE             = (jnrlistE>=0) ? f+j_coord_offsetE : scratch;
2862             fjptrF             = (jnrlistF>=0) ? f+j_coord_offsetF : scratch;
2863             fjptrG             = (jnrlistG>=0) ? f+j_coord_offsetG : scratch;
2864             fjptrH             = (jnrlistH>=0) ? f+j_coord_offsetH : scratch;
2865
2866             gmx_mm256_decrement_4rvec_8ptr_swizzle_ps(fjptrA,fjptrB,fjptrC,fjptrD,fjptrE,fjptrF,fjptrG,fjptrH,
2867                                                       fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,
2868                                                       fjx2,fjy2,fjz2,fjx3,fjy3,fjz3);
2869
2870             /* Inner loop uses 390 flops */
2871         }
2872
2873         /* End of innermost loop */
2874
2875         gmx_mm256_update_iforce_4atom_swizzle_ps(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
2876                                                  f+i_coord_offset,fshift+i_shift_offset);
2877
2878         /* Increment number of inner iterations */
2879         inneriter                  += j_index_end - j_index_start;
2880
2881         /* Outer loop uses 24 flops */
2882     }
2883
2884     /* Increment number of outer iterations */
2885     outeriter        += nri;
2886
2887     /* Update outer/inner flops */
2888
2889     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4W4_F,outeriter*24 + inneriter*390);
2890 }