src/gromacs/gmxlib/nonbonded/nb_kernel_avx_128_fma_double/nb_kernel_ElecCoul_VdwCSTab_GeomW4W4_avx_128_fma_double.c

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