K-computer specific modifications
[alexxy/gromacs.git] / src / gromacs / gmxlib / nonbonded / nb_kernel_sparc64_hpc_ace_double / nb_kernel_ElecEwSw_VdwNone_GeomW4P1_sparc64_hpc_ace_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 sparc64_hpc_ace_double kernel generator.
37  */
38 #ifdef HAVE_CONFIG_H
39 #include <config.h>
40 #endif
41
42 #include <math.h>
43
44 #include "../nb_kernel.h"
45 #include "types/simple.h"
46 #include "gromacs/legacyheaders/vec.h"
47 #include "nrnb.h"
48
49 #include "kernelutil_sparc64_hpc_ace_double.h"
50
51 /*
52  * Gromacs nonbonded kernel:   nb_kernel_ElecEwSw_VdwNone_GeomW4P1_VF_sparc64_hpc_ace_double
53  * Electrostatics interaction: Ewald
54  * VdW interaction:            None
55  * Geometry:                   Water4-Particle
56  * Calculate force/pot:        PotentialAndForce
57  */
58 void
59 nb_kernel_ElecEwSw_VdwNone_GeomW4P1_VF_sparc64_hpc_ace_double
60                     (t_nblist                    * gmx_restrict       nlist,
61                      rvec                        * gmx_restrict          xx,
62                      rvec                        * gmx_restrict          ff,
63                      t_forcerec                  * gmx_restrict          fr,
64                      t_mdatoms                   * gmx_restrict     mdatoms,
65                      nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
66                      t_nrnb                      * gmx_restrict        nrnb)
67 {
68     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
69      * just 0 for non-waters.
70      * Suffixes A,B refer to j loop unrolling done with double precision SIMD, e.g. for the two different
71      * jnr indices corresponding to data put in the four positions in the SIMD register.
72      */
73     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
74     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
75     int              jnrA,jnrB;
76     int              j_coord_offsetA,j_coord_offsetB;
77     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
78     real             rcutoff_scalar;
79     real             *shiftvec,*fshift,*x,*f;
80     _fjsp_v2r8       tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
81     int              vdwioffset1;
82     _fjsp_v2r8       ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
83     int              vdwioffset2;
84     _fjsp_v2r8       ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
85     int              vdwioffset3;
86     _fjsp_v2r8       ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
87     int              vdwjidx0A,vdwjidx0B;
88     _fjsp_v2r8       jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
89     _fjsp_v2r8       dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
90     _fjsp_v2r8       dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
91     _fjsp_v2r8       dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
92     _fjsp_v2r8       velec,felec,velecsum,facel,crf,krf,krf2;
93     real             *charge;
94     _fjsp_v2r8       ewtabscale,eweps,sh_ewald,ewrt,ewtabhalfspace,ewtabF,ewtabFn,ewtabD,ewtabV;
95     real             *ewtab;
96     _fjsp_v2r8       rswitch,swV3,swV4,swV5,swF2,swF3,swF4,d,d2,sw,dsw;
97     real             rswitch_scalar,d_scalar;
98     _fjsp_v2r8       itab_tmp;
99     _fjsp_v2r8       dummy_mask,cutoff_mask;
100     _fjsp_v2r8       one     = gmx_fjsp_set1_v2r8(1.0);
101     _fjsp_v2r8       two     = gmx_fjsp_set1_v2r8(2.0);
102     union { _fjsp_v2r8 simd; long long int i[2]; } vfconv,gbconv,ewconv;
103
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            = gmx_fjsp_set1_v2r8(fr->epsfac);
116     charge           = mdatoms->chargeA;
117
118     sh_ewald         = gmx_fjsp_set1_v2r8(fr->ic->sh_ewald);
119     ewtab            = fr->ic->tabq_coul_FDV0;
120     ewtabscale       = gmx_fjsp_set1_v2r8(fr->ic->tabq_scale);
121     ewtabhalfspace   = gmx_fjsp_set1_v2r8(0.5/fr->ic->tabq_scale);
122
123     /* Setup water-specific parameters */
124     inr              = nlist->iinr[0];
125     iq1              = _fjsp_mul_v2r8(facel,gmx_fjsp_set1_v2r8(charge[inr+1]));
126     iq2              = _fjsp_mul_v2r8(facel,gmx_fjsp_set1_v2r8(charge[inr+2]));
127     iq3              = _fjsp_mul_v2r8(facel,gmx_fjsp_set1_v2r8(charge[inr+3]));
128
129     /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
130     rcutoff_scalar   = fr->rcoulomb;
131     rcutoff          = gmx_fjsp_set1_v2r8(rcutoff_scalar);
132     rcutoff2         = _fjsp_mul_v2r8(rcutoff,rcutoff);
133
134     rswitch_scalar   = fr->rcoulomb_switch;
135     rswitch          = gmx_fjsp_set1_v2r8(rswitch_scalar);
136     /* Setup switch parameters */
137     d_scalar         = rcutoff_scalar-rswitch_scalar;
138     d                = gmx_fjsp_set1_v2r8(d_scalar);
139     swV3             = gmx_fjsp_set1_v2r8(-10.0/(d_scalar*d_scalar*d_scalar));
140     swV4             = gmx_fjsp_set1_v2r8( 15.0/(d_scalar*d_scalar*d_scalar*d_scalar));
141     swV5             = gmx_fjsp_set1_v2r8( -6.0/(d_scalar*d_scalar*d_scalar*d_scalar*d_scalar));
142     swF2             = gmx_fjsp_set1_v2r8(-30.0/(d_scalar*d_scalar*d_scalar));
143     swF3             = gmx_fjsp_set1_v2r8( 60.0/(d_scalar*d_scalar*d_scalar*d_scalar));
144     swF4             = gmx_fjsp_set1_v2r8(-30.0/(d_scalar*d_scalar*d_scalar*d_scalar*d_scalar));
145
146     /* Avoid stupid compiler warnings */
147     jnrA = jnrB = 0;
148     j_coord_offsetA = 0;
149     j_coord_offsetB = 0;
150
151     outeriter        = 0;
152     inneriter        = 0;
153
154     /* Start outer loop over neighborlists */
155     for(iidx=0; iidx<nri; iidx++)
156     {
157         /* Load shift vector for this list */
158         i_shift_offset   = DIM*shiftidx[iidx];
159
160         /* Load limits for loop over neighbors */
161         j_index_start    = jindex[iidx];
162         j_index_end      = jindex[iidx+1];
163
164         /* Get outer coordinate index */
165         inr              = iinr[iidx];
166         i_coord_offset   = DIM*inr;
167
168         /* Load i particle coords and add shift vector */
169         gmx_fjsp_load_shift_and_3rvec_broadcast_v2r8(shiftvec+i_shift_offset,x+i_coord_offset+DIM,
170                                                  &ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
171
172         fix1             = _fjsp_setzero_v2r8();
173         fiy1             = _fjsp_setzero_v2r8();
174         fiz1             = _fjsp_setzero_v2r8();
175         fix2             = _fjsp_setzero_v2r8();
176         fiy2             = _fjsp_setzero_v2r8();
177         fiz2             = _fjsp_setzero_v2r8();
178         fix3             = _fjsp_setzero_v2r8();
179         fiy3             = _fjsp_setzero_v2r8();
180         fiz3             = _fjsp_setzero_v2r8();
181
182         /* Reset potential sums */
183         velecsum         = _fjsp_setzero_v2r8();
184
185         /* Start inner kernel loop */
186         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
187         {
188
189             /* Get j neighbor index, and coordinate index */
190             jnrA             = jjnr[jidx];
191             jnrB             = jjnr[jidx+1];
192             j_coord_offsetA  = DIM*jnrA;
193             j_coord_offsetB  = DIM*jnrB;
194
195             /* load j atom coordinates */
196             gmx_fjsp_load_1rvec_2ptr_swizzle_v2r8(x+j_coord_offsetA,x+j_coord_offsetB,
197                                               &jx0,&jy0,&jz0);
198
199             /* Calculate displacement vector */
200             dx10             = _fjsp_sub_v2r8(ix1,jx0);
201             dy10             = _fjsp_sub_v2r8(iy1,jy0);
202             dz10             = _fjsp_sub_v2r8(iz1,jz0);
203             dx20             = _fjsp_sub_v2r8(ix2,jx0);
204             dy20             = _fjsp_sub_v2r8(iy2,jy0);
205             dz20             = _fjsp_sub_v2r8(iz2,jz0);
206             dx30             = _fjsp_sub_v2r8(ix3,jx0);
207             dy30             = _fjsp_sub_v2r8(iy3,jy0);
208             dz30             = _fjsp_sub_v2r8(iz3,jz0);
209
210             /* Calculate squared distance and things based on it */
211             rsq10            = gmx_fjsp_calc_rsq_v2r8(dx10,dy10,dz10);
212             rsq20            = gmx_fjsp_calc_rsq_v2r8(dx20,dy20,dz20);
213             rsq30            = gmx_fjsp_calc_rsq_v2r8(dx30,dy30,dz30);
214
215             rinv10           = gmx_fjsp_invsqrt_v2r8(rsq10);
216             rinv20           = gmx_fjsp_invsqrt_v2r8(rsq20);
217             rinv30           = gmx_fjsp_invsqrt_v2r8(rsq30);
218
219             rinvsq10         = _fjsp_mul_v2r8(rinv10,rinv10);
220             rinvsq20         = _fjsp_mul_v2r8(rinv20,rinv20);
221             rinvsq30         = _fjsp_mul_v2r8(rinv30,rinv30);
222
223             /* Load parameters for j particles */
224             jq0              = gmx_fjsp_load_2real_swizzle_v2r8(charge+jnrA+0,charge+jnrB+0);
225
226             fjx0             = _fjsp_setzero_v2r8();
227             fjy0             = _fjsp_setzero_v2r8();
228             fjz0             = _fjsp_setzero_v2r8();
229
230             /**************************
231              * CALCULATE INTERACTIONS *
232              **************************/
233
234             if (gmx_fjsp_any_lt_v2r8(rsq10,rcutoff2))
235             {
236
237             r10              = _fjsp_mul_v2r8(rsq10,rinv10);
238
239             /* Compute parameters for interactions between i and j atoms */
240             qq10             = _fjsp_mul_v2r8(iq1,jq0);
241
242             /* EWALD ELECTROSTATICS */
243
244             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
245             ewrt             = _fjsp_mul_v2r8(r10,ewtabscale);
246             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
247             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
248             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
249
250             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
251             ewtabD           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[1] );
252             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
253             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
254             ewtabFn          = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[1] +2);
255             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
256             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
257             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
258             velec            = _fjsp_mul_v2r8(qq10,_fjsp_sub_v2r8(rinv10,velec));
259             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq10,rinv10),_fjsp_sub_v2r8(rinvsq10,felec));
260
261             d                = _fjsp_sub_v2r8(r10,rswitch);
262             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
263             d2               = _fjsp_mul_v2r8(d,d);
264             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
265
266             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
267
268             /* Evaluate switch function */
269             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
270             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv10,_fjsp_mul_v2r8(velec,dsw)) );
271             velec            = _fjsp_mul_v2r8(velec,sw);
272             cutoff_mask      = _fjsp_cmplt_v2r8(rsq10,rcutoff2);
273
274             /* Update potential sum for this i atom from the interaction with this j atom. */
275             velec            = _fjsp_and_v2r8(velec,cutoff_mask);
276             velecsum         = _fjsp_add_v2r8(velecsum,velec);
277
278             fscal            = felec;
279
280             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
281
282             /* Update vectorial force */
283             fix1             = _fjsp_madd_v2r8(dx10,fscal,fix1);
284             fiy1             = _fjsp_madd_v2r8(dy10,fscal,fiy1);
285             fiz1             = _fjsp_madd_v2r8(dz10,fscal,fiz1);
286             
287             fjx0             = _fjsp_madd_v2r8(dx10,fscal,fjx0);
288             fjy0             = _fjsp_madd_v2r8(dy10,fscal,fjy0);
289             fjz0             = _fjsp_madd_v2r8(dz10,fscal,fjz0);
290
291             }
292
293             /**************************
294              * CALCULATE INTERACTIONS *
295              **************************/
296
297             if (gmx_fjsp_any_lt_v2r8(rsq20,rcutoff2))
298             {
299
300             r20              = _fjsp_mul_v2r8(rsq20,rinv20);
301
302             /* Compute parameters for interactions between i and j atoms */
303             qq20             = _fjsp_mul_v2r8(iq2,jq0);
304
305             /* EWALD ELECTROSTATICS */
306
307             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
308             ewrt             = _fjsp_mul_v2r8(r20,ewtabscale);
309             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
310             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
311             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
312
313             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
314             ewtabD           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[1] );
315             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
316             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
317             ewtabFn          = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[1] +2);
318             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
319             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
320             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
321             velec            = _fjsp_mul_v2r8(qq20,_fjsp_sub_v2r8(rinv20,velec));
322             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq20,rinv20),_fjsp_sub_v2r8(rinvsq20,felec));
323
324             d                = _fjsp_sub_v2r8(r20,rswitch);
325             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
326             d2               = _fjsp_mul_v2r8(d,d);
327             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
328
329             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
330
331             /* Evaluate switch function */
332             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
333             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv20,_fjsp_mul_v2r8(velec,dsw)) );
334             velec            = _fjsp_mul_v2r8(velec,sw);
335             cutoff_mask      = _fjsp_cmplt_v2r8(rsq20,rcutoff2);
336
337             /* Update potential sum for this i atom from the interaction with this j atom. */
338             velec            = _fjsp_and_v2r8(velec,cutoff_mask);
339             velecsum         = _fjsp_add_v2r8(velecsum,velec);
340
341             fscal            = felec;
342
343             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
344
345             /* Update vectorial force */
346             fix2             = _fjsp_madd_v2r8(dx20,fscal,fix2);
347             fiy2             = _fjsp_madd_v2r8(dy20,fscal,fiy2);
348             fiz2             = _fjsp_madd_v2r8(dz20,fscal,fiz2);
349             
350             fjx0             = _fjsp_madd_v2r8(dx20,fscal,fjx0);
351             fjy0             = _fjsp_madd_v2r8(dy20,fscal,fjy0);
352             fjz0             = _fjsp_madd_v2r8(dz20,fscal,fjz0);
353
354             }
355
356             /**************************
357              * CALCULATE INTERACTIONS *
358              **************************/
359
360             if (gmx_fjsp_any_lt_v2r8(rsq30,rcutoff2))
361             {
362
363             r30              = _fjsp_mul_v2r8(rsq30,rinv30);
364
365             /* Compute parameters for interactions between i and j atoms */
366             qq30             = _fjsp_mul_v2r8(iq3,jq0);
367
368             /* EWALD ELECTROSTATICS */
369
370             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
371             ewrt             = _fjsp_mul_v2r8(r30,ewtabscale);
372             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
373             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
374             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
375
376             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
377             ewtabD           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[1] );
378             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
379             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
380             ewtabFn          = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[1] +2);
381             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
382             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
383             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
384             velec            = _fjsp_mul_v2r8(qq30,_fjsp_sub_v2r8(rinv30,velec));
385             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq30,rinv30),_fjsp_sub_v2r8(rinvsq30,felec));
386
387             d                = _fjsp_sub_v2r8(r30,rswitch);
388             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
389             d2               = _fjsp_mul_v2r8(d,d);
390             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
391
392             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
393
394             /* Evaluate switch function */
395             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
396             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv30,_fjsp_mul_v2r8(velec,dsw)) );
397             velec            = _fjsp_mul_v2r8(velec,sw);
398             cutoff_mask      = _fjsp_cmplt_v2r8(rsq30,rcutoff2);
399
400             /* Update potential sum for this i atom from the interaction with this j atom. */
401             velec            = _fjsp_and_v2r8(velec,cutoff_mask);
402             velecsum         = _fjsp_add_v2r8(velecsum,velec);
403
404             fscal            = felec;
405
406             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
407
408             /* Update vectorial force */
409             fix3             = _fjsp_madd_v2r8(dx30,fscal,fix3);
410             fiy3             = _fjsp_madd_v2r8(dy30,fscal,fiy3);
411             fiz3             = _fjsp_madd_v2r8(dz30,fscal,fiz3);
412             
413             fjx0             = _fjsp_madd_v2r8(dx30,fscal,fjx0);
414             fjy0             = _fjsp_madd_v2r8(dy30,fscal,fjy0);
415             fjz0             = _fjsp_madd_v2r8(dz30,fscal,fjz0);
416
417             }
418
419             gmx_fjsp_decrement_1rvec_2ptr_swizzle_v2r8(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0);
420
421             /* Inner loop uses 207 flops */
422         }
423
424         if(jidx<j_index_end)
425         {
426
427             jnrA             = jjnr[jidx];
428             j_coord_offsetA  = DIM*jnrA;
429
430             /* load j atom coordinates */
431             gmx_fjsp_load_1rvec_1ptr_swizzle_v2r8(x+j_coord_offsetA,
432                                               &jx0,&jy0,&jz0);
433
434             /* Calculate displacement vector */
435             dx10             = _fjsp_sub_v2r8(ix1,jx0);
436             dy10             = _fjsp_sub_v2r8(iy1,jy0);
437             dz10             = _fjsp_sub_v2r8(iz1,jz0);
438             dx20             = _fjsp_sub_v2r8(ix2,jx0);
439             dy20             = _fjsp_sub_v2r8(iy2,jy0);
440             dz20             = _fjsp_sub_v2r8(iz2,jz0);
441             dx30             = _fjsp_sub_v2r8(ix3,jx0);
442             dy30             = _fjsp_sub_v2r8(iy3,jy0);
443             dz30             = _fjsp_sub_v2r8(iz3,jz0);
444
445             /* Calculate squared distance and things based on it */
446             rsq10            = gmx_fjsp_calc_rsq_v2r8(dx10,dy10,dz10);
447             rsq20            = gmx_fjsp_calc_rsq_v2r8(dx20,dy20,dz20);
448             rsq30            = gmx_fjsp_calc_rsq_v2r8(dx30,dy30,dz30);
449
450             rinv10           = gmx_fjsp_invsqrt_v2r8(rsq10);
451             rinv20           = gmx_fjsp_invsqrt_v2r8(rsq20);
452             rinv30           = gmx_fjsp_invsqrt_v2r8(rsq30);
453
454             rinvsq10         = _fjsp_mul_v2r8(rinv10,rinv10);
455             rinvsq20         = _fjsp_mul_v2r8(rinv20,rinv20);
456             rinvsq30         = _fjsp_mul_v2r8(rinv30,rinv30);
457
458             /* Load parameters for j particles */
459             jq0              = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(),charge+jnrA+0);
460
461             fjx0             = _fjsp_setzero_v2r8();
462             fjy0             = _fjsp_setzero_v2r8();
463             fjz0             = _fjsp_setzero_v2r8();
464
465             /**************************
466              * CALCULATE INTERACTIONS *
467              **************************/
468
469             if (gmx_fjsp_any_lt_v2r8(rsq10,rcutoff2))
470             {
471
472             r10              = _fjsp_mul_v2r8(rsq10,rinv10);
473
474             /* Compute parameters for interactions between i and j atoms */
475             qq10             = _fjsp_mul_v2r8(iq1,jq0);
476
477             /* EWALD ELECTROSTATICS */
478
479             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
480             ewrt             = _fjsp_mul_v2r8(r10,ewtabscale);
481             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
482             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
483             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
484
485             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
486             ewtabD           = _fjsp_setzero_v2r8();
487             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
488             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
489             ewtabFn          = _fjsp_setzero_v2r8();
490             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
491             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
492             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
493             velec            = _fjsp_mul_v2r8(qq10,_fjsp_sub_v2r8(rinv10,velec));
494             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq10,rinv10),_fjsp_sub_v2r8(rinvsq10,felec));
495
496             d                = _fjsp_sub_v2r8(r10,rswitch);
497             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
498             d2               = _fjsp_mul_v2r8(d,d);
499             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
500
501             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
502
503             /* Evaluate switch function */
504             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
505             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv10,_fjsp_mul_v2r8(velec,dsw)) );
506             velec            = _fjsp_mul_v2r8(velec,sw);
507             cutoff_mask      = _fjsp_cmplt_v2r8(rsq10,rcutoff2);
508
509             /* Update potential sum for this i atom from the interaction with this j atom. */
510             velec            = _fjsp_and_v2r8(velec,cutoff_mask);
511             velec            = _fjsp_unpacklo_v2r8(velec,_fjsp_setzero_v2r8());
512             velecsum         = _fjsp_add_v2r8(velecsum,velec);
513
514             fscal            = felec;
515
516             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
517
518             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
519
520             /* Update vectorial force */
521             fix1             = _fjsp_madd_v2r8(dx10,fscal,fix1);
522             fiy1             = _fjsp_madd_v2r8(dy10,fscal,fiy1);
523             fiz1             = _fjsp_madd_v2r8(dz10,fscal,fiz1);
524             
525             fjx0             = _fjsp_madd_v2r8(dx10,fscal,fjx0);
526             fjy0             = _fjsp_madd_v2r8(dy10,fscal,fjy0);
527             fjz0             = _fjsp_madd_v2r8(dz10,fscal,fjz0);
528
529             }
530
531             /**************************
532              * CALCULATE INTERACTIONS *
533              **************************/
534
535             if (gmx_fjsp_any_lt_v2r8(rsq20,rcutoff2))
536             {
537
538             r20              = _fjsp_mul_v2r8(rsq20,rinv20);
539
540             /* Compute parameters for interactions between i and j atoms */
541             qq20             = _fjsp_mul_v2r8(iq2,jq0);
542
543             /* EWALD ELECTROSTATICS */
544
545             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
546             ewrt             = _fjsp_mul_v2r8(r20,ewtabscale);
547             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
548             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
549             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
550
551             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
552             ewtabD           = _fjsp_setzero_v2r8();
553             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
554             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
555             ewtabFn          = _fjsp_setzero_v2r8();
556             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
557             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
558             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
559             velec            = _fjsp_mul_v2r8(qq20,_fjsp_sub_v2r8(rinv20,velec));
560             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq20,rinv20),_fjsp_sub_v2r8(rinvsq20,felec));
561
562             d                = _fjsp_sub_v2r8(r20,rswitch);
563             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
564             d2               = _fjsp_mul_v2r8(d,d);
565             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
566
567             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
568
569             /* Evaluate switch function */
570             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
571             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv20,_fjsp_mul_v2r8(velec,dsw)) );
572             velec            = _fjsp_mul_v2r8(velec,sw);
573             cutoff_mask      = _fjsp_cmplt_v2r8(rsq20,rcutoff2);
574
575             /* Update potential sum for this i atom from the interaction with this j atom. */
576             velec            = _fjsp_and_v2r8(velec,cutoff_mask);
577             velec            = _fjsp_unpacklo_v2r8(velec,_fjsp_setzero_v2r8());
578             velecsum         = _fjsp_add_v2r8(velecsum,velec);
579
580             fscal            = felec;
581
582             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
583
584             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
585
586             /* Update vectorial force */
587             fix2             = _fjsp_madd_v2r8(dx20,fscal,fix2);
588             fiy2             = _fjsp_madd_v2r8(dy20,fscal,fiy2);
589             fiz2             = _fjsp_madd_v2r8(dz20,fscal,fiz2);
590             
591             fjx0             = _fjsp_madd_v2r8(dx20,fscal,fjx0);
592             fjy0             = _fjsp_madd_v2r8(dy20,fscal,fjy0);
593             fjz0             = _fjsp_madd_v2r8(dz20,fscal,fjz0);
594
595             }
596
597             /**************************
598              * CALCULATE INTERACTIONS *
599              **************************/
600
601             if (gmx_fjsp_any_lt_v2r8(rsq30,rcutoff2))
602             {
603
604             r30              = _fjsp_mul_v2r8(rsq30,rinv30);
605
606             /* Compute parameters for interactions between i and j atoms */
607             qq30             = _fjsp_mul_v2r8(iq3,jq0);
608
609             /* EWALD ELECTROSTATICS */
610
611             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
612             ewrt             = _fjsp_mul_v2r8(r30,ewtabscale);
613             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
614             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
615             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
616
617             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
618             ewtabD           = _fjsp_setzero_v2r8();
619             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
620             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
621             ewtabFn          = _fjsp_setzero_v2r8();
622             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
623             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
624             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
625             velec            = _fjsp_mul_v2r8(qq30,_fjsp_sub_v2r8(rinv30,velec));
626             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq30,rinv30),_fjsp_sub_v2r8(rinvsq30,felec));
627
628             d                = _fjsp_sub_v2r8(r30,rswitch);
629             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
630             d2               = _fjsp_mul_v2r8(d,d);
631             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
632
633             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
634
635             /* Evaluate switch function */
636             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
637             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv30,_fjsp_mul_v2r8(velec,dsw)) );
638             velec            = _fjsp_mul_v2r8(velec,sw);
639             cutoff_mask      = _fjsp_cmplt_v2r8(rsq30,rcutoff2);
640
641             /* Update potential sum for this i atom from the interaction with this j atom. */
642             velec            = _fjsp_and_v2r8(velec,cutoff_mask);
643             velec            = _fjsp_unpacklo_v2r8(velec,_fjsp_setzero_v2r8());
644             velecsum         = _fjsp_add_v2r8(velecsum,velec);
645
646             fscal            = felec;
647
648             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
649
650             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
651
652             /* Update vectorial force */
653             fix3             = _fjsp_madd_v2r8(dx30,fscal,fix3);
654             fiy3             = _fjsp_madd_v2r8(dy30,fscal,fiy3);
655             fiz3             = _fjsp_madd_v2r8(dz30,fscal,fiz3);
656             
657             fjx0             = _fjsp_madd_v2r8(dx30,fscal,fjx0);
658             fjy0             = _fjsp_madd_v2r8(dy30,fscal,fjy0);
659             fjz0             = _fjsp_madd_v2r8(dz30,fscal,fjz0);
660
661             }
662
663             gmx_fjsp_decrement_1rvec_1ptr_swizzle_v2r8(f+j_coord_offsetA,fjx0,fjy0,fjz0);
664
665             /* Inner loop uses 207 flops */
666         }
667
668         /* End of innermost loop */
669
670         gmx_fjsp_update_iforce_3atom_swizzle_v2r8(fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
671                                               f+i_coord_offset+DIM,fshift+i_shift_offset);
672
673         ggid                        = gid[iidx];
674         /* Update potential energies */
675         gmx_fjsp_update_1pot_v2r8(velecsum,kernel_data->energygrp_elec+ggid);
676
677         /* Increment number of inner iterations */
678         inneriter                  += j_index_end - j_index_start;
679
680         /* Outer loop uses 19 flops */
681     }
682
683     /* Increment number of outer iterations */
684     outeriter        += nri;
685
686     /* Update outer/inner flops */
687
688     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W4_VF,outeriter*19 + inneriter*207);
689 }
690 /*
691  * Gromacs nonbonded kernel:   nb_kernel_ElecEwSw_VdwNone_GeomW4P1_F_sparc64_hpc_ace_double
692  * Electrostatics interaction: Ewald
693  * VdW interaction:            None
694  * Geometry:                   Water4-Particle
695  * Calculate force/pot:        Force
696  */
697 void
698 nb_kernel_ElecEwSw_VdwNone_GeomW4P1_F_sparc64_hpc_ace_double
699                     (t_nblist                    * gmx_restrict       nlist,
700                      rvec                        * gmx_restrict          xx,
701                      rvec                        * gmx_restrict          ff,
702                      t_forcerec                  * gmx_restrict          fr,
703                      t_mdatoms                   * gmx_restrict     mdatoms,
704                      nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
705                      t_nrnb                      * gmx_restrict        nrnb)
706 {
707     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
708      * just 0 for non-waters.
709      * Suffixes A,B refer to j loop unrolling done with double precision SIMD, e.g. for the two different
710      * jnr indices corresponding to data put in the four positions in the SIMD register.
711      */
712     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
713     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
714     int              jnrA,jnrB;
715     int              j_coord_offsetA,j_coord_offsetB;
716     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
717     real             rcutoff_scalar;
718     real             *shiftvec,*fshift,*x,*f;
719     _fjsp_v2r8       tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
720     int              vdwioffset1;
721     _fjsp_v2r8       ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
722     int              vdwioffset2;
723     _fjsp_v2r8       ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
724     int              vdwioffset3;
725     _fjsp_v2r8       ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
726     int              vdwjidx0A,vdwjidx0B;
727     _fjsp_v2r8       jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
728     _fjsp_v2r8       dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
729     _fjsp_v2r8       dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
730     _fjsp_v2r8       dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
731     _fjsp_v2r8       velec,felec,velecsum,facel,crf,krf,krf2;
732     real             *charge;
733     _fjsp_v2r8       ewtabscale,eweps,sh_ewald,ewrt,ewtabhalfspace,ewtabF,ewtabFn,ewtabD,ewtabV;
734     real             *ewtab;
735     _fjsp_v2r8       rswitch,swV3,swV4,swV5,swF2,swF3,swF4,d,d2,sw,dsw;
736     real             rswitch_scalar,d_scalar;
737     _fjsp_v2r8       itab_tmp;
738     _fjsp_v2r8       dummy_mask,cutoff_mask;
739     _fjsp_v2r8       one     = gmx_fjsp_set1_v2r8(1.0);
740     _fjsp_v2r8       two     = gmx_fjsp_set1_v2r8(2.0);
741     union { _fjsp_v2r8 simd; long long int i[2]; } vfconv,gbconv,ewconv;
742
743     x                = xx[0];
744     f                = ff[0];
745
746     nri              = nlist->nri;
747     iinr             = nlist->iinr;
748     jindex           = nlist->jindex;
749     jjnr             = nlist->jjnr;
750     shiftidx         = nlist->shift;
751     gid              = nlist->gid;
752     shiftvec         = fr->shift_vec[0];
753     fshift           = fr->fshift[0];
754     facel            = gmx_fjsp_set1_v2r8(fr->epsfac);
755     charge           = mdatoms->chargeA;
756
757     sh_ewald         = gmx_fjsp_set1_v2r8(fr->ic->sh_ewald);
758     ewtab            = fr->ic->tabq_coul_FDV0;
759     ewtabscale       = gmx_fjsp_set1_v2r8(fr->ic->tabq_scale);
760     ewtabhalfspace   = gmx_fjsp_set1_v2r8(0.5/fr->ic->tabq_scale);
761
762     /* Setup water-specific parameters */
763     inr              = nlist->iinr[0];
764     iq1              = _fjsp_mul_v2r8(facel,gmx_fjsp_set1_v2r8(charge[inr+1]));
765     iq2              = _fjsp_mul_v2r8(facel,gmx_fjsp_set1_v2r8(charge[inr+2]));
766     iq3              = _fjsp_mul_v2r8(facel,gmx_fjsp_set1_v2r8(charge[inr+3]));
767
768     /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
769     rcutoff_scalar   = fr->rcoulomb;
770     rcutoff          = gmx_fjsp_set1_v2r8(rcutoff_scalar);
771     rcutoff2         = _fjsp_mul_v2r8(rcutoff,rcutoff);
772
773     rswitch_scalar   = fr->rcoulomb_switch;
774     rswitch          = gmx_fjsp_set1_v2r8(rswitch_scalar);
775     /* Setup switch parameters */
776     d_scalar         = rcutoff_scalar-rswitch_scalar;
777     d                = gmx_fjsp_set1_v2r8(d_scalar);
778     swV3             = gmx_fjsp_set1_v2r8(-10.0/(d_scalar*d_scalar*d_scalar));
779     swV4             = gmx_fjsp_set1_v2r8( 15.0/(d_scalar*d_scalar*d_scalar*d_scalar));
780     swV5             = gmx_fjsp_set1_v2r8( -6.0/(d_scalar*d_scalar*d_scalar*d_scalar*d_scalar));
781     swF2             = gmx_fjsp_set1_v2r8(-30.0/(d_scalar*d_scalar*d_scalar));
782     swF3             = gmx_fjsp_set1_v2r8( 60.0/(d_scalar*d_scalar*d_scalar*d_scalar));
783     swF4             = gmx_fjsp_set1_v2r8(-30.0/(d_scalar*d_scalar*d_scalar*d_scalar*d_scalar));
784
785     /* Avoid stupid compiler warnings */
786     jnrA = jnrB = 0;
787     j_coord_offsetA = 0;
788     j_coord_offsetB = 0;
789
790     outeriter        = 0;
791     inneriter        = 0;
792
793     /* Start outer loop over neighborlists */
794     for(iidx=0; iidx<nri; iidx++)
795     {
796         /* Load shift vector for this list */
797         i_shift_offset   = DIM*shiftidx[iidx];
798
799         /* Load limits for loop over neighbors */
800         j_index_start    = jindex[iidx];
801         j_index_end      = jindex[iidx+1];
802
803         /* Get outer coordinate index */
804         inr              = iinr[iidx];
805         i_coord_offset   = DIM*inr;
806
807         /* Load i particle coords and add shift vector */
808         gmx_fjsp_load_shift_and_3rvec_broadcast_v2r8(shiftvec+i_shift_offset,x+i_coord_offset+DIM,
809                                                  &ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
810
811         fix1             = _fjsp_setzero_v2r8();
812         fiy1             = _fjsp_setzero_v2r8();
813         fiz1             = _fjsp_setzero_v2r8();
814         fix2             = _fjsp_setzero_v2r8();
815         fiy2             = _fjsp_setzero_v2r8();
816         fiz2             = _fjsp_setzero_v2r8();
817         fix3             = _fjsp_setzero_v2r8();
818         fiy3             = _fjsp_setzero_v2r8();
819         fiz3             = _fjsp_setzero_v2r8();
820
821         /* Start inner kernel loop */
822         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
823         {
824
825             /* Get j neighbor index, and coordinate index */
826             jnrA             = jjnr[jidx];
827             jnrB             = jjnr[jidx+1];
828             j_coord_offsetA  = DIM*jnrA;
829             j_coord_offsetB  = DIM*jnrB;
830
831             /* load j atom coordinates */
832             gmx_fjsp_load_1rvec_2ptr_swizzle_v2r8(x+j_coord_offsetA,x+j_coord_offsetB,
833                                               &jx0,&jy0,&jz0);
834
835             /* Calculate displacement vector */
836             dx10             = _fjsp_sub_v2r8(ix1,jx0);
837             dy10             = _fjsp_sub_v2r8(iy1,jy0);
838             dz10             = _fjsp_sub_v2r8(iz1,jz0);
839             dx20             = _fjsp_sub_v2r8(ix2,jx0);
840             dy20             = _fjsp_sub_v2r8(iy2,jy0);
841             dz20             = _fjsp_sub_v2r8(iz2,jz0);
842             dx30             = _fjsp_sub_v2r8(ix3,jx0);
843             dy30             = _fjsp_sub_v2r8(iy3,jy0);
844             dz30             = _fjsp_sub_v2r8(iz3,jz0);
845
846             /* Calculate squared distance and things based on it */
847             rsq10            = gmx_fjsp_calc_rsq_v2r8(dx10,dy10,dz10);
848             rsq20            = gmx_fjsp_calc_rsq_v2r8(dx20,dy20,dz20);
849             rsq30            = gmx_fjsp_calc_rsq_v2r8(dx30,dy30,dz30);
850
851             rinv10           = gmx_fjsp_invsqrt_v2r8(rsq10);
852             rinv20           = gmx_fjsp_invsqrt_v2r8(rsq20);
853             rinv30           = gmx_fjsp_invsqrt_v2r8(rsq30);
854
855             rinvsq10         = _fjsp_mul_v2r8(rinv10,rinv10);
856             rinvsq20         = _fjsp_mul_v2r8(rinv20,rinv20);
857             rinvsq30         = _fjsp_mul_v2r8(rinv30,rinv30);
858
859             /* Load parameters for j particles */
860             jq0              = gmx_fjsp_load_2real_swizzle_v2r8(charge+jnrA+0,charge+jnrB+0);
861
862             fjx0             = _fjsp_setzero_v2r8();
863             fjy0             = _fjsp_setzero_v2r8();
864             fjz0             = _fjsp_setzero_v2r8();
865
866             /**************************
867              * CALCULATE INTERACTIONS *
868              **************************/
869
870             if (gmx_fjsp_any_lt_v2r8(rsq10,rcutoff2))
871             {
872
873             r10              = _fjsp_mul_v2r8(rsq10,rinv10);
874
875             /* Compute parameters for interactions between i and j atoms */
876             qq10             = _fjsp_mul_v2r8(iq1,jq0);
877
878             /* EWALD ELECTROSTATICS */
879
880             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
881             ewrt             = _fjsp_mul_v2r8(r10,ewtabscale);
882             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
883             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
884             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
885
886             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
887             ewtabD           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[1] );
888             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
889             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
890             ewtabFn          = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[1] +2);
891             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
892             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
893             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
894             velec            = _fjsp_mul_v2r8(qq10,_fjsp_sub_v2r8(rinv10,velec));
895             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq10,rinv10),_fjsp_sub_v2r8(rinvsq10,felec));
896
897             d                = _fjsp_sub_v2r8(r10,rswitch);
898             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
899             d2               = _fjsp_mul_v2r8(d,d);
900             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
901
902             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
903
904             /* Evaluate switch function */
905             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
906             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv10,_fjsp_mul_v2r8(velec,dsw)) );
907             cutoff_mask      = _fjsp_cmplt_v2r8(rsq10,rcutoff2);
908
909             fscal            = felec;
910
911             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
912
913             /* Update vectorial force */
914             fix1             = _fjsp_madd_v2r8(dx10,fscal,fix1);
915             fiy1             = _fjsp_madd_v2r8(dy10,fscal,fiy1);
916             fiz1             = _fjsp_madd_v2r8(dz10,fscal,fiz1);
917             
918             fjx0             = _fjsp_madd_v2r8(dx10,fscal,fjx0);
919             fjy0             = _fjsp_madd_v2r8(dy10,fscal,fjy0);
920             fjz0             = _fjsp_madd_v2r8(dz10,fscal,fjz0);
921
922             }
923
924             /**************************
925              * CALCULATE INTERACTIONS *
926              **************************/
927
928             if (gmx_fjsp_any_lt_v2r8(rsq20,rcutoff2))
929             {
930
931             r20              = _fjsp_mul_v2r8(rsq20,rinv20);
932
933             /* Compute parameters for interactions between i and j atoms */
934             qq20             = _fjsp_mul_v2r8(iq2,jq0);
935
936             /* EWALD ELECTROSTATICS */
937
938             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
939             ewrt             = _fjsp_mul_v2r8(r20,ewtabscale);
940             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
941             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
942             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
943
944             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
945             ewtabD           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[1] );
946             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
947             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
948             ewtabFn          = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[1] +2);
949             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
950             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
951             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
952             velec            = _fjsp_mul_v2r8(qq20,_fjsp_sub_v2r8(rinv20,velec));
953             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq20,rinv20),_fjsp_sub_v2r8(rinvsq20,felec));
954
955             d                = _fjsp_sub_v2r8(r20,rswitch);
956             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
957             d2               = _fjsp_mul_v2r8(d,d);
958             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
959
960             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
961
962             /* Evaluate switch function */
963             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
964             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv20,_fjsp_mul_v2r8(velec,dsw)) );
965             cutoff_mask      = _fjsp_cmplt_v2r8(rsq20,rcutoff2);
966
967             fscal            = felec;
968
969             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
970
971             /* Update vectorial force */
972             fix2             = _fjsp_madd_v2r8(dx20,fscal,fix2);
973             fiy2             = _fjsp_madd_v2r8(dy20,fscal,fiy2);
974             fiz2             = _fjsp_madd_v2r8(dz20,fscal,fiz2);
975             
976             fjx0             = _fjsp_madd_v2r8(dx20,fscal,fjx0);
977             fjy0             = _fjsp_madd_v2r8(dy20,fscal,fjy0);
978             fjz0             = _fjsp_madd_v2r8(dz20,fscal,fjz0);
979
980             }
981
982             /**************************
983              * CALCULATE INTERACTIONS *
984              **************************/
985
986             if (gmx_fjsp_any_lt_v2r8(rsq30,rcutoff2))
987             {
988
989             r30              = _fjsp_mul_v2r8(rsq30,rinv30);
990
991             /* Compute parameters for interactions between i and j atoms */
992             qq30             = _fjsp_mul_v2r8(iq3,jq0);
993
994             /* EWALD ELECTROSTATICS */
995
996             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
997             ewrt             = _fjsp_mul_v2r8(r30,ewtabscale);
998             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
999             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
1000             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
1001
1002             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
1003             ewtabD           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[1] );
1004             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
1005             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
1006             ewtabFn          = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[1] +2);
1007             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
1008             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
1009             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
1010             velec            = _fjsp_mul_v2r8(qq30,_fjsp_sub_v2r8(rinv30,velec));
1011             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq30,rinv30),_fjsp_sub_v2r8(rinvsq30,felec));
1012
1013             d                = _fjsp_sub_v2r8(r30,rswitch);
1014             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
1015             d2               = _fjsp_mul_v2r8(d,d);
1016             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
1017
1018             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
1019
1020             /* Evaluate switch function */
1021             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
1022             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv30,_fjsp_mul_v2r8(velec,dsw)) );
1023             cutoff_mask      = _fjsp_cmplt_v2r8(rsq30,rcutoff2);
1024
1025             fscal            = felec;
1026
1027             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
1028
1029             /* Update vectorial force */
1030             fix3             = _fjsp_madd_v2r8(dx30,fscal,fix3);
1031             fiy3             = _fjsp_madd_v2r8(dy30,fscal,fiy3);
1032             fiz3             = _fjsp_madd_v2r8(dz30,fscal,fiz3);
1033             
1034             fjx0             = _fjsp_madd_v2r8(dx30,fscal,fjx0);
1035             fjy0             = _fjsp_madd_v2r8(dy30,fscal,fjy0);
1036             fjz0             = _fjsp_madd_v2r8(dz30,fscal,fjz0);
1037
1038             }
1039
1040             gmx_fjsp_decrement_1rvec_2ptr_swizzle_v2r8(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0);
1041
1042             /* Inner loop uses 198 flops */
1043         }
1044
1045         if(jidx<j_index_end)
1046         {
1047
1048             jnrA             = jjnr[jidx];
1049             j_coord_offsetA  = DIM*jnrA;
1050
1051             /* load j atom coordinates */
1052             gmx_fjsp_load_1rvec_1ptr_swizzle_v2r8(x+j_coord_offsetA,
1053                                               &jx0,&jy0,&jz0);
1054
1055             /* Calculate displacement vector */
1056             dx10             = _fjsp_sub_v2r8(ix1,jx0);
1057             dy10             = _fjsp_sub_v2r8(iy1,jy0);
1058             dz10             = _fjsp_sub_v2r8(iz1,jz0);
1059             dx20             = _fjsp_sub_v2r8(ix2,jx0);
1060             dy20             = _fjsp_sub_v2r8(iy2,jy0);
1061             dz20             = _fjsp_sub_v2r8(iz2,jz0);
1062             dx30             = _fjsp_sub_v2r8(ix3,jx0);
1063             dy30             = _fjsp_sub_v2r8(iy3,jy0);
1064             dz30             = _fjsp_sub_v2r8(iz3,jz0);
1065
1066             /* Calculate squared distance and things based on it */
1067             rsq10            = gmx_fjsp_calc_rsq_v2r8(dx10,dy10,dz10);
1068             rsq20            = gmx_fjsp_calc_rsq_v2r8(dx20,dy20,dz20);
1069             rsq30            = gmx_fjsp_calc_rsq_v2r8(dx30,dy30,dz30);
1070
1071             rinv10           = gmx_fjsp_invsqrt_v2r8(rsq10);
1072             rinv20           = gmx_fjsp_invsqrt_v2r8(rsq20);
1073             rinv30           = gmx_fjsp_invsqrt_v2r8(rsq30);
1074
1075             rinvsq10         = _fjsp_mul_v2r8(rinv10,rinv10);
1076             rinvsq20         = _fjsp_mul_v2r8(rinv20,rinv20);
1077             rinvsq30         = _fjsp_mul_v2r8(rinv30,rinv30);
1078
1079             /* Load parameters for j particles */
1080             jq0              = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(),charge+jnrA+0);
1081
1082             fjx0             = _fjsp_setzero_v2r8();
1083             fjy0             = _fjsp_setzero_v2r8();
1084             fjz0             = _fjsp_setzero_v2r8();
1085
1086             /**************************
1087              * CALCULATE INTERACTIONS *
1088              **************************/
1089
1090             if (gmx_fjsp_any_lt_v2r8(rsq10,rcutoff2))
1091             {
1092
1093             r10              = _fjsp_mul_v2r8(rsq10,rinv10);
1094
1095             /* Compute parameters for interactions between i and j atoms */
1096             qq10             = _fjsp_mul_v2r8(iq1,jq0);
1097
1098             /* EWALD ELECTROSTATICS */
1099
1100             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1101             ewrt             = _fjsp_mul_v2r8(r10,ewtabscale);
1102             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
1103             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
1104             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
1105
1106             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
1107             ewtabD           = _fjsp_setzero_v2r8();
1108             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
1109             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
1110             ewtabFn          = _fjsp_setzero_v2r8();
1111             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
1112             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
1113             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
1114             velec            = _fjsp_mul_v2r8(qq10,_fjsp_sub_v2r8(rinv10,velec));
1115             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq10,rinv10),_fjsp_sub_v2r8(rinvsq10,felec));
1116
1117             d                = _fjsp_sub_v2r8(r10,rswitch);
1118             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
1119             d2               = _fjsp_mul_v2r8(d,d);
1120             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
1121
1122             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
1123
1124             /* Evaluate switch function */
1125             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
1126             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv10,_fjsp_mul_v2r8(velec,dsw)) );
1127             cutoff_mask      = _fjsp_cmplt_v2r8(rsq10,rcutoff2);
1128
1129             fscal            = felec;
1130
1131             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
1132
1133             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
1134
1135             /* Update vectorial force */
1136             fix1             = _fjsp_madd_v2r8(dx10,fscal,fix1);
1137             fiy1             = _fjsp_madd_v2r8(dy10,fscal,fiy1);
1138             fiz1             = _fjsp_madd_v2r8(dz10,fscal,fiz1);
1139             
1140             fjx0             = _fjsp_madd_v2r8(dx10,fscal,fjx0);
1141             fjy0             = _fjsp_madd_v2r8(dy10,fscal,fjy0);
1142             fjz0             = _fjsp_madd_v2r8(dz10,fscal,fjz0);
1143
1144             }
1145
1146             /**************************
1147              * CALCULATE INTERACTIONS *
1148              **************************/
1149
1150             if (gmx_fjsp_any_lt_v2r8(rsq20,rcutoff2))
1151             {
1152
1153             r20              = _fjsp_mul_v2r8(rsq20,rinv20);
1154
1155             /* Compute parameters for interactions between i and j atoms */
1156             qq20             = _fjsp_mul_v2r8(iq2,jq0);
1157
1158             /* EWALD ELECTROSTATICS */
1159
1160             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1161             ewrt             = _fjsp_mul_v2r8(r20,ewtabscale);
1162             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
1163             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
1164             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
1165
1166             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
1167             ewtabD           = _fjsp_setzero_v2r8();
1168             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
1169             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
1170             ewtabFn          = _fjsp_setzero_v2r8();
1171             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
1172             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
1173             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
1174             velec            = _fjsp_mul_v2r8(qq20,_fjsp_sub_v2r8(rinv20,velec));
1175             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq20,rinv20),_fjsp_sub_v2r8(rinvsq20,felec));
1176
1177             d                = _fjsp_sub_v2r8(r20,rswitch);
1178             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
1179             d2               = _fjsp_mul_v2r8(d,d);
1180             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
1181
1182             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
1183
1184             /* Evaluate switch function */
1185             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
1186             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv20,_fjsp_mul_v2r8(velec,dsw)) );
1187             cutoff_mask      = _fjsp_cmplt_v2r8(rsq20,rcutoff2);
1188
1189             fscal            = felec;
1190
1191             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
1192
1193             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
1194
1195             /* Update vectorial force */
1196             fix2             = _fjsp_madd_v2r8(dx20,fscal,fix2);
1197             fiy2             = _fjsp_madd_v2r8(dy20,fscal,fiy2);
1198             fiz2             = _fjsp_madd_v2r8(dz20,fscal,fiz2);
1199             
1200             fjx0             = _fjsp_madd_v2r8(dx20,fscal,fjx0);
1201             fjy0             = _fjsp_madd_v2r8(dy20,fscal,fjy0);
1202             fjz0             = _fjsp_madd_v2r8(dz20,fscal,fjz0);
1203
1204             }
1205
1206             /**************************
1207              * CALCULATE INTERACTIONS *
1208              **************************/
1209
1210             if (gmx_fjsp_any_lt_v2r8(rsq30,rcutoff2))
1211             {
1212
1213             r30              = _fjsp_mul_v2r8(rsq30,rinv30);
1214
1215             /* Compute parameters for interactions between i and j atoms */
1216             qq30             = _fjsp_mul_v2r8(iq3,jq0);
1217
1218             /* EWALD ELECTROSTATICS */
1219
1220             /* Calculate Ewald table index by multiplying r with scale and truncate to integer */
1221             ewrt             = _fjsp_mul_v2r8(r30,ewtabscale);
1222             itab_tmp         = _fjsp_dtox_v2r8(ewrt);
1223             eweps            = _fjsp_sub_v2r8(ewrt,_fjsp_xtod_v2r8(itab_tmp));
1224             _fjsp_store_v2r8(&ewconv.simd,itab_tmp);
1225
1226             ewtabF           = _fjsp_load_v2r8( ewtab + 4*ewconv.i[0] );
1227             ewtabD           = _fjsp_setzero_v2r8();
1228             GMX_FJSP_TRANSPOSE2_V2R8(ewtabF,ewtabD);
1229             ewtabV           = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(), ewtab + 4*ewconv.i[0] +2);
1230             ewtabFn          = _fjsp_setzero_v2r8();
1231             GMX_FJSP_TRANSPOSE2_V2R8(ewtabV,ewtabFn);
1232             felec            = _fjsp_madd_v2r8(eweps,ewtabD,ewtabF);
1233             velec            = _fjsp_nmsub_v2r8(_fjsp_mul_v2r8(ewtabhalfspace,eweps) ,_fjsp_add_v2r8(ewtabF,felec), ewtabV);
1234             velec            = _fjsp_mul_v2r8(qq30,_fjsp_sub_v2r8(rinv30,velec));
1235             felec            = _fjsp_mul_v2r8(_fjsp_mul_v2r8(qq30,rinv30),_fjsp_sub_v2r8(rinvsq30,felec));
1236
1237             d                = _fjsp_sub_v2r8(r30,rswitch);
1238             d                = _fjsp_max_v2r8(d,_fjsp_setzero_v2r8());
1239             d2               = _fjsp_mul_v2r8(d,d);
1240             sw               = _fjsp_add_v2r8(one,_fjsp_mul_v2r8(d2,_fjsp_mul_v2r8(d,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swV5,swV4),swV3))));
1241
1242             dsw              = _fjsp_mul_v2r8(d2,_fjsp_madd_v2r8(d,_fjsp_madd_v2r8(d,swF4,swF3),swF2));
1243
1244             /* Evaluate switch function */
1245             /* fscal'=f'/r=-(v*sw)'/r=-(v'*sw+v*dsw)/r=-v'*sw/r-v*dsw/r=fscal*sw-v*dsw/r */
1246             felec            = _fjsp_msub_v2r8( felec,sw , _fjsp_mul_v2r8(rinv30,_fjsp_mul_v2r8(velec,dsw)) );
1247             cutoff_mask      = _fjsp_cmplt_v2r8(rsq30,rcutoff2);
1248
1249             fscal            = felec;
1250
1251             fscal            = _fjsp_and_v2r8(fscal,cutoff_mask);
1252
1253             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
1254
1255             /* Update vectorial force */
1256             fix3             = _fjsp_madd_v2r8(dx30,fscal,fix3);
1257             fiy3             = _fjsp_madd_v2r8(dy30,fscal,fiy3);
1258             fiz3             = _fjsp_madd_v2r8(dz30,fscal,fiz3);
1259             
1260             fjx0             = _fjsp_madd_v2r8(dx30,fscal,fjx0);
1261             fjy0             = _fjsp_madd_v2r8(dy30,fscal,fjy0);
1262             fjz0             = _fjsp_madd_v2r8(dz30,fscal,fjz0);
1263
1264             }
1265
1266             gmx_fjsp_decrement_1rvec_1ptr_swizzle_v2r8(f+j_coord_offsetA,fjx0,fjy0,fjz0);
1267
1268             /* Inner loop uses 198 flops */
1269         }
1270
1271         /* End of innermost loop */
1272
1273         gmx_fjsp_update_iforce_3atom_swizzle_v2r8(fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
1274                                               f+i_coord_offset+DIM,fshift+i_shift_offset);
1275
1276         /* Increment number of inner iterations */
1277         inneriter                  += j_index_end - j_index_start;
1278
1279         /* Outer loop uses 18 flops */
1280     }
1281
1282     /* Increment number of outer iterations */
1283     outeriter        += nri;
1284
1285     /* Update outer/inner flops */
1286
1287     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W4_F,outeriter*18 + inneriter*198);
1288 }