src/gromacs/gmxlib/nonbonded/nb_kernel_avx_256_double/nb_kernel_ElecCSTab_VdwCSTab_GeomW4W4_avx_256_double.c

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