src/gmxlib/nonbonded/nb_kernel_avx_256_double/nb_kernel_ElecCSTab_VdwCSTab_GeomW4P1_avx_256_double.c

   1 /*
   2  * Note: this file was generated by the Gromacs avx_256_double kernel generator.
   3  *
   4  *                This source code is part of
   5  *
   6  *                 G   R   O   M   A   C   S
   7  *
   8  * Copyright (c) 2001-2012, The GROMACS Development Team
   9  *
  10  * Gromacs is a library for molecular simulation and trajectory analysis,
  11  * written by Erik Lindahl, David van der Spoel, Berk Hess, and others - for
  12  * a full list of developers and information, check out http://www.gromacs.org
  13  *
  14  * This program is free software; you can redistribute it and/or modify it under
  15  * the terms of the GNU Lesser General Public License as published by the Free
  16  * Software Foundation; either version 2 of the License, or (at your option) any
  17  * later version.
  18  *
  19  * To help fund GROMACS development, we humbly ask that you cite
  20  * the papers people have written on it - you can find them on the website.
  21  */
  22 #ifdef HAVE_CONFIG_H
  23 #include <config.h>
  24 #endif
  25
  26 #include <math.h>
  27
  28 #include "../nb_kernel.h"
  29 #include "types/simple.h"
  30 #include "vec.h"
  31 #include "nrnb.h"
  32
  33 #include "gmx_math_x86_avx_256_double.h"
  34 #include "kernelutil_x86_avx_256_double.h"
  35
  36 /*
  37  * Gromacs nonbonded kernel:   nb_kernel_ElecCSTab_VdwCSTab_GeomW4P1_VF_avx_256_double
  38  * Electrostatics interaction: CubicSplineTable
  39  * VdW interaction:            CubicSplineTable
  40  * Geometry:                   Water4-Particle
  41  * Calculate force/pot:        PotentialAndForce
  42  */
  43 void
  44 nb_kernel_ElecCSTab_VdwCSTab_GeomW4P1_VF_avx_256_double
  45                     (t_nblist * gmx_restrict                nlist,
  46                      rvec * gmx_restrict                    xx,
  47                      rvec * gmx_restrict                    ff,
  48                      t_forcerec * gmx_restrict              fr,
  49                      t_mdatoms * gmx_restrict               mdatoms,
  50                      nb_kernel_data_t * gmx_restrict        kernel_data,
  51                      t_nrnb * gmx_restrict                  nrnb)
  52 {
  53     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
  54      * just 0 for non-waters.
  55      * Suffixes A,B,C,D refer to j loop unrolling done with AVX, e.g. for the four different
  56      * jnr indices corresponding to data put in the four positions in the SIMD register.
  57      */
  58     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
  59     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
  60     int              jnrA,jnrB,jnrC,jnrD;
  61     int              jnrlistA,jnrlistB,jnrlistC,jnrlistD;
  62     int              jnrlistE,jnrlistF,jnrlistG,jnrlistH;
  63     int              j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
  64     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
  65     real             rcutoff_scalar;
  66     real             *shiftvec,*fshift,*x,*f;
  67     real             *fjptrA,*fjptrB,*fjptrC,*fjptrD;
  68     real             scratch[4*DIM];
  69     __m256d          tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
  70     real *           vdwioffsetptr0;
  71     __m256d          ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
  72     real *           vdwioffsetptr1;
  73     __m256d          ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
  74     real *           vdwioffsetptr2;
  75     __m256d          ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
  76     real *           vdwioffsetptr3;
  77     __m256d          ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
  78     int              vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
  79     __m256d          jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
  80     __m256d          dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
  81     __m256d          dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
  82     __m256d          dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
  83     __m256d          dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
  84     __m256d          velec,felec,velecsum,facel,crf,krf,krf2;
  85     real             *charge;
  86     int              nvdwtype;
  87     __m256d          rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
  88     int              *vdwtype;
  89     real             *vdwparam;
  90     __m256d          one_sixth   = _mm256_set1_pd(1.0/6.0);
  91     __m256d          one_twelfth = _mm256_set1_pd(1.0/12.0);
  92     __m128i          vfitab;
  93     __m128i          ifour       = _mm_set1_epi32(4);
  94     __m256d          rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF;
  95     real             *vftab;
  96     __m256d          dummy_mask,cutoff_mask;
  97     __m128           tmpmask0,tmpmask1;
  98     __m256d          signbit = _mm256_castsi256_pd( _mm256_set1_epi32(0x80000000) );
  99     __m256d          one     = _mm256_set1_pd(1.0);
 100     __m256d          two     = _mm256_set1_pd(2.0);
 101     x                = xx[0];
 102     f                = ff[0];
 103
 104     nri              = nlist->nri;
 105     iinr             = nlist->iinr;
 106     jindex           = nlist->jindex;
 107     jjnr             = nlist->jjnr;
 108     shiftidx         = nlist->shift;
 109     gid              = nlist->gid;
 110     shiftvec         = fr->shift_vec[0];
 111     fshift           = fr->fshift[0];
 112     facel            = _mm256_set1_pd(fr->epsfac);
 113     charge           = mdatoms->chargeA;
 114     nvdwtype         = fr->ntype;
 115     vdwparam         = fr->nbfp;
 116     vdwtype          = mdatoms->typeA;
 117
 118     vftab            = kernel_data->table_elec_vdw->data;
 119     vftabscale       = _mm256_set1_pd(kernel_data->table_elec_vdw->scale);
 120
 121     /* Setup water-specific parameters */
 122     inr              = nlist->iinr[0];
 123     iq1              = _mm256_mul_pd(facel,_mm256_set1_pd(charge[inr+1]));
 124     iq2              = _mm256_mul_pd(facel,_mm256_set1_pd(charge[inr+2]));
 125     iq3              = _mm256_mul_pd(facel,_mm256_set1_pd(charge[inr+3]));
 126     vdwioffsetptr0   = vdwparam+2*nvdwtype*vdwtype[inr+0];
 127
 128     /* Avoid stupid compiler warnings */
 129     jnrA = jnrB = jnrC = jnrD = 0;
 130     j_coord_offsetA = 0;
 131     j_coord_offsetB = 0;
 132     j_coord_offsetC = 0;
 133     j_coord_offsetD = 0;
 134
 135     outeriter        = 0;
 136     inneriter        = 0;
 137
 138     for(iidx=0;iidx<4*DIM;iidx++)
 139     {
 140         scratch[iidx] = 0.0;
 141     }
 142
 143     /* Start outer loop over neighborlists */
 144     for(iidx=0; iidx<nri; iidx++)
 145     {
 146         /* Load shift vector for this list */
 147         i_shift_offset   = DIM*shiftidx[iidx];
 148
 149         /* Load limits for loop over neighbors */
 150         j_index_start    = jindex[iidx];
 151         j_index_end      = jindex[iidx+1];
 152
 153         /* Get outer coordinate index */
 154         inr              = iinr[iidx];
 155         i_coord_offset   = DIM*inr;
 156
 157         /* Load i particle coords and add shift vector */
 158         gmx_mm256_load_shift_and_4rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
 159                                                     &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
 160
 161         fix0             = _mm256_setzero_pd();
 162         fiy0             = _mm256_setzero_pd();
 163         fiz0             = _mm256_setzero_pd();
 164         fix1             = _mm256_setzero_pd();
 165         fiy1             = _mm256_setzero_pd();
 166         fiz1             = _mm256_setzero_pd();
 167         fix2             = _mm256_setzero_pd();
 168         fiy2             = _mm256_setzero_pd();
 169         fiz2             = _mm256_setzero_pd();
 170         fix3             = _mm256_setzero_pd();
 171         fiy3             = _mm256_setzero_pd();
 172         fiz3             = _mm256_setzero_pd();
 173
 174         /* Reset potential sums */
 175         velecsum         = _mm256_setzero_pd();
 176         vvdwsum          = _mm256_setzero_pd();
 177
 178         /* Start inner kernel loop */
 179         for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
 180         {
 181
 182             /* Get j neighbor index, and coordinate index */
 183             jnrA             = jjnr[jidx];
 184             jnrB             = jjnr[jidx+1];
 185             jnrC             = jjnr[jidx+2];
 186             jnrD             = jjnr[jidx+3];
 187             j_coord_offsetA  = DIM*jnrA;
 188             j_coord_offsetB  = DIM*jnrB;
 189             j_coord_offsetC  = DIM*jnrC;
 190             j_coord_offsetD  = DIM*jnrD;
 191
 192             /* load j atom coordinates */
 193             gmx_mm256_load_1rvec_4ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
 194                                                  x+j_coord_offsetC,x+j_coord_offsetD,
 195                                                  &jx0,&jy0,&jz0);
 196
 197             /* Calculate displacement vector */
 198             dx00             = _mm256_sub_pd(ix0,jx0);
 199             dy00             = _mm256_sub_pd(iy0,jy0);
 200             dz00             = _mm256_sub_pd(iz0,jz0);
 201             dx10             = _mm256_sub_pd(ix1,jx0);
 202             dy10             = _mm256_sub_pd(iy1,jy0);
 203             dz10             = _mm256_sub_pd(iz1,jz0);
 204             dx20             = _mm256_sub_pd(ix2,jx0);
 205             dy20             = _mm256_sub_pd(iy2,jy0);
 206             dz20             = _mm256_sub_pd(iz2,jz0);
 207             dx30             = _mm256_sub_pd(ix3,jx0);
 208             dy30             = _mm256_sub_pd(iy3,jy0);
 209             dz30             = _mm256_sub_pd(iz3,jz0);
 210
 211             /* Calculate squared distance and things based on it */
 212             rsq00            = gmx_mm256_calc_rsq_pd(dx00,dy00,dz00);
 213             rsq10            = gmx_mm256_calc_rsq_pd(dx10,dy10,dz10);
 214             rsq20            = gmx_mm256_calc_rsq_pd(dx20,dy20,dz20);
 215             rsq30            = gmx_mm256_calc_rsq_pd(dx30,dy30,dz30);
 216
 217             rinv00           = gmx_mm256_invsqrt_pd(rsq00);
 218             rinv10           = gmx_mm256_invsqrt_pd(rsq10);
 219             rinv20           = gmx_mm256_invsqrt_pd(rsq20);
 220             rinv30           = gmx_mm256_invsqrt_pd(rsq30);
 221
 222             /* Load parameters for j particles */
 223             jq0              = gmx_mm256_load_4real_swizzle_pd(charge+jnrA+0,charge+jnrB+0,
 224                                                                  charge+jnrC+0,charge+jnrD+0);
 225             vdwjidx0A        = 2*vdwtype[jnrA+0];
 226             vdwjidx0B        = 2*vdwtype[jnrB+0];
 227             vdwjidx0C        = 2*vdwtype[jnrC+0];
 228             vdwjidx0D        = 2*vdwtype[jnrD+0];
 229
 230             fjx0             = _mm256_setzero_pd();
 231             fjy0             = _mm256_setzero_pd();
 232             fjz0             = _mm256_setzero_pd();
 233
 234             /**************************
 235              * CALCULATE INTERACTIONS *
 236              **************************/
 237
 238             r00              = _mm256_mul_pd(rsq00,rinv00);
 239
 240             /* Compute parameters for interactions between i and j atoms */
 241             gmx_mm256_load_4pair_swizzle_pd(vdwioffsetptr0+vdwjidx0A,
 242                                             vdwioffsetptr0+vdwjidx0B,
 243                                             vdwioffsetptr0+vdwjidx0C,
 244                                             vdwioffsetptr0+vdwjidx0D,
 245                                             &c6_00,&c12_00);
 246
 247             /* Calculate table index by multiplying r with table scale and truncate to integer */
 248             rt               = _mm256_mul_pd(r00,vftabscale);
 249             vfitab           = _mm256_cvttpd_epi32(rt);
 250             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 251             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 252
 253             /* CUBIC SPLINE TABLE DISPERSION */
 254             vfitab           = _mm_add_epi32(vfitab,ifour);
 255             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 256             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 257             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 258             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 259             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 260             Heps             = _mm256_mul_pd(vfeps,H);
 261             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 262             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 263             vvdw6            = _mm256_mul_pd(c6_00,VV);
 264             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 265             fvdw6            = _mm256_mul_pd(c6_00,FF);
 266
 267             /* CUBIC SPLINE TABLE REPULSION */
 268             vfitab           = _mm_add_epi32(vfitab,ifour);
 269             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 270             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 271             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 272             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 273             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 274             Heps             = _mm256_mul_pd(vfeps,H);
 275             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 276             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 277             vvdw12           = _mm256_mul_pd(c12_00,VV);
 278             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 279             fvdw12           = _mm256_mul_pd(c12_00,FF);
 280             vvdw             = _mm256_add_pd(vvdw12,vvdw6);
 281             fvdw             = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_add_pd(fvdw6,fvdw12),_mm256_mul_pd(vftabscale,rinv00)));
 282
 283             /* Update potential sum for this i atom from the interaction with this j atom. */
 284             vvdwsum          = _mm256_add_pd(vvdwsum,vvdw);
 285
 286             fscal            = fvdw;
 287
 288             /* Calculate temporary vectorial force */
 289             tx               = _mm256_mul_pd(fscal,dx00);
 290             ty               = _mm256_mul_pd(fscal,dy00);
 291             tz               = _mm256_mul_pd(fscal,dz00);
 292
 293             /* Update vectorial force */
 294             fix0             = _mm256_add_pd(fix0,tx);
 295             fiy0             = _mm256_add_pd(fiy0,ty);
 296             fiz0             = _mm256_add_pd(fiz0,tz);
 297
 298             fjx0             = _mm256_add_pd(fjx0,tx);
 299             fjy0             = _mm256_add_pd(fjy0,ty);
 300             fjz0             = _mm256_add_pd(fjz0,tz);
 301
 302             /**************************
 303              * CALCULATE INTERACTIONS *
 304              **************************/
 305
 306             r10              = _mm256_mul_pd(rsq10,rinv10);
 307
 308             /* Compute parameters for interactions between i and j atoms */
 309             qq10             = _mm256_mul_pd(iq1,jq0);
 310
 311             /* Calculate table index by multiplying r with table scale and truncate to integer */
 312             rt               = _mm256_mul_pd(r10,vftabscale);
 313             vfitab           = _mm256_cvttpd_epi32(rt);
 314             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 315             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 316
 317             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 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             velec            = _mm256_mul_pd(qq10,VV);
 327             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 328             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq10,FF),_mm256_mul_pd(vftabscale,rinv10)));
 329
 330             /* Update potential sum for this i atom from the interaction with this j atom. */
 331             velecsum         = _mm256_add_pd(velecsum,velec);
 332
 333             fscal            = felec;
 334
 335             /* Calculate temporary vectorial force */
 336             tx               = _mm256_mul_pd(fscal,dx10);
 337             ty               = _mm256_mul_pd(fscal,dy10);
 338             tz               = _mm256_mul_pd(fscal,dz10);
 339
 340             /* Update vectorial force */
 341             fix1             = _mm256_add_pd(fix1,tx);
 342             fiy1             = _mm256_add_pd(fiy1,ty);
 343             fiz1             = _mm256_add_pd(fiz1,tz);
 344
 345             fjx0             = _mm256_add_pd(fjx0,tx);
 346             fjy0             = _mm256_add_pd(fjy0,ty);
 347             fjz0             = _mm256_add_pd(fjz0,tz);
 348
 349             /**************************
 350              * CALCULATE INTERACTIONS *
 351              **************************/
 352
 353             r20              = _mm256_mul_pd(rsq20,rinv20);
 354
 355             /* Compute parameters for interactions between i and j atoms */
 356             qq20             = _mm256_mul_pd(iq2,jq0);
 357
 358             /* Calculate table index by multiplying r with table scale and truncate to integer */
 359             rt               = _mm256_mul_pd(r20,vftabscale);
 360             vfitab           = _mm256_cvttpd_epi32(rt);
 361             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 362             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 363
 364             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 365             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 366             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 367             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 368             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 369             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 370             Heps             = _mm256_mul_pd(vfeps,H);
 371             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 372             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 373             velec            = _mm256_mul_pd(qq20,VV);
 374             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 375             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq20,FF),_mm256_mul_pd(vftabscale,rinv20)));
 376
 377             /* Update potential sum for this i atom from the interaction with this j atom. */
 378             velecsum         = _mm256_add_pd(velecsum,velec);
 379
 380             fscal            = felec;
 381
 382             /* Calculate temporary vectorial force */
 383             tx               = _mm256_mul_pd(fscal,dx20);
 384             ty               = _mm256_mul_pd(fscal,dy20);
 385             tz               = _mm256_mul_pd(fscal,dz20);
 386
 387             /* Update vectorial force */
 388             fix2             = _mm256_add_pd(fix2,tx);
 389             fiy2             = _mm256_add_pd(fiy2,ty);
 390             fiz2             = _mm256_add_pd(fiz2,tz);
 391
 392             fjx0             = _mm256_add_pd(fjx0,tx);
 393             fjy0             = _mm256_add_pd(fjy0,ty);
 394             fjz0             = _mm256_add_pd(fjz0,tz);
 395
 396             /**************************
 397              * CALCULATE INTERACTIONS *
 398              **************************/
 399
 400             r30              = _mm256_mul_pd(rsq30,rinv30);
 401
 402             /* Compute parameters for interactions between i and j atoms */
 403             qq30             = _mm256_mul_pd(iq3,jq0);
 404
 405             /* Calculate table index by multiplying r with table scale and truncate to integer */
 406             rt               = _mm256_mul_pd(r30,vftabscale);
 407             vfitab           = _mm256_cvttpd_epi32(rt);
 408             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 409             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 410
 411             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 412             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 413             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 414             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 415             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 416             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 417             Heps             = _mm256_mul_pd(vfeps,H);
 418             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 419             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 420             velec            = _mm256_mul_pd(qq30,VV);
 421             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 422             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq30,FF),_mm256_mul_pd(vftabscale,rinv30)));
 423
 424             /* Update potential sum for this i atom from the interaction with this j atom. */
 425             velecsum         = _mm256_add_pd(velecsum,velec);
 426
 427             fscal            = felec;
 428
 429             /* Calculate temporary vectorial force */
 430             tx               = _mm256_mul_pd(fscal,dx30);
 431             ty               = _mm256_mul_pd(fscal,dy30);
 432             tz               = _mm256_mul_pd(fscal,dz30);
 433
 434             /* Update vectorial force */
 435             fix3             = _mm256_add_pd(fix3,tx);
 436             fiy3             = _mm256_add_pd(fiy3,ty);
 437             fiz3             = _mm256_add_pd(fiz3,tz);
 438
 439             fjx0             = _mm256_add_pd(fjx0,tx);
 440             fjy0             = _mm256_add_pd(fjy0,ty);
 441             fjz0             = _mm256_add_pd(fjz0,tz);
 442
 443             fjptrA             = f+j_coord_offsetA;
 444             fjptrB             = f+j_coord_offsetB;
 445             fjptrC             = f+j_coord_offsetC;
 446             fjptrD             = f+j_coord_offsetD;
 447
 448             gmx_mm256_decrement_1rvec_4ptr_swizzle_pd(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
 449
 450             /* Inner loop uses 188 flops */
 451         }
 452
 453         if(jidx<j_index_end)
 454         {
 455
 456             /* Get j neighbor index, and coordinate index */
 457             jnrlistA         = jjnr[jidx];
 458             jnrlistB         = jjnr[jidx+1];
 459             jnrlistC         = jjnr[jidx+2];
 460             jnrlistD         = jjnr[jidx+3];
 461             /* Sign of each element will be negative for non-real atoms.
 462              * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
 463              * so use it as val = _mm_andnot_pd(mask,val) to clear dummy entries.
 464              */
 465             tmpmask0 = gmx_mm_castsi128_pd(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
 466
 467             tmpmask1 = _mm_permute_ps(tmpmask0,_GMX_MM_PERMUTE(3,3,2,2));
 468             tmpmask0 = _mm_permute_ps(tmpmask0,_GMX_MM_PERMUTE(1,1,0,0));
 469             dummy_mask = _mm256_castps_pd(gmx_mm256_set_m128(tmpmask1,tmpmask0));
 470
 471             jnrA       = (jnrlistA>=0) ? jnrlistA : 0;
 472             jnrB       = (jnrlistB>=0) ? jnrlistB : 0;
 473             jnrC       = (jnrlistC>=0) ? jnrlistC : 0;
 474             jnrD       = (jnrlistD>=0) ? jnrlistD : 0;
 475             j_coord_offsetA  = DIM*jnrA;
 476             j_coord_offsetB  = DIM*jnrB;
 477             j_coord_offsetC  = DIM*jnrC;
 478             j_coord_offsetD  = DIM*jnrD;
 479
 480             /* load j atom coordinates */
 481             gmx_mm256_load_1rvec_4ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
 482                                                  x+j_coord_offsetC,x+j_coord_offsetD,
 483                                                  &jx0,&jy0,&jz0);
 484
 485             /* Calculate displacement vector */
 486             dx00             = _mm256_sub_pd(ix0,jx0);
 487             dy00             = _mm256_sub_pd(iy0,jy0);
 488             dz00             = _mm256_sub_pd(iz0,jz0);
 489             dx10             = _mm256_sub_pd(ix1,jx0);
 490             dy10             = _mm256_sub_pd(iy1,jy0);
 491             dz10             = _mm256_sub_pd(iz1,jz0);
 492             dx20             = _mm256_sub_pd(ix2,jx0);
 493             dy20             = _mm256_sub_pd(iy2,jy0);
 494             dz20             = _mm256_sub_pd(iz2,jz0);
 495             dx30             = _mm256_sub_pd(ix3,jx0);
 496             dy30             = _mm256_sub_pd(iy3,jy0);
 497             dz30             = _mm256_sub_pd(iz3,jz0);
 498
 499             /* Calculate squared distance and things based on it */
 500             rsq00            = gmx_mm256_calc_rsq_pd(dx00,dy00,dz00);
 501             rsq10            = gmx_mm256_calc_rsq_pd(dx10,dy10,dz10);
 502             rsq20            = gmx_mm256_calc_rsq_pd(dx20,dy20,dz20);
 503             rsq30            = gmx_mm256_calc_rsq_pd(dx30,dy30,dz30);
 504
 505             rinv00           = gmx_mm256_invsqrt_pd(rsq00);
 506             rinv10           = gmx_mm256_invsqrt_pd(rsq10);
 507             rinv20           = gmx_mm256_invsqrt_pd(rsq20);
 508             rinv30           = gmx_mm256_invsqrt_pd(rsq30);
 509
 510             /* Load parameters for j particles */
 511             jq0              = gmx_mm256_load_4real_swizzle_pd(charge+jnrA+0,charge+jnrB+0,
 512                                                                  charge+jnrC+0,charge+jnrD+0);
 513             vdwjidx0A        = 2*vdwtype[jnrA+0];
 514             vdwjidx0B        = 2*vdwtype[jnrB+0];
 515             vdwjidx0C        = 2*vdwtype[jnrC+0];
 516             vdwjidx0D        = 2*vdwtype[jnrD+0];
 517
 518             fjx0             = _mm256_setzero_pd();
 519             fjy0             = _mm256_setzero_pd();
 520             fjz0             = _mm256_setzero_pd();
 521
 522             /**************************
 523              * CALCULATE INTERACTIONS *
 524              **************************/
 525
 526             r00              = _mm256_mul_pd(rsq00,rinv00);
 527             r00              = _mm256_andnot_pd(dummy_mask,r00);
 528
 529             /* Compute parameters for interactions between i and j atoms */
 530             gmx_mm256_load_4pair_swizzle_pd(vdwioffsetptr0+vdwjidx0A,
 531                                             vdwioffsetptr0+vdwjidx0B,
 532                                             vdwioffsetptr0+vdwjidx0C,
 533                                             vdwioffsetptr0+vdwjidx0D,
 534                                             &c6_00,&c12_00);
 535
 536             /* Calculate table index by multiplying r with table scale and truncate to integer */
 537             rt               = _mm256_mul_pd(r00,vftabscale);
 538             vfitab           = _mm256_cvttpd_epi32(rt);
 539             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 540             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 541
 542             /* CUBIC SPLINE TABLE DISPERSION */
 543             vfitab           = _mm_add_epi32(vfitab,ifour);
 544             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 545             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 546             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 547             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 548             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 549             Heps             = _mm256_mul_pd(vfeps,H);
 550             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 551             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 552             vvdw6            = _mm256_mul_pd(c6_00,VV);
 553             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 554             fvdw6            = _mm256_mul_pd(c6_00,FF);
 555
 556             /* CUBIC SPLINE TABLE REPULSION */
 557             vfitab           = _mm_add_epi32(vfitab,ifour);
 558             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 559             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 560             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 561             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 562             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 563             Heps             = _mm256_mul_pd(vfeps,H);
 564             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 565             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 566             vvdw12           = _mm256_mul_pd(c12_00,VV);
 567             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 568             fvdw12           = _mm256_mul_pd(c12_00,FF);
 569             vvdw             = _mm256_add_pd(vvdw12,vvdw6);
 570             fvdw             = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_add_pd(fvdw6,fvdw12),_mm256_mul_pd(vftabscale,rinv00)));
 571
 572             /* Update potential sum for this i atom from the interaction with this j atom. */
 573             vvdw             = _mm256_andnot_pd(dummy_mask,vvdw);
 574             vvdwsum          = _mm256_add_pd(vvdwsum,vvdw);
 575
 576             fscal            = fvdw;
 577
 578             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
 579
 580             /* Calculate temporary vectorial force */
 581             tx               = _mm256_mul_pd(fscal,dx00);
 582             ty               = _mm256_mul_pd(fscal,dy00);
 583             tz               = _mm256_mul_pd(fscal,dz00);
 584
 585             /* Update vectorial force */
 586             fix0             = _mm256_add_pd(fix0,tx);
 587             fiy0             = _mm256_add_pd(fiy0,ty);
 588             fiz0             = _mm256_add_pd(fiz0,tz);
 589
 590             fjx0             = _mm256_add_pd(fjx0,tx);
 591             fjy0             = _mm256_add_pd(fjy0,ty);
 592             fjz0             = _mm256_add_pd(fjz0,tz);
 593
 594             /**************************
 595              * CALCULATE INTERACTIONS *
 596              **************************/
 597
 598             r10              = _mm256_mul_pd(rsq10,rinv10);
 599             r10              = _mm256_andnot_pd(dummy_mask,r10);
 600
 601             /* Compute parameters for interactions between i and j atoms */
 602             qq10             = _mm256_mul_pd(iq1,jq0);
 603
 604             /* Calculate table index by multiplying r with table scale and truncate to integer */
 605             rt               = _mm256_mul_pd(r10,vftabscale);
 606             vfitab           = _mm256_cvttpd_epi32(rt);
 607             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 608             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 609
 610             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 611             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 612             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 613             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 614             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 615             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 616             Heps             = _mm256_mul_pd(vfeps,H);
 617             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 618             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 619             velec            = _mm256_mul_pd(qq10,VV);
 620             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 621             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq10,FF),_mm256_mul_pd(vftabscale,rinv10)));
 622
 623             /* Update potential sum for this i atom from the interaction with this j atom. */
 624             velec            = _mm256_andnot_pd(dummy_mask,velec);
 625             velecsum         = _mm256_add_pd(velecsum,velec);
 626
 627             fscal            = felec;
 628
 629             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
 630
 631             /* Calculate temporary vectorial force */
 632             tx               = _mm256_mul_pd(fscal,dx10);
 633             ty               = _mm256_mul_pd(fscal,dy10);
 634             tz               = _mm256_mul_pd(fscal,dz10);
 635
 636             /* Update vectorial force */
 637             fix1             = _mm256_add_pd(fix1,tx);
 638             fiy1             = _mm256_add_pd(fiy1,ty);
 639             fiz1             = _mm256_add_pd(fiz1,tz);
 640
 641             fjx0             = _mm256_add_pd(fjx0,tx);
 642             fjy0             = _mm256_add_pd(fjy0,ty);
 643             fjz0             = _mm256_add_pd(fjz0,tz);
 644
 645             /**************************
 646              * CALCULATE INTERACTIONS *
 647              **************************/
 648
 649             r20              = _mm256_mul_pd(rsq20,rinv20);
 650             r20              = _mm256_andnot_pd(dummy_mask,r20);
 651
 652             /* Compute parameters for interactions between i and j atoms */
 653             qq20             = _mm256_mul_pd(iq2,jq0);
 654
 655             /* Calculate table index by multiplying r with table scale and truncate to integer */
 656             rt               = _mm256_mul_pd(r20,vftabscale);
 657             vfitab           = _mm256_cvttpd_epi32(rt);
 658             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 659             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 660
 661             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 662             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 663             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 664             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 665             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 666             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 667             Heps             = _mm256_mul_pd(vfeps,H);
 668             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 669             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 670             velec            = _mm256_mul_pd(qq20,VV);
 671             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 672             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq20,FF),_mm256_mul_pd(vftabscale,rinv20)));
 673
 674             /* Update potential sum for this i atom from the interaction with this j atom. */
 675             velec            = _mm256_andnot_pd(dummy_mask,velec);
 676             velecsum         = _mm256_add_pd(velecsum,velec);
 677
 678             fscal            = felec;
 679
 680             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
 681
 682             /* Calculate temporary vectorial force */
 683             tx               = _mm256_mul_pd(fscal,dx20);
 684             ty               = _mm256_mul_pd(fscal,dy20);
 685             tz               = _mm256_mul_pd(fscal,dz20);
 686
 687             /* Update vectorial force */
 688             fix2             = _mm256_add_pd(fix2,tx);
 689             fiy2             = _mm256_add_pd(fiy2,ty);
 690             fiz2             = _mm256_add_pd(fiz2,tz);
 691
 692             fjx0             = _mm256_add_pd(fjx0,tx);
 693             fjy0             = _mm256_add_pd(fjy0,ty);
 694             fjz0             = _mm256_add_pd(fjz0,tz);
 695
 696             /**************************
 697              * CALCULATE INTERACTIONS *
 698              **************************/
 699
 700             r30              = _mm256_mul_pd(rsq30,rinv30);
 701             r30              = _mm256_andnot_pd(dummy_mask,r30);
 702
 703             /* Compute parameters for interactions between i and j atoms */
 704             qq30             = _mm256_mul_pd(iq3,jq0);
 705
 706             /* Calculate table index by multiplying r with table scale and truncate to integer */
 707             rt               = _mm256_mul_pd(r30,vftabscale);
 708             vfitab           = _mm256_cvttpd_epi32(rt);
 709             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 710             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 711
 712             /* CUBIC SPLINE TABLE ELECTROSTATICS */
 713             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 714             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 715             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 716             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 717             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
 718             Heps             = _mm256_mul_pd(vfeps,H);
 719             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
 720             VV               = _mm256_add_pd(Y,_mm256_mul_pd(vfeps,Fp));
 721             velec            = _mm256_mul_pd(qq30,VV);
 722             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
 723             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq30,FF),_mm256_mul_pd(vftabscale,rinv30)));
 724
 725             /* Update potential sum for this i atom from the interaction with this j atom. */
 726             velec            = _mm256_andnot_pd(dummy_mask,velec);
 727             velecsum         = _mm256_add_pd(velecsum,velec);
 728
 729             fscal            = felec;
 730
 731             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
 732
 733             /* Calculate temporary vectorial force */
 734             tx               = _mm256_mul_pd(fscal,dx30);
 735             ty               = _mm256_mul_pd(fscal,dy30);
 736             tz               = _mm256_mul_pd(fscal,dz30);
 737
 738             /* Update vectorial force */
 739             fix3             = _mm256_add_pd(fix3,tx);
 740             fiy3             = _mm256_add_pd(fiy3,ty);
 741             fiz3             = _mm256_add_pd(fiz3,tz);
 742
 743             fjx0             = _mm256_add_pd(fjx0,tx);
 744             fjy0             = _mm256_add_pd(fjy0,ty);
 745             fjz0             = _mm256_add_pd(fjz0,tz);
 746
 747             fjptrA             = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
 748             fjptrB             = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
 749             fjptrC             = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
 750             fjptrD             = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
 751
 752             gmx_mm256_decrement_1rvec_4ptr_swizzle_pd(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
 753
 754             /* Inner loop uses 192 flops */
 755         }
 756
 757         /* End of innermost loop */
 758
 759         gmx_mm256_update_iforce_4atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
 760                                                  f+i_coord_offset,fshift+i_shift_offset);
 761
 762         ggid                        = gid[iidx];
 763         /* Update potential energies */
 764         gmx_mm256_update_1pot_pd(velecsum,kernel_data->energygrp_elec+ggid);
 765         gmx_mm256_update_1pot_pd(vvdwsum,kernel_data->energygrp_vdw+ggid);
 766
 767         /* Increment number of inner iterations */
 768         inneriter                  += j_index_end - j_index_start;
 769
 770         /* Outer loop uses 26 flops */
 771     }
 772
 773     /* Increment number of outer iterations */
 774     outeriter        += nri;
 775
 776     /* Update outer/inner flops */
 777
 778     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4_VF,outeriter*26 + inneriter*192);
 779 }
 780 /*
 781  * Gromacs nonbonded kernel:   nb_kernel_ElecCSTab_VdwCSTab_GeomW4P1_F_avx_256_double
 782  * Electrostatics interaction: CubicSplineTable
 783  * VdW interaction:            CubicSplineTable
 784  * Geometry:                   Water4-Particle
 785  * Calculate force/pot:        Force
 786  */
 787 void
 788 nb_kernel_ElecCSTab_VdwCSTab_GeomW4P1_F_avx_256_double
 789                     (t_nblist * gmx_restrict                nlist,
 790                      rvec * gmx_restrict                    xx,
 791                      rvec * gmx_restrict                    ff,
 792                      t_forcerec * gmx_restrict              fr,
 793                      t_mdatoms * gmx_restrict               mdatoms,
 794                      nb_kernel_data_t * gmx_restrict        kernel_data,
 795                      t_nrnb * gmx_restrict                  nrnb)
 796 {
 797     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
 798      * just 0 for non-waters.
 799      * Suffixes A,B,C,D refer to j loop unrolling done with AVX, e.g. for the four different
 800      * jnr indices corresponding to data put in the four positions in the SIMD register.
 801      */
 802     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
 803     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
 804     int              jnrA,jnrB,jnrC,jnrD;
 805     int              jnrlistA,jnrlistB,jnrlistC,jnrlistD;
 806     int              jnrlistE,jnrlistF,jnrlistG,jnrlistH;
 807     int              j_coord_offsetA,j_coord_offsetB,j_coord_offsetC,j_coord_offsetD;
 808     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
 809     real             rcutoff_scalar;
 810     real             *shiftvec,*fshift,*x,*f;
 811     real             *fjptrA,*fjptrB,*fjptrC,*fjptrD;
 812     real             scratch[4*DIM];
 813     __m256d          tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
 814     real *           vdwioffsetptr0;
 815     __m256d          ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
 816     real *           vdwioffsetptr1;
 817     __m256d          ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
 818     real *           vdwioffsetptr2;
 819     __m256d          ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
 820     real *           vdwioffsetptr3;
 821     __m256d          ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
 822     int              vdwjidx0A,vdwjidx0B,vdwjidx0C,vdwjidx0D;
 823     __m256d          jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
 824     __m256d          dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
 825     __m256d          dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
 826     __m256d          dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
 827     __m256d          dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
 828     __m256d          velec,felec,velecsum,facel,crf,krf,krf2;
 829     real             *charge;
 830     int              nvdwtype;
 831     __m256d          rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
 832     int              *vdwtype;
 833     real             *vdwparam;
 834     __m256d          one_sixth   = _mm256_set1_pd(1.0/6.0);
 835     __m256d          one_twelfth = _mm256_set1_pd(1.0/12.0);
 836     __m128i          vfitab;
 837     __m128i          ifour       = _mm_set1_epi32(4);
 838     __m256d          rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF;
 839     real             *vftab;
 840     __m256d          dummy_mask,cutoff_mask;
 841     __m128           tmpmask0,tmpmask1;
 842     __m256d          signbit = _mm256_castsi256_pd( _mm256_set1_epi32(0x80000000) );
 843     __m256d          one     = _mm256_set1_pd(1.0);
 844     __m256d          two     = _mm256_set1_pd(2.0);
 845     x                = xx[0];
 846     f                = ff[0];
 847
 848     nri              = nlist->nri;
 849     iinr             = nlist->iinr;
 850     jindex           = nlist->jindex;
 851     jjnr             = nlist->jjnr;
 852     shiftidx         = nlist->shift;
 853     gid              = nlist->gid;
 854     shiftvec         = fr->shift_vec[0];
 855     fshift           = fr->fshift[0];
 856     facel            = _mm256_set1_pd(fr->epsfac);
 857     charge           = mdatoms->chargeA;
 858     nvdwtype         = fr->ntype;
 859     vdwparam         = fr->nbfp;
 860     vdwtype          = mdatoms->typeA;
 861
 862     vftab            = kernel_data->table_elec_vdw->data;
 863     vftabscale       = _mm256_set1_pd(kernel_data->table_elec_vdw->scale);
 864
 865     /* Setup water-specific parameters */
 866     inr              = nlist->iinr[0];
 867     iq1              = _mm256_mul_pd(facel,_mm256_set1_pd(charge[inr+1]));
 868     iq2              = _mm256_mul_pd(facel,_mm256_set1_pd(charge[inr+2]));
 869     iq3              = _mm256_mul_pd(facel,_mm256_set1_pd(charge[inr+3]));
 870     vdwioffsetptr0   = vdwparam+2*nvdwtype*vdwtype[inr+0];
 871
 872     /* Avoid stupid compiler warnings */
 873     jnrA = jnrB = jnrC = jnrD = 0;
 874     j_coord_offsetA = 0;
 875     j_coord_offsetB = 0;
 876     j_coord_offsetC = 0;
 877     j_coord_offsetD = 0;
 878
 879     outeriter        = 0;
 880     inneriter        = 0;
 881
 882     for(iidx=0;iidx<4*DIM;iidx++)
 883     {
 884         scratch[iidx] = 0.0;
 885     }
 886
 887     /* Start outer loop over neighborlists */
 888     for(iidx=0; iidx<nri; iidx++)
 889     {
 890         /* Load shift vector for this list */
 891         i_shift_offset   = DIM*shiftidx[iidx];
 892
 893         /* Load limits for loop over neighbors */
 894         j_index_start    = jindex[iidx];
 895         j_index_end      = jindex[iidx+1];
 896
 897         /* Get outer coordinate index */
 898         inr              = iinr[iidx];
 899         i_coord_offset   = DIM*inr;
 900
 901         /* Load i particle coords and add shift vector */
 902         gmx_mm256_load_shift_and_4rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
 903                                                     &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
 904
 905         fix0             = _mm256_setzero_pd();
 906         fiy0             = _mm256_setzero_pd();
 907         fiz0             = _mm256_setzero_pd();
 908         fix1             = _mm256_setzero_pd();
 909         fiy1             = _mm256_setzero_pd();
 910         fiz1             = _mm256_setzero_pd();
 911         fix2             = _mm256_setzero_pd();
 912         fiy2             = _mm256_setzero_pd();
 913         fiz2             = _mm256_setzero_pd();
 914         fix3             = _mm256_setzero_pd();
 915         fiy3             = _mm256_setzero_pd();
 916         fiz3             = _mm256_setzero_pd();
 917
 918         /* Start inner kernel loop */
 919         for(jidx=j_index_start; jidx<j_index_end && jjnr[jidx+3]>=0; jidx+=4)
 920         {
 921
 922             /* Get j neighbor index, and coordinate index */
 923             jnrA             = jjnr[jidx];
 924             jnrB             = jjnr[jidx+1];
 925             jnrC             = jjnr[jidx+2];
 926             jnrD             = jjnr[jidx+3];
 927             j_coord_offsetA  = DIM*jnrA;
 928             j_coord_offsetB  = DIM*jnrB;
 929             j_coord_offsetC  = DIM*jnrC;
 930             j_coord_offsetD  = DIM*jnrD;
 931
 932             /* load j atom coordinates */
 933             gmx_mm256_load_1rvec_4ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
 934                                                  x+j_coord_offsetC,x+j_coord_offsetD,
 935                                                  &jx0,&jy0,&jz0);
 936
 937             /* Calculate displacement vector */
 938             dx00             = _mm256_sub_pd(ix0,jx0);
 939             dy00             = _mm256_sub_pd(iy0,jy0);
 940             dz00             = _mm256_sub_pd(iz0,jz0);
 941             dx10             = _mm256_sub_pd(ix1,jx0);
 942             dy10             = _mm256_sub_pd(iy1,jy0);
 943             dz10             = _mm256_sub_pd(iz1,jz0);
 944             dx20             = _mm256_sub_pd(ix2,jx0);
 945             dy20             = _mm256_sub_pd(iy2,jy0);
 946             dz20             = _mm256_sub_pd(iz2,jz0);
 947             dx30             = _mm256_sub_pd(ix3,jx0);
 948             dy30             = _mm256_sub_pd(iy3,jy0);
 949             dz30             = _mm256_sub_pd(iz3,jz0);
 950
 951             /* Calculate squared distance and things based on it */
 952             rsq00            = gmx_mm256_calc_rsq_pd(dx00,dy00,dz00);
 953             rsq10            = gmx_mm256_calc_rsq_pd(dx10,dy10,dz10);
 954             rsq20            = gmx_mm256_calc_rsq_pd(dx20,dy20,dz20);
 955             rsq30            = gmx_mm256_calc_rsq_pd(dx30,dy30,dz30);
 956
 957             rinv00           = gmx_mm256_invsqrt_pd(rsq00);
 958             rinv10           = gmx_mm256_invsqrt_pd(rsq10);
 959             rinv20           = gmx_mm256_invsqrt_pd(rsq20);
 960             rinv30           = gmx_mm256_invsqrt_pd(rsq30);
 961
 962             /* Load parameters for j particles */
 963             jq0              = gmx_mm256_load_4real_swizzle_pd(charge+jnrA+0,charge+jnrB+0,
 964                                                                  charge+jnrC+0,charge+jnrD+0);
 965             vdwjidx0A        = 2*vdwtype[jnrA+0];
 966             vdwjidx0B        = 2*vdwtype[jnrB+0];
 967             vdwjidx0C        = 2*vdwtype[jnrC+0];
 968             vdwjidx0D        = 2*vdwtype[jnrD+0];
 969
 970             fjx0             = _mm256_setzero_pd();
 971             fjy0             = _mm256_setzero_pd();
 972             fjz0             = _mm256_setzero_pd();
 973
 974             /**************************
 975              * CALCULATE INTERACTIONS *
 976              **************************/
 977
 978             r00              = _mm256_mul_pd(rsq00,rinv00);
 979
 980             /* Compute parameters for interactions between i and j atoms */
 981             gmx_mm256_load_4pair_swizzle_pd(vdwioffsetptr0+vdwjidx0A,
 982                                             vdwioffsetptr0+vdwjidx0B,
 983                                             vdwioffsetptr0+vdwjidx0C,
 984                                             vdwioffsetptr0+vdwjidx0D,
 985                                             &c6_00,&c12_00);
 986
 987             /* Calculate table index by multiplying r with table scale and truncate to integer */
 988             rt               = _mm256_mul_pd(r00,vftabscale);
 989             vfitab           = _mm256_cvttpd_epi32(rt);
 990             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
 991             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
 992
 993             /* CUBIC SPLINE TABLE DISPERSION */
 994             vfitab           = _mm_add_epi32(vfitab,ifour);
 995             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
 996             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
 997             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
 998             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
 999             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1000             Heps             = _mm256_mul_pd(vfeps,H);
1001             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1002             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1003             fvdw6            = _mm256_mul_pd(c6_00,FF);
1004
1005             /* CUBIC SPLINE TABLE REPULSION */
1006             vfitab           = _mm_add_epi32(vfitab,ifour);
1007             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1008             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1009             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1010             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1011             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1012             Heps             = _mm256_mul_pd(vfeps,H);
1013             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1014             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1015             fvdw12           = _mm256_mul_pd(c12_00,FF);
1016             fvdw             = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_add_pd(fvdw6,fvdw12),_mm256_mul_pd(vftabscale,rinv00)));
1017
1018             fscal            = fvdw;
1019
1020             /* Calculate temporary vectorial force */
1021             tx               = _mm256_mul_pd(fscal,dx00);
1022             ty               = _mm256_mul_pd(fscal,dy00);
1023             tz               = _mm256_mul_pd(fscal,dz00);
1024
1025             /* Update vectorial force */
1026             fix0             = _mm256_add_pd(fix0,tx);
1027             fiy0             = _mm256_add_pd(fiy0,ty);
1028             fiz0             = _mm256_add_pd(fiz0,tz);
1029
1030             fjx0             = _mm256_add_pd(fjx0,tx);
1031             fjy0             = _mm256_add_pd(fjy0,ty);
1032             fjz0             = _mm256_add_pd(fjz0,tz);
1033
1034             /**************************
1035              * CALCULATE INTERACTIONS *
1036              **************************/
1037
1038             r10              = _mm256_mul_pd(rsq10,rinv10);
1039
1040             /* Compute parameters for interactions between i and j atoms */
1041             qq10             = _mm256_mul_pd(iq1,jq0);
1042
1043             /* Calculate table index by multiplying r with table scale and truncate to integer */
1044             rt               = _mm256_mul_pd(r10,vftabscale);
1045             vfitab           = _mm256_cvttpd_epi32(rt);
1046             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
1047             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
1048
1049             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1050             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1051             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1052             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1053             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1054             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1055             Heps             = _mm256_mul_pd(vfeps,H);
1056             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1057             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1058             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq10,FF),_mm256_mul_pd(vftabscale,rinv10)));
1059
1060             fscal            = felec;
1061
1062             /* Calculate temporary vectorial force */
1063             tx               = _mm256_mul_pd(fscal,dx10);
1064             ty               = _mm256_mul_pd(fscal,dy10);
1065             tz               = _mm256_mul_pd(fscal,dz10);
1066
1067             /* Update vectorial force */
1068             fix1             = _mm256_add_pd(fix1,tx);
1069             fiy1             = _mm256_add_pd(fiy1,ty);
1070             fiz1             = _mm256_add_pd(fiz1,tz);
1071
1072             fjx0             = _mm256_add_pd(fjx0,tx);
1073             fjy0             = _mm256_add_pd(fjy0,ty);
1074             fjz0             = _mm256_add_pd(fjz0,tz);
1075
1076             /**************************
1077              * CALCULATE INTERACTIONS *
1078              **************************/
1079
1080             r20              = _mm256_mul_pd(rsq20,rinv20);
1081
1082             /* Compute parameters for interactions between i and j atoms */
1083             qq20             = _mm256_mul_pd(iq2,jq0);
1084
1085             /* Calculate table index by multiplying r with table scale and truncate to integer */
1086             rt               = _mm256_mul_pd(r20,vftabscale);
1087             vfitab           = _mm256_cvttpd_epi32(rt);
1088             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
1089             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
1090
1091             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1092             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1093             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1094             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1095             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1096             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1097             Heps             = _mm256_mul_pd(vfeps,H);
1098             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1099             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1100             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq20,FF),_mm256_mul_pd(vftabscale,rinv20)));
1101
1102             fscal            = felec;
1103
1104             /* Calculate temporary vectorial force */
1105             tx               = _mm256_mul_pd(fscal,dx20);
1106             ty               = _mm256_mul_pd(fscal,dy20);
1107             tz               = _mm256_mul_pd(fscal,dz20);
1108
1109             /* Update vectorial force */
1110             fix2             = _mm256_add_pd(fix2,tx);
1111             fiy2             = _mm256_add_pd(fiy2,ty);
1112             fiz2             = _mm256_add_pd(fiz2,tz);
1113
1114             fjx0             = _mm256_add_pd(fjx0,tx);
1115             fjy0             = _mm256_add_pd(fjy0,ty);
1116             fjz0             = _mm256_add_pd(fjz0,tz);
1117
1118             /**************************
1119              * CALCULATE INTERACTIONS *
1120              **************************/
1121
1122             r30              = _mm256_mul_pd(rsq30,rinv30);
1123
1124             /* Compute parameters for interactions between i and j atoms */
1125             qq30             = _mm256_mul_pd(iq3,jq0);
1126
1127             /* Calculate table index by multiplying r with table scale and truncate to integer */
1128             rt               = _mm256_mul_pd(r30,vftabscale);
1129             vfitab           = _mm256_cvttpd_epi32(rt);
1130             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
1131             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
1132
1133             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1134             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1135             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1136             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1137             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1138             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1139             Heps             = _mm256_mul_pd(vfeps,H);
1140             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1141             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1142             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq30,FF),_mm256_mul_pd(vftabscale,rinv30)));
1143
1144             fscal            = felec;
1145
1146             /* Calculate temporary vectorial force */
1147             tx               = _mm256_mul_pd(fscal,dx30);
1148             ty               = _mm256_mul_pd(fscal,dy30);
1149             tz               = _mm256_mul_pd(fscal,dz30);
1150
1151             /* Update vectorial force */
1152             fix3             = _mm256_add_pd(fix3,tx);
1153             fiy3             = _mm256_add_pd(fiy3,ty);
1154             fiz3             = _mm256_add_pd(fiz3,tz);
1155
1156             fjx0             = _mm256_add_pd(fjx0,tx);
1157             fjy0             = _mm256_add_pd(fjy0,ty);
1158             fjz0             = _mm256_add_pd(fjz0,tz);
1159
1160             fjptrA             = f+j_coord_offsetA;
1161             fjptrB             = f+j_coord_offsetB;
1162             fjptrC             = f+j_coord_offsetC;
1163             fjptrD             = f+j_coord_offsetD;
1164
1165             gmx_mm256_decrement_1rvec_4ptr_swizzle_pd(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
1166
1167             /* Inner loop uses 168 flops */
1168         }
1169
1170         if(jidx<j_index_end)
1171         {
1172
1173             /* Get j neighbor index, and coordinate index */
1174             jnrlistA         = jjnr[jidx];
1175             jnrlistB         = jjnr[jidx+1];
1176             jnrlistC         = jjnr[jidx+2];
1177             jnrlistD         = jjnr[jidx+3];
1178             /* Sign of each element will be negative for non-real atoms.
1179              * This mask will be 0xFFFFFFFF for dummy entries and 0x0 for real ones,
1180              * so use it as val = _mm_andnot_pd(mask,val) to clear dummy entries.
1181              */
1182             tmpmask0 = gmx_mm_castsi128_pd(_mm_cmplt_epi32(_mm_loadu_si128((const __m128i *)(jjnr+jidx)),_mm_setzero_si128()));
1183
1184             tmpmask1 = _mm_permute_ps(tmpmask0,_GMX_MM_PERMUTE(3,3,2,2));
1185             tmpmask0 = _mm_permute_ps(tmpmask0,_GMX_MM_PERMUTE(1,1,0,0));
1186             dummy_mask = _mm256_castps_pd(gmx_mm256_set_m128(tmpmask1,tmpmask0));
1187
1188             jnrA       = (jnrlistA>=0) ? jnrlistA : 0;
1189             jnrB       = (jnrlistB>=0) ? jnrlistB : 0;
1190             jnrC       = (jnrlistC>=0) ? jnrlistC : 0;
1191             jnrD       = (jnrlistD>=0) ? jnrlistD : 0;
1192             j_coord_offsetA  = DIM*jnrA;
1193             j_coord_offsetB  = DIM*jnrB;
1194             j_coord_offsetC  = DIM*jnrC;
1195             j_coord_offsetD  = DIM*jnrD;
1196
1197             /* load j atom coordinates */
1198             gmx_mm256_load_1rvec_4ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
1199                                                  x+j_coord_offsetC,x+j_coord_offsetD,
1200                                                  &jx0,&jy0,&jz0);
1201
1202             /* Calculate displacement vector */
1203             dx00             = _mm256_sub_pd(ix0,jx0);
1204             dy00             = _mm256_sub_pd(iy0,jy0);
1205             dz00             = _mm256_sub_pd(iz0,jz0);
1206             dx10             = _mm256_sub_pd(ix1,jx0);
1207             dy10             = _mm256_sub_pd(iy1,jy0);
1208             dz10             = _mm256_sub_pd(iz1,jz0);
1209             dx20             = _mm256_sub_pd(ix2,jx0);
1210             dy20             = _mm256_sub_pd(iy2,jy0);
1211             dz20             = _mm256_sub_pd(iz2,jz0);
1212             dx30             = _mm256_sub_pd(ix3,jx0);
1213             dy30             = _mm256_sub_pd(iy3,jy0);
1214             dz30             = _mm256_sub_pd(iz3,jz0);
1215
1216             /* Calculate squared distance and things based on it */
1217             rsq00            = gmx_mm256_calc_rsq_pd(dx00,dy00,dz00);
1218             rsq10            = gmx_mm256_calc_rsq_pd(dx10,dy10,dz10);
1219             rsq20            = gmx_mm256_calc_rsq_pd(dx20,dy20,dz20);
1220             rsq30            = gmx_mm256_calc_rsq_pd(dx30,dy30,dz30);
1221
1222             rinv00           = gmx_mm256_invsqrt_pd(rsq00);
1223             rinv10           = gmx_mm256_invsqrt_pd(rsq10);
1224             rinv20           = gmx_mm256_invsqrt_pd(rsq20);
1225             rinv30           = gmx_mm256_invsqrt_pd(rsq30);
1226
1227             /* Load parameters for j particles */
1228             jq0              = gmx_mm256_load_4real_swizzle_pd(charge+jnrA+0,charge+jnrB+0,
1229                                                                  charge+jnrC+0,charge+jnrD+0);
1230             vdwjidx0A        = 2*vdwtype[jnrA+0];
1231             vdwjidx0B        = 2*vdwtype[jnrB+0];
1232             vdwjidx0C        = 2*vdwtype[jnrC+0];
1233             vdwjidx0D        = 2*vdwtype[jnrD+0];
1234
1235             fjx0             = _mm256_setzero_pd();
1236             fjy0             = _mm256_setzero_pd();
1237             fjz0             = _mm256_setzero_pd();
1238
1239             /**************************
1240              * CALCULATE INTERACTIONS *
1241              **************************/
1242
1243             r00              = _mm256_mul_pd(rsq00,rinv00);
1244             r00              = _mm256_andnot_pd(dummy_mask,r00);
1245
1246             /* Compute parameters for interactions between i and j atoms */
1247             gmx_mm256_load_4pair_swizzle_pd(vdwioffsetptr0+vdwjidx0A,
1248                                             vdwioffsetptr0+vdwjidx0B,
1249                                             vdwioffsetptr0+vdwjidx0C,
1250                                             vdwioffsetptr0+vdwjidx0D,
1251                                             &c6_00,&c12_00);
1252
1253             /* Calculate table index by multiplying r with table scale and truncate to integer */
1254             rt               = _mm256_mul_pd(r00,vftabscale);
1255             vfitab           = _mm256_cvttpd_epi32(rt);
1256             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
1257             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
1258
1259             /* CUBIC SPLINE TABLE DISPERSION */
1260             vfitab           = _mm_add_epi32(vfitab,ifour);
1261             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1262             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1263             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1264             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1265             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1266             Heps             = _mm256_mul_pd(vfeps,H);
1267             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1268             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1269             fvdw6            = _mm256_mul_pd(c6_00,FF);
1270
1271             /* CUBIC SPLINE TABLE REPULSION */
1272             vfitab           = _mm_add_epi32(vfitab,ifour);
1273             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1274             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1275             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1276             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1277             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1278             Heps             = _mm256_mul_pd(vfeps,H);
1279             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1280             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1281             fvdw12           = _mm256_mul_pd(c12_00,FF);
1282             fvdw             = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_add_pd(fvdw6,fvdw12),_mm256_mul_pd(vftabscale,rinv00)));
1283
1284             fscal            = fvdw;
1285
1286             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
1287
1288             /* Calculate temporary vectorial force */
1289             tx               = _mm256_mul_pd(fscal,dx00);
1290             ty               = _mm256_mul_pd(fscal,dy00);
1291             tz               = _mm256_mul_pd(fscal,dz00);
1292
1293             /* Update vectorial force */
1294             fix0             = _mm256_add_pd(fix0,tx);
1295             fiy0             = _mm256_add_pd(fiy0,ty);
1296             fiz0             = _mm256_add_pd(fiz0,tz);
1297
1298             fjx0             = _mm256_add_pd(fjx0,tx);
1299             fjy0             = _mm256_add_pd(fjy0,ty);
1300             fjz0             = _mm256_add_pd(fjz0,tz);
1301
1302             /**************************
1303              * CALCULATE INTERACTIONS *
1304              **************************/
1305
1306             r10              = _mm256_mul_pd(rsq10,rinv10);
1307             r10              = _mm256_andnot_pd(dummy_mask,r10);
1308
1309             /* Compute parameters for interactions between i and j atoms */
1310             qq10             = _mm256_mul_pd(iq1,jq0);
1311
1312             /* Calculate table index by multiplying r with table scale and truncate to integer */
1313             rt               = _mm256_mul_pd(r10,vftabscale);
1314             vfitab           = _mm256_cvttpd_epi32(rt);
1315             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
1316             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
1317
1318             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1319             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1320             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1321             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1322             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1323             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1324             Heps             = _mm256_mul_pd(vfeps,H);
1325             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1326             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1327             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq10,FF),_mm256_mul_pd(vftabscale,rinv10)));
1328
1329             fscal            = felec;
1330
1331             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
1332
1333             /* Calculate temporary vectorial force */
1334             tx               = _mm256_mul_pd(fscal,dx10);
1335             ty               = _mm256_mul_pd(fscal,dy10);
1336             tz               = _mm256_mul_pd(fscal,dz10);
1337
1338             /* Update vectorial force */
1339             fix1             = _mm256_add_pd(fix1,tx);
1340             fiy1             = _mm256_add_pd(fiy1,ty);
1341             fiz1             = _mm256_add_pd(fiz1,tz);
1342
1343             fjx0             = _mm256_add_pd(fjx0,tx);
1344             fjy0             = _mm256_add_pd(fjy0,ty);
1345             fjz0             = _mm256_add_pd(fjz0,tz);
1346
1347             /**************************
1348              * CALCULATE INTERACTIONS *
1349              **************************/
1350
1351             r20              = _mm256_mul_pd(rsq20,rinv20);
1352             r20              = _mm256_andnot_pd(dummy_mask,r20);
1353
1354             /* Compute parameters for interactions between i and j atoms */
1355             qq20             = _mm256_mul_pd(iq2,jq0);
1356
1357             /* Calculate table index by multiplying r with table scale and truncate to integer */
1358             rt               = _mm256_mul_pd(r20,vftabscale);
1359             vfitab           = _mm256_cvttpd_epi32(rt);
1360             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
1361             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
1362
1363             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1364             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1365             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1366             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1367             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1368             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1369             Heps             = _mm256_mul_pd(vfeps,H);
1370             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1371             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1372             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq20,FF),_mm256_mul_pd(vftabscale,rinv20)));
1373
1374             fscal            = felec;
1375
1376             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
1377
1378             /* Calculate temporary vectorial force */
1379             tx               = _mm256_mul_pd(fscal,dx20);
1380             ty               = _mm256_mul_pd(fscal,dy20);
1381             tz               = _mm256_mul_pd(fscal,dz20);
1382
1383             /* Update vectorial force */
1384             fix2             = _mm256_add_pd(fix2,tx);
1385             fiy2             = _mm256_add_pd(fiy2,ty);
1386             fiz2             = _mm256_add_pd(fiz2,tz);
1387
1388             fjx0             = _mm256_add_pd(fjx0,tx);
1389             fjy0             = _mm256_add_pd(fjy0,ty);
1390             fjz0             = _mm256_add_pd(fjz0,tz);
1391
1392             /**************************
1393              * CALCULATE INTERACTIONS *
1394              **************************/
1395
1396             r30              = _mm256_mul_pd(rsq30,rinv30);
1397             r30              = _mm256_andnot_pd(dummy_mask,r30);
1398
1399             /* Compute parameters for interactions between i and j atoms */
1400             qq30             = _mm256_mul_pd(iq3,jq0);
1401
1402             /* Calculate table index by multiplying r with table scale and truncate to integer */
1403             rt               = _mm256_mul_pd(r30,vftabscale);
1404             vfitab           = _mm256_cvttpd_epi32(rt);
1405             vfeps            = _mm256_sub_pd(rt,_mm256_round_pd(rt, _MM_FROUND_FLOOR));
1406             vfitab           = _mm_slli_epi32(_mm_add_epi32(vfitab,_mm_slli_epi32(vfitab,1)),2);
1407
1408             /* CUBIC SPLINE TABLE ELECTROSTATICS */
1409             Y                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,0) );
1410             F                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,1) );
1411             G                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,2) );
1412             H                = _mm256_load_pd( vftab + _mm_extract_epi32(vfitab,3) );
1413             GMX_MM256_FULLTRANSPOSE4_PD(Y,F,G,H);
1414             Heps             = _mm256_mul_pd(vfeps,H);
1415             Fp               = _mm256_add_pd(F,_mm256_mul_pd(vfeps,_mm256_add_pd(G,Heps)));
1416             FF               = _mm256_add_pd(Fp,_mm256_mul_pd(vfeps,_mm256_add_pd(G,_mm256_add_pd(Heps,Heps))));
1417             felec            = _mm256_xor_pd(signbit,_mm256_mul_pd(_mm256_mul_pd(qq30,FF),_mm256_mul_pd(vftabscale,rinv30)));
1418
1419             fscal            = felec;
1420
1421             fscal            = _mm256_andnot_pd(dummy_mask,fscal);
1422
1423             /* Calculate temporary vectorial force */
1424             tx               = _mm256_mul_pd(fscal,dx30);
1425             ty               = _mm256_mul_pd(fscal,dy30);
1426             tz               = _mm256_mul_pd(fscal,dz30);
1427
1428             /* Update vectorial force */
1429             fix3             = _mm256_add_pd(fix3,tx);
1430             fiy3             = _mm256_add_pd(fiy3,ty);
1431             fiz3             = _mm256_add_pd(fiz3,tz);
1432
1433             fjx0             = _mm256_add_pd(fjx0,tx);
1434             fjy0             = _mm256_add_pd(fjy0,ty);
1435             fjz0             = _mm256_add_pd(fjz0,tz);
1436
1437             fjptrA             = (jnrlistA>=0) ? f+j_coord_offsetA : scratch;
1438             fjptrB             = (jnrlistB>=0) ? f+j_coord_offsetB : scratch;
1439             fjptrC             = (jnrlistC>=0) ? f+j_coord_offsetC : scratch;
1440             fjptrD             = (jnrlistD>=0) ? f+j_coord_offsetD : scratch;
1441
1442             gmx_mm256_decrement_1rvec_4ptr_swizzle_pd(fjptrA,fjptrB,fjptrC,fjptrD,fjx0,fjy0,fjz0);
1443
1444             /* Inner loop uses 172 flops */
1445         }
1446
1447         /* End of innermost loop */
1448
1449         gmx_mm256_update_iforce_4atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
1450                                                  f+i_coord_offset,fshift+i_shift_offset);
1451
1452         /* Increment number of inner iterations */
1453         inneriter                  += j_index_end - j_index_start;
1454
1455         /* Outer loop uses 24 flops */
1456     }
1457
1458     /* Increment number of outer iterations */
1459     outeriter        += nri;
1460
1461     /* Update outer/inner flops */
1462
1463     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W4_F,outeriter*24 + inneriter*172);
1464 }