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

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2012,2013,2014, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 /*
  36  * Note: this file was generated by the GROMACS avx_128_fma_double kernel generator.
  37  */
  38 #include "gmxpre.h"
  39
  40 #include "config.h"
  41
  42 #include <math.h>
  43
  44 #include "../nb_kernel.h"
  45 #include "gromacs/legacyheaders/types/simple.h"
  46 #include "gromacs/math/vec.h"
  47 #include "gromacs/legacyheaders/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_ElecRF_VdwLJ_GeomW3W3_VF_avx_128_fma_double
  54  * Electrostatics interaction: ReactionField
  55  * VdW interaction:            LennardJones
  56  * Geometry:                   Water3-Water3
  57  * Calculate force/pot:        PotentialAndForce
  58  */
  59 void
  60 nb_kernel_ElecRF_VdwLJ_GeomW3W3_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              vdwjidx0A,vdwjidx0B;
  89     __m128d          jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
  90     int              vdwjidx1A,vdwjidx1B;
  91     __m128d          jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
  92     int              vdwjidx2A,vdwjidx2B;
  93     __m128d          jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
  94     __m128d          dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
  95     __m128d          dx01,dy01,dz01,rsq01,rinv01,rinvsq01,r01,qq01,c6_01,c12_01;
  96     __m128d          dx02,dy02,dz02,rsq02,rinv02,rinvsq02,r02,qq02,c6_02,c12_02;
  97     __m128d          dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
  98     __m128d          dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
  99     __m128d          dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
 100     __m128d          dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
 101     __m128d          dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
 102     __m128d          dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
 103     __m128d          velec,felec,velecsum,facel,crf,krf,krf2;
 104     real             *charge;
 105     int              nvdwtype;
 106     __m128d          rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
 107     int              *vdwtype;
 108     real             *vdwparam;
 109     __m128d          one_sixth   = _mm_set1_pd(1.0/6.0);
 110     __m128d          one_twelfth = _mm_set1_pd(1.0/12.0);
 111     __m128d          dummy_mask,cutoff_mask;
 112     __m128d          signbit   = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
 113     __m128d          one     = _mm_set1_pd(1.0);
 114     __m128d          two     = _mm_set1_pd(2.0);
 115     x                = xx[0];
 116     f                = ff[0];
 117
 118     nri              = nlist->nri;
 119     iinr             = nlist->iinr;
 120     jindex           = nlist->jindex;
 121     jjnr             = nlist->jjnr;
 122     shiftidx         = nlist->shift;
 123     gid              = nlist->gid;
 124     shiftvec         = fr->shift_vec[0];
 125     fshift           = fr->fshift[0];
 126     facel            = _mm_set1_pd(fr->epsfac);
 127     charge           = mdatoms->chargeA;
 128     krf              = _mm_set1_pd(fr->ic->k_rf);
 129     krf2             = _mm_set1_pd(fr->ic->k_rf*2.0);
 130     crf              = _mm_set1_pd(fr->ic->c_rf);
 131     nvdwtype         = fr->ntype;
 132     vdwparam         = fr->nbfp;
 133     vdwtype          = mdatoms->typeA;
 134
 135     /* Setup water-specific parameters */
 136     inr              = nlist->iinr[0];
 137     iq0              = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+0]));
 138     iq1              = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
 139     iq2              = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
 140     vdwioffset0      = 2*nvdwtype*vdwtype[inr+0];
 141
 142     jq0              = _mm_set1_pd(charge[inr+0]);
 143     jq1              = _mm_set1_pd(charge[inr+1]);
 144     jq2              = _mm_set1_pd(charge[inr+2]);
 145     vdwjidx0A        = 2*vdwtype[inr+0];
 146     qq00             = _mm_mul_pd(iq0,jq0);
 147     c6_00            = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A]);
 148     c12_00           = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A+1]);
 149     qq01             = _mm_mul_pd(iq0,jq1);
 150     qq02             = _mm_mul_pd(iq0,jq2);
 151     qq10             = _mm_mul_pd(iq1,jq0);
 152     qq11             = _mm_mul_pd(iq1,jq1);
 153     qq12             = _mm_mul_pd(iq1,jq2);
 154     qq20             = _mm_mul_pd(iq2,jq0);
 155     qq21             = _mm_mul_pd(iq2,jq1);
 156     qq22             = _mm_mul_pd(iq2,jq2);
 157
 158     /* Avoid stupid compiler warnings */
 159     jnrA = jnrB = 0;
 160     j_coord_offsetA = 0;
 161     j_coord_offsetB = 0;
 162
 163     outeriter        = 0;
 164     inneriter        = 0;
 165
 166     /* Start outer loop over neighborlists */
 167     for(iidx=0; iidx<nri; iidx++)
 168     {
 169         /* Load shift vector for this list */
 170         i_shift_offset   = DIM*shiftidx[iidx];
 171
 172         /* Load limits for loop over neighbors */
 173         j_index_start    = jindex[iidx];
 174         j_index_end      = jindex[iidx+1];
 175
 176         /* Get outer coordinate index */
 177         inr              = iinr[iidx];
 178         i_coord_offset   = DIM*inr;
 179
 180         /* Load i particle coords and add shift vector */
 181         gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
 182                                                  &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2);
 183
 184         fix0             = _mm_setzero_pd();
 185         fiy0             = _mm_setzero_pd();
 186         fiz0             = _mm_setzero_pd();
 187         fix1             = _mm_setzero_pd();
 188         fiy1             = _mm_setzero_pd();
 189         fiz1             = _mm_setzero_pd();
 190         fix2             = _mm_setzero_pd();
 191         fiy2             = _mm_setzero_pd();
 192         fiz2             = _mm_setzero_pd();
 193
 194         /* Reset potential sums */
 195         velecsum         = _mm_setzero_pd();
 196         vvdwsum          = _mm_setzero_pd();
 197
 198         /* Start inner kernel loop */
 199         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
 200         {
 201
 202             /* Get j neighbor index, and coordinate index */
 203             jnrA             = jjnr[jidx];
 204             jnrB             = jjnr[jidx+1];
 205             j_coord_offsetA  = DIM*jnrA;
 206             j_coord_offsetB  = DIM*jnrB;
 207
 208             /* load j atom coordinates */
 209             gmx_mm_load_3rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
 210                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
 211
 212             /* Calculate displacement vector */
 213             dx00             = _mm_sub_pd(ix0,jx0);
 214             dy00             = _mm_sub_pd(iy0,jy0);
 215             dz00             = _mm_sub_pd(iz0,jz0);
 216             dx01             = _mm_sub_pd(ix0,jx1);
 217             dy01             = _mm_sub_pd(iy0,jy1);
 218             dz01             = _mm_sub_pd(iz0,jz1);
 219             dx02             = _mm_sub_pd(ix0,jx2);
 220             dy02             = _mm_sub_pd(iy0,jy2);
 221             dz02             = _mm_sub_pd(iz0,jz2);
 222             dx10             = _mm_sub_pd(ix1,jx0);
 223             dy10             = _mm_sub_pd(iy1,jy0);
 224             dz10             = _mm_sub_pd(iz1,jz0);
 225             dx11             = _mm_sub_pd(ix1,jx1);
 226             dy11             = _mm_sub_pd(iy1,jy1);
 227             dz11             = _mm_sub_pd(iz1,jz1);
 228             dx12             = _mm_sub_pd(ix1,jx2);
 229             dy12             = _mm_sub_pd(iy1,jy2);
 230             dz12             = _mm_sub_pd(iz1,jz2);
 231             dx20             = _mm_sub_pd(ix2,jx0);
 232             dy20             = _mm_sub_pd(iy2,jy0);
 233             dz20             = _mm_sub_pd(iz2,jz0);
 234             dx21             = _mm_sub_pd(ix2,jx1);
 235             dy21             = _mm_sub_pd(iy2,jy1);
 236             dz21             = _mm_sub_pd(iz2,jz1);
 237             dx22             = _mm_sub_pd(ix2,jx2);
 238             dy22             = _mm_sub_pd(iy2,jy2);
 239             dz22             = _mm_sub_pd(iz2,jz2);
 240
 241             /* Calculate squared distance and things based on it */
 242             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
 243             rsq01            = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
 244             rsq02            = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
 245             rsq10            = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
 246             rsq11            = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
 247             rsq12            = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
 248             rsq20            = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
 249             rsq21            = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
 250             rsq22            = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
 251
 252             rinv00           = gmx_mm_invsqrt_pd(rsq00);
 253             rinv01           = gmx_mm_invsqrt_pd(rsq01);
 254             rinv02           = gmx_mm_invsqrt_pd(rsq02);
 255             rinv10           = gmx_mm_invsqrt_pd(rsq10);
 256             rinv11           = gmx_mm_invsqrt_pd(rsq11);
 257             rinv12           = gmx_mm_invsqrt_pd(rsq12);
 258             rinv20           = gmx_mm_invsqrt_pd(rsq20);
 259             rinv21           = gmx_mm_invsqrt_pd(rsq21);
 260             rinv22           = gmx_mm_invsqrt_pd(rsq22);
 261
 262             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
 263             rinvsq01         = _mm_mul_pd(rinv01,rinv01);
 264             rinvsq02         = _mm_mul_pd(rinv02,rinv02);
 265             rinvsq10         = _mm_mul_pd(rinv10,rinv10);
 266             rinvsq11         = _mm_mul_pd(rinv11,rinv11);
 267             rinvsq12         = _mm_mul_pd(rinv12,rinv12);
 268             rinvsq20         = _mm_mul_pd(rinv20,rinv20);
 269             rinvsq21         = _mm_mul_pd(rinv21,rinv21);
 270             rinvsq22         = _mm_mul_pd(rinv22,rinv22);
 271
 272             fjx0             = _mm_setzero_pd();
 273             fjy0             = _mm_setzero_pd();
 274             fjz0             = _mm_setzero_pd();
 275             fjx1             = _mm_setzero_pd();
 276             fjy1             = _mm_setzero_pd();
 277             fjz1             = _mm_setzero_pd();
 278             fjx2             = _mm_setzero_pd();
 279             fjy2             = _mm_setzero_pd();
 280             fjz2             = _mm_setzero_pd();
 281
 282             /**************************
 283              * CALCULATE INTERACTIONS *
 284              **************************/
 285
 286             /* REACTION-FIELD ELECTROSTATICS */
 287             velec            = _mm_mul_pd(qq00,_mm_sub_pd(_mm_macc_pd(krf,rsq00,rinv00),crf));
 288             felec            = _mm_mul_pd(qq00,_mm_msub_pd(rinv00,rinvsq00,krf2));
 289
 290             /* LENNARD-JONES DISPERSION/REPULSION */
 291
 292             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
 293             vvdw6            = _mm_mul_pd(c6_00,rinvsix);
 294             vvdw12           = _mm_mul_pd(c12_00,_mm_mul_pd(rinvsix,rinvsix));
 295             vvdw             = _mm_msub_pd( vvdw12,one_twelfth, _mm_mul_pd(vvdw6,one_sixth) );
 296             fvdw             = _mm_mul_pd(_mm_sub_pd(vvdw12,vvdw6),rinvsq00);
 297
 298             /* Update potential sum for this i atom from the interaction with this j atom. */
 299             velecsum         = _mm_add_pd(velecsum,velec);
 300             vvdwsum          = _mm_add_pd(vvdwsum,vvdw);
 301
 302             fscal            = _mm_add_pd(felec,fvdw);
 303
 304             /* Update vectorial force */
 305             fix0             = _mm_macc_pd(dx00,fscal,fix0);
 306             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
 307             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
 308
 309             fjx0             = _mm_macc_pd(dx00,fscal,fjx0);
 310             fjy0             = _mm_macc_pd(dy00,fscal,fjy0);
 311             fjz0             = _mm_macc_pd(dz00,fscal,fjz0);
 312
 313             /**************************
 314              * CALCULATE INTERACTIONS *
 315              **************************/
 316
 317             /* REACTION-FIELD ELECTROSTATICS */
 318             velec            = _mm_mul_pd(qq01,_mm_sub_pd(_mm_macc_pd(krf,rsq01,rinv01),crf));
 319             felec            = _mm_mul_pd(qq01,_mm_msub_pd(rinv01,rinvsq01,krf2));
 320
 321             /* Update potential sum for this i atom from the interaction with this j atom. */
 322             velecsum         = _mm_add_pd(velecsum,velec);
 323
 324             fscal            = felec;
 325
 326             /* Update vectorial force */
 327             fix0             = _mm_macc_pd(dx01,fscal,fix0);
 328             fiy0             = _mm_macc_pd(dy01,fscal,fiy0);
 329             fiz0             = _mm_macc_pd(dz01,fscal,fiz0);
 330
 331             fjx1             = _mm_macc_pd(dx01,fscal,fjx1);
 332             fjy1             = _mm_macc_pd(dy01,fscal,fjy1);
 333             fjz1             = _mm_macc_pd(dz01,fscal,fjz1);
 334
 335             /**************************
 336              * CALCULATE INTERACTIONS *
 337              **************************/
 338
 339             /* REACTION-FIELD ELECTROSTATICS */
 340             velec            = _mm_mul_pd(qq02,_mm_sub_pd(_mm_macc_pd(krf,rsq02,rinv02),crf));
 341             felec            = _mm_mul_pd(qq02,_mm_msub_pd(rinv02,rinvsq02,krf2));
 342
 343             /* Update potential sum for this i atom from the interaction with this j atom. */
 344             velecsum         = _mm_add_pd(velecsum,velec);
 345
 346             fscal            = felec;
 347
 348             /* Update vectorial force */
 349             fix0             = _mm_macc_pd(dx02,fscal,fix0);
 350             fiy0             = _mm_macc_pd(dy02,fscal,fiy0);
 351             fiz0             = _mm_macc_pd(dz02,fscal,fiz0);
 352
 353             fjx2             = _mm_macc_pd(dx02,fscal,fjx2);
 354             fjy2             = _mm_macc_pd(dy02,fscal,fjy2);
 355             fjz2             = _mm_macc_pd(dz02,fscal,fjz2);
 356
 357             /**************************
 358              * CALCULATE INTERACTIONS *
 359              **************************/
 360
 361             /* REACTION-FIELD ELECTROSTATICS */
 362             velec            = _mm_mul_pd(qq10,_mm_sub_pd(_mm_macc_pd(krf,rsq10,rinv10),crf));
 363             felec            = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
 364
 365             /* Update potential sum for this i atom from the interaction with this j atom. */
 366             velecsum         = _mm_add_pd(velecsum,velec);
 367
 368             fscal            = felec;
 369
 370             /* Update vectorial force */
 371             fix1             = _mm_macc_pd(dx10,fscal,fix1);
 372             fiy1             = _mm_macc_pd(dy10,fscal,fiy1);
 373             fiz1             = _mm_macc_pd(dz10,fscal,fiz1);
 374
 375             fjx0             = _mm_macc_pd(dx10,fscal,fjx0);
 376             fjy0             = _mm_macc_pd(dy10,fscal,fjy0);
 377             fjz0             = _mm_macc_pd(dz10,fscal,fjz0);
 378
 379             /**************************
 380              * CALCULATE INTERACTIONS *
 381              **************************/
 382
 383             /* REACTION-FIELD ELECTROSTATICS */
 384             velec            = _mm_mul_pd(qq11,_mm_sub_pd(_mm_macc_pd(krf,rsq11,rinv11),crf));
 385             felec            = _mm_mul_pd(qq11,_mm_msub_pd(rinv11,rinvsq11,krf2));
 386
 387             /* Update potential sum for this i atom from the interaction with this j atom. */
 388             velecsum         = _mm_add_pd(velecsum,velec);
 389
 390             fscal            = felec;
 391
 392             /* Update vectorial force */
 393             fix1             = _mm_macc_pd(dx11,fscal,fix1);
 394             fiy1             = _mm_macc_pd(dy11,fscal,fiy1);
 395             fiz1             = _mm_macc_pd(dz11,fscal,fiz1);
 396
 397             fjx1             = _mm_macc_pd(dx11,fscal,fjx1);
 398             fjy1             = _mm_macc_pd(dy11,fscal,fjy1);
 399             fjz1             = _mm_macc_pd(dz11,fscal,fjz1);
 400
 401             /**************************
 402              * CALCULATE INTERACTIONS *
 403              **************************/
 404
 405             /* REACTION-FIELD ELECTROSTATICS */
 406             velec            = _mm_mul_pd(qq12,_mm_sub_pd(_mm_macc_pd(krf,rsq12,rinv12),crf));
 407             felec            = _mm_mul_pd(qq12,_mm_msub_pd(rinv12,rinvsq12,krf2));
 408
 409             /* Update potential sum for this i atom from the interaction with this j atom. */
 410             velecsum         = _mm_add_pd(velecsum,velec);
 411
 412             fscal            = felec;
 413
 414             /* Update vectorial force */
 415             fix1             = _mm_macc_pd(dx12,fscal,fix1);
 416             fiy1             = _mm_macc_pd(dy12,fscal,fiy1);
 417             fiz1             = _mm_macc_pd(dz12,fscal,fiz1);
 418
 419             fjx2             = _mm_macc_pd(dx12,fscal,fjx2);
 420             fjy2             = _mm_macc_pd(dy12,fscal,fjy2);
 421             fjz2             = _mm_macc_pd(dz12,fscal,fjz2);
 422
 423             /**************************
 424              * CALCULATE INTERACTIONS *
 425              **************************/
 426
 427             /* REACTION-FIELD ELECTROSTATICS */
 428             velec            = _mm_mul_pd(qq20,_mm_sub_pd(_mm_macc_pd(krf,rsq20,rinv20),crf));
 429             felec            = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
 430
 431             /* Update potential sum for this i atom from the interaction with this j atom. */
 432             velecsum         = _mm_add_pd(velecsum,velec);
 433
 434             fscal            = felec;
 435
 436             /* Update vectorial force */
 437             fix2             = _mm_macc_pd(dx20,fscal,fix2);
 438             fiy2             = _mm_macc_pd(dy20,fscal,fiy2);
 439             fiz2             = _mm_macc_pd(dz20,fscal,fiz2);
 440
 441             fjx0             = _mm_macc_pd(dx20,fscal,fjx0);
 442             fjy0             = _mm_macc_pd(dy20,fscal,fjy0);
 443             fjz0             = _mm_macc_pd(dz20,fscal,fjz0);
 444
 445             /**************************
 446              * CALCULATE INTERACTIONS *
 447              **************************/
 448
 449             /* REACTION-FIELD ELECTROSTATICS */
 450             velec            = _mm_mul_pd(qq21,_mm_sub_pd(_mm_macc_pd(krf,rsq21,rinv21),crf));
 451             felec            = _mm_mul_pd(qq21,_mm_msub_pd(rinv21,rinvsq21,krf2));
 452
 453             /* Update potential sum for this i atom from the interaction with this j atom. */
 454             velecsum         = _mm_add_pd(velecsum,velec);
 455
 456             fscal            = felec;
 457
 458             /* Update vectorial force */
 459             fix2             = _mm_macc_pd(dx21,fscal,fix2);
 460             fiy2             = _mm_macc_pd(dy21,fscal,fiy2);
 461             fiz2             = _mm_macc_pd(dz21,fscal,fiz2);
 462
 463             fjx1             = _mm_macc_pd(dx21,fscal,fjx1);
 464             fjy1             = _mm_macc_pd(dy21,fscal,fjy1);
 465             fjz1             = _mm_macc_pd(dz21,fscal,fjz1);
 466
 467             /**************************
 468              * CALCULATE INTERACTIONS *
 469              **************************/
 470
 471             /* REACTION-FIELD ELECTROSTATICS */
 472             velec            = _mm_mul_pd(qq22,_mm_sub_pd(_mm_macc_pd(krf,rsq22,rinv22),crf));
 473             felec            = _mm_mul_pd(qq22,_mm_msub_pd(rinv22,rinvsq22,krf2));
 474
 475             /* Update potential sum for this i atom from the interaction with this j atom. */
 476             velecsum         = _mm_add_pd(velecsum,velec);
 477
 478             fscal            = felec;
 479
 480             /* Update vectorial force */
 481             fix2             = _mm_macc_pd(dx22,fscal,fix2);
 482             fiy2             = _mm_macc_pd(dy22,fscal,fiy2);
 483             fiz2             = _mm_macc_pd(dz22,fscal,fiz2);
 484
 485             fjx2             = _mm_macc_pd(dx22,fscal,fjx2);
 486             fjy2             = _mm_macc_pd(dy22,fscal,fjy2);
 487             fjz2             = _mm_macc_pd(dz22,fscal,fjz2);
 488
 489             gmx_mm_decrement_3rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
 490
 491             /* Inner loop uses 327 flops */
 492         }
 493
 494         if(jidx<j_index_end)
 495         {
 496
 497             jnrA             = jjnr[jidx];
 498             j_coord_offsetA  = DIM*jnrA;
 499
 500             /* load j atom coordinates */
 501             gmx_mm_load_3rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
 502                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
 503
 504             /* Calculate displacement vector */
 505             dx00             = _mm_sub_pd(ix0,jx0);
 506             dy00             = _mm_sub_pd(iy0,jy0);
 507             dz00             = _mm_sub_pd(iz0,jz0);
 508             dx01             = _mm_sub_pd(ix0,jx1);
 509             dy01             = _mm_sub_pd(iy0,jy1);
 510             dz01             = _mm_sub_pd(iz0,jz1);
 511             dx02             = _mm_sub_pd(ix0,jx2);
 512             dy02             = _mm_sub_pd(iy0,jy2);
 513             dz02             = _mm_sub_pd(iz0,jz2);
 514             dx10             = _mm_sub_pd(ix1,jx0);
 515             dy10             = _mm_sub_pd(iy1,jy0);
 516             dz10             = _mm_sub_pd(iz1,jz0);
 517             dx11             = _mm_sub_pd(ix1,jx1);
 518             dy11             = _mm_sub_pd(iy1,jy1);
 519             dz11             = _mm_sub_pd(iz1,jz1);
 520             dx12             = _mm_sub_pd(ix1,jx2);
 521             dy12             = _mm_sub_pd(iy1,jy2);
 522             dz12             = _mm_sub_pd(iz1,jz2);
 523             dx20             = _mm_sub_pd(ix2,jx0);
 524             dy20             = _mm_sub_pd(iy2,jy0);
 525             dz20             = _mm_sub_pd(iz2,jz0);
 526             dx21             = _mm_sub_pd(ix2,jx1);
 527             dy21             = _mm_sub_pd(iy2,jy1);
 528             dz21             = _mm_sub_pd(iz2,jz1);
 529             dx22             = _mm_sub_pd(ix2,jx2);
 530             dy22             = _mm_sub_pd(iy2,jy2);
 531             dz22             = _mm_sub_pd(iz2,jz2);
 532
 533             /* Calculate squared distance and things based on it */
 534             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
 535             rsq01            = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
 536             rsq02            = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
 537             rsq10            = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
 538             rsq11            = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
 539             rsq12            = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
 540             rsq20            = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
 541             rsq21            = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
 542             rsq22            = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
 543
 544             rinv00           = gmx_mm_invsqrt_pd(rsq00);
 545             rinv01           = gmx_mm_invsqrt_pd(rsq01);
 546             rinv02           = gmx_mm_invsqrt_pd(rsq02);
 547             rinv10           = gmx_mm_invsqrt_pd(rsq10);
 548             rinv11           = gmx_mm_invsqrt_pd(rsq11);
 549             rinv12           = gmx_mm_invsqrt_pd(rsq12);
 550             rinv20           = gmx_mm_invsqrt_pd(rsq20);
 551             rinv21           = gmx_mm_invsqrt_pd(rsq21);
 552             rinv22           = gmx_mm_invsqrt_pd(rsq22);
 553
 554             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
 555             rinvsq01         = _mm_mul_pd(rinv01,rinv01);
 556             rinvsq02         = _mm_mul_pd(rinv02,rinv02);
 557             rinvsq10         = _mm_mul_pd(rinv10,rinv10);
 558             rinvsq11         = _mm_mul_pd(rinv11,rinv11);
 559             rinvsq12         = _mm_mul_pd(rinv12,rinv12);
 560             rinvsq20         = _mm_mul_pd(rinv20,rinv20);
 561             rinvsq21         = _mm_mul_pd(rinv21,rinv21);
 562             rinvsq22         = _mm_mul_pd(rinv22,rinv22);
 563
 564             fjx0             = _mm_setzero_pd();
 565             fjy0             = _mm_setzero_pd();
 566             fjz0             = _mm_setzero_pd();
 567             fjx1             = _mm_setzero_pd();
 568             fjy1             = _mm_setzero_pd();
 569             fjz1             = _mm_setzero_pd();
 570             fjx2             = _mm_setzero_pd();
 571             fjy2             = _mm_setzero_pd();
 572             fjz2             = _mm_setzero_pd();
 573
 574             /**************************
 575              * CALCULATE INTERACTIONS *
 576              **************************/
 577
 578             /* REACTION-FIELD ELECTROSTATICS */
 579             velec            = _mm_mul_pd(qq00,_mm_sub_pd(_mm_macc_pd(krf,rsq00,rinv00),crf));
 580             felec            = _mm_mul_pd(qq00,_mm_msub_pd(rinv00,rinvsq00,krf2));
 581
 582             /* LENNARD-JONES DISPERSION/REPULSION */
 583
 584             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
 585             vvdw6            = _mm_mul_pd(c6_00,rinvsix);
 586             vvdw12           = _mm_mul_pd(c12_00,_mm_mul_pd(rinvsix,rinvsix));
 587             vvdw             = _mm_msub_pd( vvdw12,one_twelfth, _mm_mul_pd(vvdw6,one_sixth) );
 588             fvdw             = _mm_mul_pd(_mm_sub_pd(vvdw12,vvdw6),rinvsq00);
 589
 590             /* Update potential sum for this i atom from the interaction with this j atom. */
 591             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 592             velecsum         = _mm_add_pd(velecsum,velec);
 593             vvdw             = _mm_unpacklo_pd(vvdw,_mm_setzero_pd());
 594             vvdwsum          = _mm_add_pd(vvdwsum,vvdw);
 595
 596             fscal            = _mm_add_pd(felec,fvdw);
 597
 598             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 599
 600             /* Update vectorial force */
 601             fix0             = _mm_macc_pd(dx00,fscal,fix0);
 602             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
 603             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
 604
 605             fjx0             = _mm_macc_pd(dx00,fscal,fjx0);
 606             fjy0             = _mm_macc_pd(dy00,fscal,fjy0);
 607             fjz0             = _mm_macc_pd(dz00,fscal,fjz0);
 608
 609             /**************************
 610              * CALCULATE INTERACTIONS *
 611              **************************/
 612
 613             /* REACTION-FIELD ELECTROSTATICS */
 614             velec            = _mm_mul_pd(qq01,_mm_sub_pd(_mm_macc_pd(krf,rsq01,rinv01),crf));
 615             felec            = _mm_mul_pd(qq01,_mm_msub_pd(rinv01,rinvsq01,krf2));
 616
 617             /* Update potential sum for this i atom from the interaction with this j atom. */
 618             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 619             velecsum         = _mm_add_pd(velecsum,velec);
 620
 621             fscal            = felec;
 622
 623             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 624
 625             /* Update vectorial force */
 626             fix0             = _mm_macc_pd(dx01,fscal,fix0);
 627             fiy0             = _mm_macc_pd(dy01,fscal,fiy0);
 628             fiz0             = _mm_macc_pd(dz01,fscal,fiz0);
 629
 630             fjx1             = _mm_macc_pd(dx01,fscal,fjx1);
 631             fjy1             = _mm_macc_pd(dy01,fscal,fjy1);
 632             fjz1             = _mm_macc_pd(dz01,fscal,fjz1);
 633
 634             /**************************
 635              * CALCULATE INTERACTIONS *
 636              **************************/
 637
 638             /* REACTION-FIELD ELECTROSTATICS */
 639             velec            = _mm_mul_pd(qq02,_mm_sub_pd(_mm_macc_pd(krf,rsq02,rinv02),crf));
 640             felec            = _mm_mul_pd(qq02,_mm_msub_pd(rinv02,rinvsq02,krf2));
 641
 642             /* Update potential sum for this i atom from the interaction with this j atom. */
 643             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 644             velecsum         = _mm_add_pd(velecsum,velec);
 645
 646             fscal            = felec;
 647
 648             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 649
 650             /* Update vectorial force */
 651             fix0             = _mm_macc_pd(dx02,fscal,fix0);
 652             fiy0             = _mm_macc_pd(dy02,fscal,fiy0);
 653             fiz0             = _mm_macc_pd(dz02,fscal,fiz0);
 654
 655             fjx2             = _mm_macc_pd(dx02,fscal,fjx2);
 656             fjy2             = _mm_macc_pd(dy02,fscal,fjy2);
 657             fjz2             = _mm_macc_pd(dz02,fscal,fjz2);
 658
 659             /**************************
 660              * CALCULATE INTERACTIONS *
 661              **************************/
 662
 663             /* REACTION-FIELD ELECTROSTATICS */
 664             velec            = _mm_mul_pd(qq10,_mm_sub_pd(_mm_macc_pd(krf,rsq10,rinv10),crf));
 665             felec            = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
 666
 667             /* Update potential sum for this i atom from the interaction with this j atom. */
 668             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 669             velecsum         = _mm_add_pd(velecsum,velec);
 670
 671             fscal            = felec;
 672
 673             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 674
 675             /* Update vectorial force */
 676             fix1             = _mm_macc_pd(dx10,fscal,fix1);
 677             fiy1             = _mm_macc_pd(dy10,fscal,fiy1);
 678             fiz1             = _mm_macc_pd(dz10,fscal,fiz1);
 679
 680             fjx0             = _mm_macc_pd(dx10,fscal,fjx0);
 681             fjy0             = _mm_macc_pd(dy10,fscal,fjy0);
 682             fjz0             = _mm_macc_pd(dz10,fscal,fjz0);
 683
 684             /**************************
 685              * CALCULATE INTERACTIONS *
 686              **************************/
 687
 688             /* REACTION-FIELD ELECTROSTATICS */
 689             velec            = _mm_mul_pd(qq11,_mm_sub_pd(_mm_macc_pd(krf,rsq11,rinv11),crf));
 690             felec            = _mm_mul_pd(qq11,_mm_msub_pd(rinv11,rinvsq11,krf2));
 691
 692             /* Update potential sum for this i atom from the interaction with this j atom. */
 693             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 694             velecsum         = _mm_add_pd(velecsum,velec);
 695
 696             fscal            = felec;
 697
 698             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 699
 700             /* Update vectorial force */
 701             fix1             = _mm_macc_pd(dx11,fscal,fix1);
 702             fiy1             = _mm_macc_pd(dy11,fscal,fiy1);
 703             fiz1             = _mm_macc_pd(dz11,fscal,fiz1);
 704
 705             fjx1             = _mm_macc_pd(dx11,fscal,fjx1);
 706             fjy1             = _mm_macc_pd(dy11,fscal,fjy1);
 707             fjz1             = _mm_macc_pd(dz11,fscal,fjz1);
 708
 709             /**************************
 710              * CALCULATE INTERACTIONS *
 711              **************************/
 712
 713             /* REACTION-FIELD ELECTROSTATICS */
 714             velec            = _mm_mul_pd(qq12,_mm_sub_pd(_mm_macc_pd(krf,rsq12,rinv12),crf));
 715             felec            = _mm_mul_pd(qq12,_mm_msub_pd(rinv12,rinvsq12,krf2));
 716
 717             /* Update potential sum for this i atom from the interaction with this j atom. */
 718             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 719             velecsum         = _mm_add_pd(velecsum,velec);
 720
 721             fscal            = felec;
 722
 723             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 724
 725             /* Update vectorial force */
 726             fix1             = _mm_macc_pd(dx12,fscal,fix1);
 727             fiy1             = _mm_macc_pd(dy12,fscal,fiy1);
 728             fiz1             = _mm_macc_pd(dz12,fscal,fiz1);
 729
 730             fjx2             = _mm_macc_pd(dx12,fscal,fjx2);
 731             fjy2             = _mm_macc_pd(dy12,fscal,fjy2);
 732             fjz2             = _mm_macc_pd(dz12,fscal,fjz2);
 733
 734             /**************************
 735              * CALCULATE INTERACTIONS *
 736              **************************/
 737
 738             /* REACTION-FIELD ELECTROSTATICS */
 739             velec            = _mm_mul_pd(qq20,_mm_sub_pd(_mm_macc_pd(krf,rsq20,rinv20),crf));
 740             felec            = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
 741
 742             /* Update potential sum for this i atom from the interaction with this j atom. */
 743             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 744             velecsum         = _mm_add_pd(velecsum,velec);
 745
 746             fscal            = felec;
 747
 748             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 749
 750             /* Update vectorial force */
 751             fix2             = _mm_macc_pd(dx20,fscal,fix2);
 752             fiy2             = _mm_macc_pd(dy20,fscal,fiy2);
 753             fiz2             = _mm_macc_pd(dz20,fscal,fiz2);
 754
 755             fjx0             = _mm_macc_pd(dx20,fscal,fjx0);
 756             fjy0             = _mm_macc_pd(dy20,fscal,fjy0);
 757             fjz0             = _mm_macc_pd(dz20,fscal,fjz0);
 758
 759             /**************************
 760              * CALCULATE INTERACTIONS *
 761              **************************/
 762
 763             /* REACTION-FIELD ELECTROSTATICS */
 764             velec            = _mm_mul_pd(qq21,_mm_sub_pd(_mm_macc_pd(krf,rsq21,rinv21),crf));
 765             felec            = _mm_mul_pd(qq21,_mm_msub_pd(rinv21,rinvsq21,krf2));
 766
 767             /* Update potential sum for this i atom from the interaction with this j atom. */
 768             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 769             velecsum         = _mm_add_pd(velecsum,velec);
 770
 771             fscal            = felec;
 772
 773             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 774
 775             /* Update vectorial force */
 776             fix2             = _mm_macc_pd(dx21,fscal,fix2);
 777             fiy2             = _mm_macc_pd(dy21,fscal,fiy2);
 778             fiz2             = _mm_macc_pd(dz21,fscal,fiz2);
 779
 780             fjx1             = _mm_macc_pd(dx21,fscal,fjx1);
 781             fjy1             = _mm_macc_pd(dy21,fscal,fjy1);
 782             fjz1             = _mm_macc_pd(dz21,fscal,fjz1);
 783
 784             /**************************
 785              * CALCULATE INTERACTIONS *
 786              **************************/
 787
 788             /* REACTION-FIELD ELECTROSTATICS */
 789             velec            = _mm_mul_pd(qq22,_mm_sub_pd(_mm_macc_pd(krf,rsq22,rinv22),crf));
 790             felec            = _mm_mul_pd(qq22,_mm_msub_pd(rinv22,rinvsq22,krf2));
 791
 792             /* Update potential sum for this i atom from the interaction with this j atom. */
 793             velec            = _mm_unpacklo_pd(velec,_mm_setzero_pd());
 794             velecsum         = _mm_add_pd(velecsum,velec);
 795
 796             fscal            = felec;
 797
 798             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
 799
 800             /* Update vectorial force */
 801             fix2             = _mm_macc_pd(dx22,fscal,fix2);
 802             fiy2             = _mm_macc_pd(dy22,fscal,fiy2);
 803             fiz2             = _mm_macc_pd(dz22,fscal,fiz2);
 804
 805             fjx2             = _mm_macc_pd(dx22,fscal,fjx2);
 806             fjy2             = _mm_macc_pd(dy22,fscal,fjy2);
 807             fjz2             = _mm_macc_pd(dz22,fscal,fjz2);
 808
 809             gmx_mm_decrement_3rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
 810
 811             /* Inner loop uses 327 flops */
 812         }
 813
 814         /* End of innermost loop */
 815
 816         gmx_mm_update_iforce_3atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,
 817                                               f+i_coord_offset,fshift+i_shift_offset);
 818
 819         ggid                        = gid[iidx];
 820         /* Update potential energies */
 821         gmx_mm_update_1pot_pd(velecsum,kernel_data->energygrp_elec+ggid);
 822         gmx_mm_update_1pot_pd(vvdwsum,kernel_data->energygrp_vdw+ggid);
 823
 824         /* Increment number of inner iterations */
 825         inneriter                  += j_index_end - j_index_start;
 826
 827         /* Outer loop uses 20 flops */
 828     }
 829
 830     /* Increment number of outer iterations */
 831     outeriter        += nri;
 832
 833     /* Update outer/inner flops */
 834
 835     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W3W3_VF,outeriter*20 + inneriter*327);
 836 }
 837 /*
 838  * Gromacs nonbonded kernel:   nb_kernel_ElecRF_VdwLJ_GeomW3W3_F_avx_128_fma_double
 839  * Electrostatics interaction: ReactionField
 840  * VdW interaction:            LennardJones
 841  * Geometry:                   Water3-Water3
 842  * Calculate force/pot:        Force
 843  */
 844 void
 845 nb_kernel_ElecRF_VdwLJ_GeomW3W3_F_avx_128_fma_double
 846                     (t_nblist                    * gmx_restrict       nlist,
 847                      rvec                        * gmx_restrict          xx,
 848                      rvec                        * gmx_restrict          ff,
 849                      t_forcerec                  * gmx_restrict          fr,
 850                      t_mdatoms                   * gmx_restrict     mdatoms,
 851                      nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
 852                      t_nrnb                      * gmx_restrict        nrnb)
 853 {
 854     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
 855      * just 0 for non-waters.
 856      * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
 857      * jnr indices corresponding to data put in the four positions in the SIMD register.
 858      */
 859     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
 860     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
 861     int              jnrA,jnrB;
 862     int              j_coord_offsetA,j_coord_offsetB;
 863     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
 864     real             rcutoff_scalar;
 865     real             *shiftvec,*fshift,*x,*f;
 866     __m128d          tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
 867     int              vdwioffset0;
 868     __m128d          ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
 869     int              vdwioffset1;
 870     __m128d          ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
 871     int              vdwioffset2;
 872     __m128d          ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
 873     int              vdwjidx0A,vdwjidx0B;
 874     __m128d          jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
 875     int              vdwjidx1A,vdwjidx1B;
 876     __m128d          jx1,jy1,jz1,fjx1,fjy1,fjz1,jq1,isaj1;
 877     int              vdwjidx2A,vdwjidx2B;
 878     __m128d          jx2,jy2,jz2,fjx2,fjy2,fjz2,jq2,isaj2;
 879     __m128d          dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
 880     __m128d          dx01,dy01,dz01,rsq01,rinv01,rinvsq01,r01,qq01,c6_01,c12_01;
 881     __m128d          dx02,dy02,dz02,rsq02,rinv02,rinvsq02,r02,qq02,c6_02,c12_02;
 882     __m128d          dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
 883     __m128d          dx11,dy11,dz11,rsq11,rinv11,rinvsq11,r11,qq11,c6_11,c12_11;
 884     __m128d          dx12,dy12,dz12,rsq12,rinv12,rinvsq12,r12,qq12,c6_12,c12_12;
 885     __m128d          dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
 886     __m128d          dx21,dy21,dz21,rsq21,rinv21,rinvsq21,r21,qq21,c6_21,c12_21;
 887     __m128d          dx22,dy22,dz22,rsq22,rinv22,rinvsq22,r22,qq22,c6_22,c12_22;
 888     __m128d          velec,felec,velecsum,facel,crf,krf,krf2;
 889     real             *charge;
 890     int              nvdwtype;
 891     __m128d          rinvsix,rvdw,vvdw,vvdw6,vvdw12,fvdw,fvdw6,fvdw12,vvdwsum,sh_vdw_invrcut6;
 892     int              *vdwtype;
 893     real             *vdwparam;
 894     __m128d          one_sixth   = _mm_set1_pd(1.0/6.0);
 895     __m128d          one_twelfth = _mm_set1_pd(1.0/12.0);
 896     __m128d          dummy_mask,cutoff_mask;
 897     __m128d          signbit   = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
 898     __m128d          one     = _mm_set1_pd(1.0);
 899     __m128d          two     = _mm_set1_pd(2.0);
 900     x                = xx[0];
 901     f                = ff[0];
 902
 903     nri              = nlist->nri;
 904     iinr             = nlist->iinr;
 905     jindex           = nlist->jindex;
 906     jjnr             = nlist->jjnr;
 907     shiftidx         = nlist->shift;
 908     gid              = nlist->gid;
 909     shiftvec         = fr->shift_vec[0];
 910     fshift           = fr->fshift[0];
 911     facel            = _mm_set1_pd(fr->epsfac);
 912     charge           = mdatoms->chargeA;
 913     krf              = _mm_set1_pd(fr->ic->k_rf);
 914     krf2             = _mm_set1_pd(fr->ic->k_rf*2.0);
 915     crf              = _mm_set1_pd(fr->ic->c_rf);
 916     nvdwtype         = fr->ntype;
 917     vdwparam         = fr->nbfp;
 918     vdwtype          = mdatoms->typeA;
 919
 920     /* Setup water-specific parameters */
 921     inr              = nlist->iinr[0];
 922     iq0              = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+0]));
 923     iq1              = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
 924     iq2              = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
 925     vdwioffset0      = 2*nvdwtype*vdwtype[inr+0];
 926
 927     jq0              = _mm_set1_pd(charge[inr+0]);
 928     jq1              = _mm_set1_pd(charge[inr+1]);
 929     jq2              = _mm_set1_pd(charge[inr+2]);
 930     vdwjidx0A        = 2*vdwtype[inr+0];
 931     qq00             = _mm_mul_pd(iq0,jq0);
 932     c6_00            = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A]);
 933     c12_00           = _mm_set1_pd(vdwparam[vdwioffset0+vdwjidx0A+1]);
 934     qq01             = _mm_mul_pd(iq0,jq1);
 935     qq02             = _mm_mul_pd(iq0,jq2);
 936     qq10             = _mm_mul_pd(iq1,jq0);
 937     qq11             = _mm_mul_pd(iq1,jq1);
 938     qq12             = _mm_mul_pd(iq1,jq2);
 939     qq20             = _mm_mul_pd(iq2,jq0);
 940     qq21             = _mm_mul_pd(iq2,jq1);
 941     qq22             = _mm_mul_pd(iq2,jq2);
 942
 943     /* Avoid stupid compiler warnings */
 944     jnrA = jnrB = 0;
 945     j_coord_offsetA = 0;
 946     j_coord_offsetB = 0;
 947
 948     outeriter        = 0;
 949     inneriter        = 0;
 950
 951     /* Start outer loop over neighborlists */
 952     for(iidx=0; iidx<nri; iidx++)
 953     {
 954         /* Load shift vector for this list */
 955         i_shift_offset   = DIM*shiftidx[iidx];
 956
 957         /* Load limits for loop over neighbors */
 958         j_index_start    = jindex[iidx];
 959         j_index_end      = jindex[iidx+1];
 960
 961         /* Get outer coordinate index */
 962         inr              = iinr[iidx];
 963         i_coord_offset   = DIM*inr;
 964
 965         /* Load i particle coords and add shift vector */
 966         gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset,
 967                                                  &ix0,&iy0,&iz0,&ix1,&iy1,&iz1,&ix2,&iy2,&iz2);
 968
 969         fix0             = _mm_setzero_pd();
 970         fiy0             = _mm_setzero_pd();
 971         fiz0             = _mm_setzero_pd();
 972         fix1             = _mm_setzero_pd();
 973         fiy1             = _mm_setzero_pd();
 974         fiz1             = _mm_setzero_pd();
 975         fix2             = _mm_setzero_pd();
 976         fiy2             = _mm_setzero_pd();
 977         fiz2             = _mm_setzero_pd();
 978
 979         /* Start inner kernel loop */
 980         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
 981         {
 982
 983             /* Get j neighbor index, and coordinate index */
 984             jnrA             = jjnr[jidx];
 985             jnrB             = jjnr[jidx+1];
 986             j_coord_offsetA  = DIM*jnrA;
 987             j_coord_offsetB  = DIM*jnrB;
 988
 989             /* load j atom coordinates */
 990             gmx_mm_load_3rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
 991                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
 992
 993             /* Calculate displacement vector */
 994             dx00             = _mm_sub_pd(ix0,jx0);
 995             dy00             = _mm_sub_pd(iy0,jy0);
 996             dz00             = _mm_sub_pd(iz0,jz0);
 997             dx01             = _mm_sub_pd(ix0,jx1);
 998             dy01             = _mm_sub_pd(iy0,jy1);
 999             dz01             = _mm_sub_pd(iz0,jz1);
1000             dx02             = _mm_sub_pd(ix0,jx2);
1001             dy02             = _mm_sub_pd(iy0,jy2);
1002             dz02             = _mm_sub_pd(iz0,jz2);
1003             dx10             = _mm_sub_pd(ix1,jx0);
1004             dy10             = _mm_sub_pd(iy1,jy0);
1005             dz10             = _mm_sub_pd(iz1,jz0);
1006             dx11             = _mm_sub_pd(ix1,jx1);
1007             dy11             = _mm_sub_pd(iy1,jy1);
1008             dz11             = _mm_sub_pd(iz1,jz1);
1009             dx12             = _mm_sub_pd(ix1,jx2);
1010             dy12             = _mm_sub_pd(iy1,jy2);
1011             dz12             = _mm_sub_pd(iz1,jz2);
1012             dx20             = _mm_sub_pd(ix2,jx0);
1013             dy20             = _mm_sub_pd(iy2,jy0);
1014             dz20             = _mm_sub_pd(iz2,jz0);
1015             dx21             = _mm_sub_pd(ix2,jx1);
1016             dy21             = _mm_sub_pd(iy2,jy1);
1017             dz21             = _mm_sub_pd(iz2,jz1);
1018             dx22             = _mm_sub_pd(ix2,jx2);
1019             dy22             = _mm_sub_pd(iy2,jy2);
1020             dz22             = _mm_sub_pd(iz2,jz2);
1021
1022             /* Calculate squared distance and things based on it */
1023             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
1024             rsq01            = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
1025             rsq02            = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
1026             rsq10            = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
1027             rsq11            = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
1028             rsq12            = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
1029             rsq20            = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
1030             rsq21            = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
1031             rsq22            = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
1032
1033             rinv00           = gmx_mm_invsqrt_pd(rsq00);
1034             rinv01           = gmx_mm_invsqrt_pd(rsq01);
1035             rinv02           = gmx_mm_invsqrt_pd(rsq02);
1036             rinv10           = gmx_mm_invsqrt_pd(rsq10);
1037             rinv11           = gmx_mm_invsqrt_pd(rsq11);
1038             rinv12           = gmx_mm_invsqrt_pd(rsq12);
1039             rinv20           = gmx_mm_invsqrt_pd(rsq20);
1040             rinv21           = gmx_mm_invsqrt_pd(rsq21);
1041             rinv22           = gmx_mm_invsqrt_pd(rsq22);
1042
1043             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
1044             rinvsq01         = _mm_mul_pd(rinv01,rinv01);
1045             rinvsq02         = _mm_mul_pd(rinv02,rinv02);
1046             rinvsq10         = _mm_mul_pd(rinv10,rinv10);
1047             rinvsq11         = _mm_mul_pd(rinv11,rinv11);
1048             rinvsq12         = _mm_mul_pd(rinv12,rinv12);
1049             rinvsq20         = _mm_mul_pd(rinv20,rinv20);
1050             rinvsq21         = _mm_mul_pd(rinv21,rinv21);
1051             rinvsq22         = _mm_mul_pd(rinv22,rinv22);
1052
1053             fjx0             = _mm_setzero_pd();
1054             fjy0             = _mm_setzero_pd();
1055             fjz0             = _mm_setzero_pd();
1056             fjx1             = _mm_setzero_pd();
1057             fjy1             = _mm_setzero_pd();
1058             fjz1             = _mm_setzero_pd();
1059             fjx2             = _mm_setzero_pd();
1060             fjy2             = _mm_setzero_pd();
1061             fjz2             = _mm_setzero_pd();
1062
1063             /**************************
1064              * CALCULATE INTERACTIONS *
1065              **************************/
1066
1067             /* REACTION-FIELD ELECTROSTATICS */
1068             felec            = _mm_mul_pd(qq00,_mm_msub_pd(rinv00,rinvsq00,krf2));
1069
1070             /* LENNARD-JONES DISPERSION/REPULSION */
1071
1072             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
1073             fvdw             = _mm_mul_pd(_mm_msub_pd(c12_00,rinvsix,c6_00),_mm_mul_pd(rinvsix,rinvsq00));
1074
1075             fscal            = _mm_add_pd(felec,fvdw);
1076
1077             /* Update vectorial force */
1078             fix0             = _mm_macc_pd(dx00,fscal,fix0);
1079             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
1080             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
1081
1082             fjx0             = _mm_macc_pd(dx00,fscal,fjx0);
1083             fjy0             = _mm_macc_pd(dy00,fscal,fjy0);
1084             fjz0             = _mm_macc_pd(dz00,fscal,fjz0);
1085
1086             /**************************
1087              * CALCULATE INTERACTIONS *
1088              **************************/
1089
1090             /* REACTION-FIELD ELECTROSTATICS */
1091             felec            = _mm_mul_pd(qq01,_mm_msub_pd(rinv01,rinvsq01,krf2));
1092
1093             fscal            = felec;
1094
1095             /* Update vectorial force */
1096             fix0             = _mm_macc_pd(dx01,fscal,fix0);
1097             fiy0             = _mm_macc_pd(dy01,fscal,fiy0);
1098             fiz0             = _mm_macc_pd(dz01,fscal,fiz0);
1099
1100             fjx1             = _mm_macc_pd(dx01,fscal,fjx1);
1101             fjy1             = _mm_macc_pd(dy01,fscal,fjy1);
1102             fjz1             = _mm_macc_pd(dz01,fscal,fjz1);
1103
1104             /**************************
1105              * CALCULATE INTERACTIONS *
1106              **************************/
1107
1108             /* REACTION-FIELD ELECTROSTATICS */
1109             felec            = _mm_mul_pd(qq02,_mm_msub_pd(rinv02,rinvsq02,krf2));
1110
1111             fscal            = felec;
1112
1113             /* Update vectorial force */
1114             fix0             = _mm_macc_pd(dx02,fscal,fix0);
1115             fiy0             = _mm_macc_pd(dy02,fscal,fiy0);
1116             fiz0             = _mm_macc_pd(dz02,fscal,fiz0);
1117
1118             fjx2             = _mm_macc_pd(dx02,fscal,fjx2);
1119             fjy2             = _mm_macc_pd(dy02,fscal,fjy2);
1120             fjz2             = _mm_macc_pd(dz02,fscal,fjz2);
1121
1122             /**************************
1123              * CALCULATE INTERACTIONS *
1124              **************************/
1125
1126             /* REACTION-FIELD ELECTROSTATICS */
1127             felec            = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
1128
1129             fscal            = felec;
1130
1131             /* Update vectorial force */
1132             fix1             = _mm_macc_pd(dx10,fscal,fix1);
1133             fiy1             = _mm_macc_pd(dy10,fscal,fiy1);
1134             fiz1             = _mm_macc_pd(dz10,fscal,fiz1);
1135
1136             fjx0             = _mm_macc_pd(dx10,fscal,fjx0);
1137             fjy0             = _mm_macc_pd(dy10,fscal,fjy0);
1138             fjz0             = _mm_macc_pd(dz10,fscal,fjz0);
1139
1140             /**************************
1141              * CALCULATE INTERACTIONS *
1142              **************************/
1143
1144             /* REACTION-FIELD ELECTROSTATICS */
1145             felec            = _mm_mul_pd(qq11,_mm_msub_pd(rinv11,rinvsq11,krf2));
1146
1147             fscal            = felec;
1148
1149             /* Update vectorial force */
1150             fix1             = _mm_macc_pd(dx11,fscal,fix1);
1151             fiy1             = _mm_macc_pd(dy11,fscal,fiy1);
1152             fiz1             = _mm_macc_pd(dz11,fscal,fiz1);
1153
1154             fjx1             = _mm_macc_pd(dx11,fscal,fjx1);
1155             fjy1             = _mm_macc_pd(dy11,fscal,fjy1);
1156             fjz1             = _mm_macc_pd(dz11,fscal,fjz1);
1157
1158             /**************************
1159              * CALCULATE INTERACTIONS *
1160              **************************/
1161
1162             /* REACTION-FIELD ELECTROSTATICS */
1163             felec            = _mm_mul_pd(qq12,_mm_msub_pd(rinv12,rinvsq12,krf2));
1164
1165             fscal            = felec;
1166
1167             /* Update vectorial force */
1168             fix1             = _mm_macc_pd(dx12,fscal,fix1);
1169             fiy1             = _mm_macc_pd(dy12,fscal,fiy1);
1170             fiz1             = _mm_macc_pd(dz12,fscal,fiz1);
1171
1172             fjx2             = _mm_macc_pd(dx12,fscal,fjx2);
1173             fjy2             = _mm_macc_pd(dy12,fscal,fjy2);
1174             fjz2             = _mm_macc_pd(dz12,fscal,fjz2);
1175
1176             /**************************
1177              * CALCULATE INTERACTIONS *
1178              **************************/
1179
1180             /* REACTION-FIELD ELECTROSTATICS */
1181             felec            = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
1182
1183             fscal            = felec;
1184
1185             /* Update vectorial force */
1186             fix2             = _mm_macc_pd(dx20,fscal,fix2);
1187             fiy2             = _mm_macc_pd(dy20,fscal,fiy2);
1188             fiz2             = _mm_macc_pd(dz20,fscal,fiz2);
1189
1190             fjx0             = _mm_macc_pd(dx20,fscal,fjx0);
1191             fjy0             = _mm_macc_pd(dy20,fscal,fjy0);
1192             fjz0             = _mm_macc_pd(dz20,fscal,fjz0);
1193
1194             /**************************
1195              * CALCULATE INTERACTIONS *
1196              **************************/
1197
1198             /* REACTION-FIELD ELECTROSTATICS */
1199             felec            = _mm_mul_pd(qq21,_mm_msub_pd(rinv21,rinvsq21,krf2));
1200
1201             fscal            = felec;
1202
1203             /* Update vectorial force */
1204             fix2             = _mm_macc_pd(dx21,fscal,fix2);
1205             fiy2             = _mm_macc_pd(dy21,fscal,fiy2);
1206             fiz2             = _mm_macc_pd(dz21,fscal,fiz2);
1207
1208             fjx1             = _mm_macc_pd(dx21,fscal,fjx1);
1209             fjy1             = _mm_macc_pd(dy21,fscal,fjy1);
1210             fjz1             = _mm_macc_pd(dz21,fscal,fjz1);
1211
1212             /**************************
1213              * CALCULATE INTERACTIONS *
1214              **************************/
1215
1216             /* REACTION-FIELD ELECTROSTATICS */
1217             felec            = _mm_mul_pd(qq22,_mm_msub_pd(rinv22,rinvsq22,krf2));
1218
1219             fscal            = felec;
1220
1221             /* Update vectorial force */
1222             fix2             = _mm_macc_pd(dx22,fscal,fix2);
1223             fiy2             = _mm_macc_pd(dy22,fscal,fiy2);
1224             fiz2             = _mm_macc_pd(dz22,fscal,fiz2);
1225
1226             fjx2             = _mm_macc_pd(dx22,fscal,fjx2);
1227             fjy2             = _mm_macc_pd(dy22,fscal,fjy2);
1228             fjz2             = _mm_macc_pd(dz22,fscal,fjz2);
1229
1230             gmx_mm_decrement_3rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
1231
1232             /* Inner loop uses 277 flops */
1233         }
1234
1235         if(jidx<j_index_end)
1236         {
1237
1238             jnrA             = jjnr[jidx];
1239             j_coord_offsetA  = DIM*jnrA;
1240
1241             /* load j atom coordinates */
1242             gmx_mm_load_3rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
1243                                               &jx0,&jy0,&jz0,&jx1,&jy1,&jz1,&jx2,&jy2,&jz2);
1244
1245             /* Calculate displacement vector */
1246             dx00             = _mm_sub_pd(ix0,jx0);
1247             dy00             = _mm_sub_pd(iy0,jy0);
1248             dz00             = _mm_sub_pd(iz0,jz0);
1249             dx01             = _mm_sub_pd(ix0,jx1);
1250             dy01             = _mm_sub_pd(iy0,jy1);
1251             dz01             = _mm_sub_pd(iz0,jz1);
1252             dx02             = _mm_sub_pd(ix0,jx2);
1253             dy02             = _mm_sub_pd(iy0,jy2);
1254             dz02             = _mm_sub_pd(iz0,jz2);
1255             dx10             = _mm_sub_pd(ix1,jx0);
1256             dy10             = _mm_sub_pd(iy1,jy0);
1257             dz10             = _mm_sub_pd(iz1,jz0);
1258             dx11             = _mm_sub_pd(ix1,jx1);
1259             dy11             = _mm_sub_pd(iy1,jy1);
1260             dz11             = _mm_sub_pd(iz1,jz1);
1261             dx12             = _mm_sub_pd(ix1,jx2);
1262             dy12             = _mm_sub_pd(iy1,jy2);
1263             dz12             = _mm_sub_pd(iz1,jz2);
1264             dx20             = _mm_sub_pd(ix2,jx0);
1265             dy20             = _mm_sub_pd(iy2,jy0);
1266             dz20             = _mm_sub_pd(iz2,jz0);
1267             dx21             = _mm_sub_pd(ix2,jx1);
1268             dy21             = _mm_sub_pd(iy2,jy1);
1269             dz21             = _mm_sub_pd(iz2,jz1);
1270             dx22             = _mm_sub_pd(ix2,jx2);
1271             dy22             = _mm_sub_pd(iy2,jy2);
1272             dz22             = _mm_sub_pd(iz2,jz2);
1273
1274             /* Calculate squared distance and things based on it */
1275             rsq00            = gmx_mm_calc_rsq_pd(dx00,dy00,dz00);
1276             rsq01            = gmx_mm_calc_rsq_pd(dx01,dy01,dz01);
1277             rsq02            = gmx_mm_calc_rsq_pd(dx02,dy02,dz02);
1278             rsq10            = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
1279             rsq11            = gmx_mm_calc_rsq_pd(dx11,dy11,dz11);
1280             rsq12            = gmx_mm_calc_rsq_pd(dx12,dy12,dz12);
1281             rsq20            = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
1282             rsq21            = gmx_mm_calc_rsq_pd(dx21,dy21,dz21);
1283             rsq22            = gmx_mm_calc_rsq_pd(dx22,dy22,dz22);
1284
1285             rinv00           = gmx_mm_invsqrt_pd(rsq00);
1286             rinv01           = gmx_mm_invsqrt_pd(rsq01);
1287             rinv02           = gmx_mm_invsqrt_pd(rsq02);
1288             rinv10           = gmx_mm_invsqrt_pd(rsq10);
1289             rinv11           = gmx_mm_invsqrt_pd(rsq11);
1290             rinv12           = gmx_mm_invsqrt_pd(rsq12);
1291             rinv20           = gmx_mm_invsqrt_pd(rsq20);
1292             rinv21           = gmx_mm_invsqrt_pd(rsq21);
1293             rinv22           = gmx_mm_invsqrt_pd(rsq22);
1294
1295             rinvsq00         = _mm_mul_pd(rinv00,rinv00);
1296             rinvsq01         = _mm_mul_pd(rinv01,rinv01);
1297             rinvsq02         = _mm_mul_pd(rinv02,rinv02);
1298             rinvsq10         = _mm_mul_pd(rinv10,rinv10);
1299             rinvsq11         = _mm_mul_pd(rinv11,rinv11);
1300             rinvsq12         = _mm_mul_pd(rinv12,rinv12);
1301             rinvsq20         = _mm_mul_pd(rinv20,rinv20);
1302             rinvsq21         = _mm_mul_pd(rinv21,rinv21);
1303             rinvsq22         = _mm_mul_pd(rinv22,rinv22);
1304
1305             fjx0             = _mm_setzero_pd();
1306             fjy0             = _mm_setzero_pd();
1307             fjz0             = _mm_setzero_pd();
1308             fjx1             = _mm_setzero_pd();
1309             fjy1             = _mm_setzero_pd();
1310             fjz1             = _mm_setzero_pd();
1311             fjx2             = _mm_setzero_pd();
1312             fjy2             = _mm_setzero_pd();
1313             fjz2             = _mm_setzero_pd();
1314
1315             /**************************
1316              * CALCULATE INTERACTIONS *
1317              **************************/
1318
1319             /* REACTION-FIELD ELECTROSTATICS */
1320             felec            = _mm_mul_pd(qq00,_mm_msub_pd(rinv00,rinvsq00,krf2));
1321
1322             /* LENNARD-JONES DISPERSION/REPULSION */
1323
1324             rinvsix          = _mm_mul_pd(_mm_mul_pd(rinvsq00,rinvsq00),rinvsq00);
1325             fvdw             = _mm_mul_pd(_mm_msub_pd(c12_00,rinvsix,c6_00),_mm_mul_pd(rinvsix,rinvsq00));
1326
1327             fscal            = _mm_add_pd(felec,fvdw);
1328
1329             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1330
1331             /* Update vectorial force */
1332             fix0             = _mm_macc_pd(dx00,fscal,fix0);
1333             fiy0             = _mm_macc_pd(dy00,fscal,fiy0);
1334             fiz0             = _mm_macc_pd(dz00,fscal,fiz0);
1335
1336             fjx0             = _mm_macc_pd(dx00,fscal,fjx0);
1337             fjy0             = _mm_macc_pd(dy00,fscal,fjy0);
1338             fjz0             = _mm_macc_pd(dz00,fscal,fjz0);
1339
1340             /**************************
1341              * CALCULATE INTERACTIONS *
1342              **************************/
1343
1344             /* REACTION-FIELD ELECTROSTATICS */
1345             felec            = _mm_mul_pd(qq01,_mm_msub_pd(rinv01,rinvsq01,krf2));
1346
1347             fscal            = felec;
1348
1349             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1350
1351             /* Update vectorial force */
1352             fix0             = _mm_macc_pd(dx01,fscal,fix0);
1353             fiy0             = _mm_macc_pd(dy01,fscal,fiy0);
1354             fiz0             = _mm_macc_pd(dz01,fscal,fiz0);
1355
1356             fjx1             = _mm_macc_pd(dx01,fscal,fjx1);
1357             fjy1             = _mm_macc_pd(dy01,fscal,fjy1);
1358             fjz1             = _mm_macc_pd(dz01,fscal,fjz1);
1359
1360             /**************************
1361              * CALCULATE INTERACTIONS *
1362              **************************/
1363
1364             /* REACTION-FIELD ELECTROSTATICS */
1365             felec            = _mm_mul_pd(qq02,_mm_msub_pd(rinv02,rinvsq02,krf2));
1366
1367             fscal            = felec;
1368
1369             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1370
1371             /* Update vectorial force */
1372             fix0             = _mm_macc_pd(dx02,fscal,fix0);
1373             fiy0             = _mm_macc_pd(dy02,fscal,fiy0);
1374             fiz0             = _mm_macc_pd(dz02,fscal,fiz0);
1375
1376             fjx2             = _mm_macc_pd(dx02,fscal,fjx2);
1377             fjy2             = _mm_macc_pd(dy02,fscal,fjy2);
1378             fjz2             = _mm_macc_pd(dz02,fscal,fjz2);
1379
1380             /**************************
1381              * CALCULATE INTERACTIONS *
1382              **************************/
1383
1384             /* REACTION-FIELD ELECTROSTATICS */
1385             felec            = _mm_mul_pd(qq10,_mm_msub_pd(rinv10,rinvsq10,krf2));
1386
1387             fscal            = felec;
1388
1389             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1390
1391             /* Update vectorial force */
1392             fix1             = _mm_macc_pd(dx10,fscal,fix1);
1393             fiy1             = _mm_macc_pd(dy10,fscal,fiy1);
1394             fiz1             = _mm_macc_pd(dz10,fscal,fiz1);
1395
1396             fjx0             = _mm_macc_pd(dx10,fscal,fjx0);
1397             fjy0             = _mm_macc_pd(dy10,fscal,fjy0);
1398             fjz0             = _mm_macc_pd(dz10,fscal,fjz0);
1399
1400             /**************************
1401              * CALCULATE INTERACTIONS *
1402              **************************/
1403
1404             /* REACTION-FIELD ELECTROSTATICS */
1405             felec            = _mm_mul_pd(qq11,_mm_msub_pd(rinv11,rinvsq11,krf2));
1406
1407             fscal            = felec;
1408
1409             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1410
1411             /* Update vectorial force */
1412             fix1             = _mm_macc_pd(dx11,fscal,fix1);
1413             fiy1             = _mm_macc_pd(dy11,fscal,fiy1);
1414             fiz1             = _mm_macc_pd(dz11,fscal,fiz1);
1415
1416             fjx1             = _mm_macc_pd(dx11,fscal,fjx1);
1417             fjy1             = _mm_macc_pd(dy11,fscal,fjy1);
1418             fjz1             = _mm_macc_pd(dz11,fscal,fjz1);
1419
1420             /**************************
1421              * CALCULATE INTERACTIONS *
1422              **************************/
1423
1424             /* REACTION-FIELD ELECTROSTATICS */
1425             felec            = _mm_mul_pd(qq12,_mm_msub_pd(rinv12,rinvsq12,krf2));
1426
1427             fscal            = felec;
1428
1429             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1430
1431             /* Update vectorial force */
1432             fix1             = _mm_macc_pd(dx12,fscal,fix1);
1433             fiy1             = _mm_macc_pd(dy12,fscal,fiy1);
1434             fiz1             = _mm_macc_pd(dz12,fscal,fiz1);
1435
1436             fjx2             = _mm_macc_pd(dx12,fscal,fjx2);
1437             fjy2             = _mm_macc_pd(dy12,fscal,fjy2);
1438             fjz2             = _mm_macc_pd(dz12,fscal,fjz2);
1439
1440             /**************************
1441              * CALCULATE INTERACTIONS *
1442              **************************/
1443
1444             /* REACTION-FIELD ELECTROSTATICS */
1445             felec            = _mm_mul_pd(qq20,_mm_msub_pd(rinv20,rinvsq20,krf2));
1446
1447             fscal            = felec;
1448
1449             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1450
1451             /* Update vectorial force */
1452             fix2             = _mm_macc_pd(dx20,fscal,fix2);
1453             fiy2             = _mm_macc_pd(dy20,fscal,fiy2);
1454             fiz2             = _mm_macc_pd(dz20,fscal,fiz2);
1455
1456             fjx0             = _mm_macc_pd(dx20,fscal,fjx0);
1457             fjy0             = _mm_macc_pd(dy20,fscal,fjy0);
1458             fjz0             = _mm_macc_pd(dz20,fscal,fjz0);
1459
1460             /**************************
1461              * CALCULATE INTERACTIONS *
1462              **************************/
1463
1464             /* REACTION-FIELD ELECTROSTATICS */
1465             felec            = _mm_mul_pd(qq21,_mm_msub_pd(rinv21,rinvsq21,krf2));
1466
1467             fscal            = felec;
1468
1469             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1470
1471             /* Update vectorial force */
1472             fix2             = _mm_macc_pd(dx21,fscal,fix2);
1473             fiy2             = _mm_macc_pd(dy21,fscal,fiy2);
1474             fiz2             = _mm_macc_pd(dz21,fscal,fiz2);
1475
1476             fjx1             = _mm_macc_pd(dx21,fscal,fjx1);
1477             fjy1             = _mm_macc_pd(dy21,fscal,fjy1);
1478             fjz1             = _mm_macc_pd(dz21,fscal,fjz1);
1479
1480             /**************************
1481              * CALCULATE INTERACTIONS *
1482              **************************/
1483
1484             /* REACTION-FIELD ELECTROSTATICS */
1485             felec            = _mm_mul_pd(qq22,_mm_msub_pd(rinv22,rinvsq22,krf2));
1486
1487             fscal            = felec;
1488
1489             fscal            = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
1490
1491             /* Update vectorial force */
1492             fix2             = _mm_macc_pd(dx22,fscal,fix2);
1493             fiy2             = _mm_macc_pd(dy22,fscal,fiy2);
1494             fiz2             = _mm_macc_pd(dz22,fscal,fiz2);
1495
1496             fjx2             = _mm_macc_pd(dx22,fscal,fjx2);
1497             fjy2             = _mm_macc_pd(dy22,fscal,fjy2);
1498             fjz2             = _mm_macc_pd(dz22,fscal,fjz2);
1499
1500             gmx_mm_decrement_3rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0,fjx1,fjy1,fjz1,fjx2,fjy2,fjz2);
1501
1502             /* Inner loop uses 277 flops */
1503         }
1504
1505         /* End of innermost loop */
1506
1507         gmx_mm_update_iforce_3atom_swizzle_pd(fix0,fiy0,fiz0,fix1,fiy1,fiz1,fix2,fiy2,fiz2,
1508                                               f+i_coord_offset,fshift+i_shift_offset);
1509
1510         /* Increment number of inner iterations */
1511         inneriter                  += j_index_end - j_index_start;
1512
1513         /* Outer loop uses 18 flops */
1514     }
1515
1516     /* Increment number of outer iterations */
1517     outeriter        += nri;
1518
1519     /* Update outer/inner flops */
1520
1521     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VDW_W3W3_F,outeriter*18 + inneriter*277);
1522 }