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

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