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

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