a4746031a8fd9c9a62f991e941d3da51f74e2907
[alexxy/gromacs.git] / src / gromacs / gmxlib / nonbonded / nb_kernel_sparc64_hpc_ace_double / nb_kernel_ElecGB_VdwNone_GeomP1P1_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 #include "config.h"
39
40 #include <math.h>
41
42 #include "../nb_kernel.h"
43 #include "types/simple.h"
44 #include "gromacs/math/vec.h"
45 #include "nrnb.h"
46
47 #include "kernelutil_sparc64_hpc_ace_double.h"
48
49 /*
50  * Gromacs nonbonded kernel:   nb_kernel_ElecGB_VdwNone_GeomP1P1_VF_sparc64_hpc_ace_double
51  * Electrostatics interaction: GeneralizedBorn
52  * VdW interaction:            None
53  * Geometry:                   Particle-Particle
54  * Calculate force/pot:        PotentialAndForce
55  */
56 void
57 nb_kernel_ElecGB_VdwNone_GeomP1P1_VF_sparc64_hpc_ace_double
58                     (t_nblist                    * gmx_restrict       nlist,
59                      rvec                        * gmx_restrict          xx,
60                      rvec                        * gmx_restrict          ff,
61                      t_forcerec                  * gmx_restrict          fr,
62                      t_mdatoms                   * gmx_restrict     mdatoms,
63                      nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
64                      t_nrnb                      * gmx_restrict        nrnb)
65 {
66     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
67      * just 0 for non-waters.
68      * Suffixes A,B refer to j loop unrolling done with double precision SIMD, e.g. for the two different
69      * jnr indices corresponding to data put in the four positions in the SIMD register.
70      */
71     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
72     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
73     int              jnrA,jnrB;
74     int              j_coord_offsetA,j_coord_offsetB;
75     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
76     real             rcutoff_scalar;
77     real             *shiftvec,*fshift,*x,*f;
78     _fjsp_v2r8       tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
79     int              vdwioffset0;
80     _fjsp_v2r8       ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
81     int              vdwjidx0A,vdwjidx0B;
82     _fjsp_v2r8       jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
83     _fjsp_v2r8       dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
84     _fjsp_v2r8       velec,felec,velecsum,facel,crf,krf,krf2;
85     real             *charge;
86     _fjsp_v2r8       vgb,fgb,vgbsum,dvdasum,gbscale,gbtabscale,isaprod,gbqqfactor,gbinvepsdiff,dvdaj,gbeps,twogbeps,dvdatmp;
87     _fjsp_v2r8       minushalf = gmx_fjsp_set1_v2r8(-0.5);
88     real             *invsqrta,*dvda,*gbtab;
89     _fjsp_v2r8       rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF,twovfeps;
90     real             *vftab;
91     _fjsp_v2r8       itab_tmp;
92     _fjsp_v2r8       dummy_mask,cutoff_mask;
93     _fjsp_v2r8       one     = gmx_fjsp_set1_v2r8(1.0);
94     _fjsp_v2r8       two     = gmx_fjsp_set1_v2r8(2.0);
95     union { _fjsp_v2r8 simd; long long int i[2]; } vfconv,gbconv,ewconv;
96
97     x                = xx[0];
98     f                = ff[0];
99
100     nri              = nlist->nri;
101     iinr             = nlist->iinr;
102     jindex           = nlist->jindex;
103     jjnr             = nlist->jjnr;
104     shiftidx         = nlist->shift;
105     gid              = nlist->gid;
106     shiftvec         = fr->shift_vec[0];
107     fshift           = fr->fshift[0];
108     facel            = gmx_fjsp_set1_v2r8(fr->epsfac);
109     charge           = mdatoms->chargeA;
110
111     invsqrta         = fr->invsqrta;
112     dvda             = fr->dvda;
113     gbtabscale       = gmx_fjsp_set1_v2r8(fr->gbtab.scale);
114     gbtab            = fr->gbtab.data;
115     gbinvepsdiff     = gmx_fjsp_set1_v2r8((1.0/fr->epsilon_r) - (1.0/fr->gb_epsilon_solvent));
116
117     /* Avoid stupid compiler warnings */
118     jnrA = jnrB = 0;
119     j_coord_offsetA = 0;
120     j_coord_offsetB = 0;
121
122     outeriter        = 0;
123     inneriter        = 0;
124
125     /* Start outer loop over neighborlists */
126     for(iidx=0; iidx<nri; iidx++)
127     {
128         /* Load shift vector for this list */
129         i_shift_offset   = DIM*shiftidx[iidx];
130
131         /* Load limits for loop over neighbors */
132         j_index_start    = jindex[iidx];
133         j_index_end      = jindex[iidx+1];
134
135         /* Get outer coordinate index */
136         inr              = iinr[iidx];
137         i_coord_offset   = DIM*inr;
138
139         /* Load i particle coords and add shift vector */
140         gmx_fjsp_load_shift_and_1rvec_broadcast_v2r8(shiftvec+i_shift_offset,x+i_coord_offset,&ix0,&iy0,&iz0);
141
142         fix0             = _fjsp_setzero_v2r8();
143         fiy0             = _fjsp_setzero_v2r8();
144         fiz0             = _fjsp_setzero_v2r8();
145
146         /* Load parameters for i particles */
147         iq0              = _fjsp_mul_v2r8(facel,gmx_fjsp_load1_v2r8(charge+inr+0));
148         isai0            = gmx_fjsp_load1_v2r8(invsqrta+inr+0);
149
150         /* Reset potential sums */
151         velecsum         = _fjsp_setzero_v2r8();
152         vgbsum           = _fjsp_setzero_v2r8();
153         dvdasum          = _fjsp_setzero_v2r8();
154
155         /* Start inner kernel loop */
156         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
157         {
158
159             /* Get j neighbor index, and coordinate index */
160             jnrA             = jjnr[jidx];
161             jnrB             = jjnr[jidx+1];
162             j_coord_offsetA  = DIM*jnrA;
163             j_coord_offsetB  = DIM*jnrB;
164
165             /* load j atom coordinates */
166             gmx_fjsp_load_1rvec_2ptr_swizzle_v2r8(x+j_coord_offsetA,x+j_coord_offsetB,
167                                               &jx0,&jy0,&jz0);
168
169             /* Calculate displacement vector */
170             dx00             = _fjsp_sub_v2r8(ix0,jx0);
171             dy00             = _fjsp_sub_v2r8(iy0,jy0);
172             dz00             = _fjsp_sub_v2r8(iz0,jz0);
173
174             /* Calculate squared distance and things based on it */
175             rsq00            = gmx_fjsp_calc_rsq_v2r8(dx00,dy00,dz00);
176
177             rinv00           = gmx_fjsp_invsqrt_v2r8(rsq00);
178
179             /* Load parameters for j particles */
180             jq0              = gmx_fjsp_load_2real_swizzle_v2r8(charge+jnrA+0,charge+jnrB+0);
181             isaj0            = gmx_fjsp_load_2real_swizzle_v2r8(invsqrta+jnrA+0,invsqrta+jnrB+0);
182
183             /**************************
184              * CALCULATE INTERACTIONS *
185              **************************/
186
187             r00              = _fjsp_mul_v2r8(rsq00,rinv00);
188
189             /* Compute parameters for interactions between i and j atoms */
190             qq00             = _fjsp_mul_v2r8(iq0,jq0);
191
192             /* GENERALIZED BORN AND COULOMB ELECTROSTATICS */
193             isaprod          = _fjsp_mul_v2r8(isai0,isaj0);
194             gbqqfactor       = _fjsp_neg_v2r8(_fjsp_mul_v2r8(qq00,_fjsp_mul_v2r8(isaprod,gbinvepsdiff)));
195             gbscale          = _fjsp_mul_v2r8(isaprod,gbtabscale);
196
197             /* Calculate generalized born table index - this is a separate table from the normal one,
198              * but we use the same procedure by multiplying r with scale and truncating to integer.
199              */
200             rt               = _fjsp_mul_v2r8(r00,gbscale);
201             itab_tmp         = _fjsp_dtox_v2r8(rt);
202             gbeps            = _fjsp_sub_v2r8(rt,_fjsp_xtod_v2r8(itab_tmp));
203             _fjsp_store_v2r8(&gbconv.simd,itab_tmp);
204
205             Y                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] );
206             F                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[1] );
207             GMX_FJSP_TRANSPOSE2_V2R8(Y,F);
208             G                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] +2);
209             H                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[1] +2);
210             GMX_FJSP_TRANSPOSE2_V2R8(G,H);
211             Fp               = _fjsp_madd_v2r8(gbeps,_fjsp_madd_v2r8(gbeps,H,G),F);
212             VV               = _fjsp_madd_v2r8(gbeps,Fp,Y);
213             vgb              = _fjsp_mul_v2r8(gbqqfactor,VV);
214
215             twogbeps         = _fjsp_add_v2r8(gbeps,gbeps);
216             FF               = _fjsp_madd_v2r8(_fjsp_madd_v2r8(twogbeps,H,G),gbeps,Fp);
217             fgb              = _fjsp_mul_v2r8(gbqqfactor,_fjsp_mul_v2r8(FF,gbscale));
218             dvdatmp          = _fjsp_mul_v2r8(minushalf,_fjsp_madd_v2r8(fgb,r00,vgb));
219             dvdasum          = _fjsp_add_v2r8(dvdasum,dvdatmp);
220             gmx_fjsp_increment_2real_swizzle_v2r8(dvda+jnrA,dvda+jnrB,_fjsp_mul_v2r8(dvdatmp,_fjsp_mul_v2r8(isaj0,isaj0)));
221             velec            = _fjsp_mul_v2r8(qq00,rinv00);
222             felec            = _fjsp_mul_v2r8(_fjsp_msub_v2r8(velec,rinv00,fgb),rinv00);
223
224             /* Update potential sum for this i atom from the interaction with this j atom. */
225             velecsum         = _fjsp_add_v2r8(velecsum,velec);
226             vgbsum           = _fjsp_add_v2r8(vgbsum,vgb);
227
228             fscal            = felec;
229
230             /* Update vectorial force */
231             fix0             = _fjsp_madd_v2r8(dx00,fscal,fix0);
232             fiy0             = _fjsp_madd_v2r8(dy00,fscal,fiy0);
233             fiz0             = _fjsp_madd_v2r8(dz00,fscal,fiz0);
234             
235             gmx_fjsp_decrement_fma_1rvec_2ptr_swizzle_v2r8(f+j_coord_offsetA,f+j_coord_offsetB,fscal,dx00,dy00,dz00);
236
237             /* Inner loop uses 61 flops */
238         }
239
240         if(jidx<j_index_end)
241         {
242
243             jnrA             = jjnr[jidx];
244             j_coord_offsetA  = DIM*jnrA;
245
246             /* load j atom coordinates */
247             gmx_fjsp_load_1rvec_1ptr_swizzle_v2r8(x+j_coord_offsetA,
248                                               &jx0,&jy0,&jz0);
249
250             /* Calculate displacement vector */
251             dx00             = _fjsp_sub_v2r8(ix0,jx0);
252             dy00             = _fjsp_sub_v2r8(iy0,jy0);
253             dz00             = _fjsp_sub_v2r8(iz0,jz0);
254
255             /* Calculate squared distance and things based on it */
256             rsq00            = gmx_fjsp_calc_rsq_v2r8(dx00,dy00,dz00);
257
258             rinv00           = gmx_fjsp_invsqrt_v2r8(rsq00);
259
260             /* Load parameters for j particles */
261             jq0              = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(),charge+jnrA+0);
262             isaj0            = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(),invsqrta+jnrA+0);
263
264             /**************************
265              * CALCULATE INTERACTIONS *
266              **************************/
267
268             r00              = _fjsp_mul_v2r8(rsq00,rinv00);
269
270             /* Compute parameters for interactions between i and j atoms */
271             qq00             = _fjsp_mul_v2r8(iq0,jq0);
272
273             /* GENERALIZED BORN AND COULOMB ELECTROSTATICS */
274             isaprod          = _fjsp_mul_v2r8(isai0,isaj0);
275             gbqqfactor       = _fjsp_neg_v2r8(_fjsp_mul_v2r8(qq00,_fjsp_mul_v2r8(isaprod,gbinvepsdiff)));
276             gbscale          = _fjsp_mul_v2r8(isaprod,gbtabscale);
277
278             /* Calculate generalized born table index - this is a separate table from the normal one,
279              * but we use the same procedure by multiplying r with scale and truncating to integer.
280              */
281             rt               = _fjsp_mul_v2r8(r00,gbscale);
282             itab_tmp         = _fjsp_dtox_v2r8(rt);
283             gbeps            = _fjsp_sub_v2r8(rt,_fjsp_xtod_v2r8(itab_tmp));
284             _fjsp_store_v2r8(&gbconv.simd,itab_tmp);
285
286             Y                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] );
287             F                = _fjsp_setzero_v2r8();
288             GMX_FJSP_TRANSPOSE2_V2R8(Y,F);
289             G                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] +2);
290             H                = _fjsp_setzero_v2r8();
291             GMX_FJSP_TRANSPOSE2_V2R8(G,H);
292             Fp               = _fjsp_madd_v2r8(gbeps,_fjsp_madd_v2r8(gbeps,H,G),F);
293             VV               = _fjsp_madd_v2r8(gbeps,Fp,Y);
294             vgb              = _fjsp_mul_v2r8(gbqqfactor,VV);
295
296             twogbeps         = _fjsp_add_v2r8(gbeps,gbeps);
297             FF               = _fjsp_madd_v2r8(_fjsp_madd_v2r8(twogbeps,H,G),gbeps,Fp);
298             fgb              = _fjsp_mul_v2r8(gbqqfactor,_fjsp_mul_v2r8(FF,gbscale));
299             dvdatmp          = _fjsp_mul_v2r8(minushalf,_fjsp_madd_v2r8(fgb,r00,vgb));
300             dvdasum          = _fjsp_add_v2r8(dvdasum,dvdatmp);
301             gmx_fjsp_increment_1real_v2r8(dvda+jnrA,_fjsp_mul_v2r8(dvdatmp,_fjsp_mul_v2r8(isaj0,isaj0)));
302             velec            = _fjsp_mul_v2r8(qq00,rinv00);
303             felec            = _fjsp_mul_v2r8(_fjsp_msub_v2r8(velec,rinv00,fgb),rinv00);
304
305             /* Update potential sum for this i atom from the interaction with this j atom. */
306             velec            = _fjsp_unpacklo_v2r8(velec,_fjsp_setzero_v2r8());
307             velecsum         = _fjsp_add_v2r8(velecsum,velec);
308             vgb              = _fjsp_unpacklo_v2r8(vgb,_fjsp_setzero_v2r8());
309             vgbsum           = _fjsp_add_v2r8(vgbsum,vgb);
310
311             fscal            = felec;
312
313             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
314
315             /* Update vectorial force */
316             fix0             = _fjsp_madd_v2r8(dx00,fscal,fix0);
317             fiy0             = _fjsp_madd_v2r8(dy00,fscal,fiy0);
318             fiz0             = _fjsp_madd_v2r8(dz00,fscal,fiz0);
319             
320             gmx_fjsp_decrement_fma_1rvec_1ptr_swizzle_v2r8(f+j_coord_offsetA,fscal,dx00,dy00,dz00);
321
322             /* Inner loop uses 61 flops */
323         }
324
325         /* End of innermost loop */
326
327         gmx_fjsp_update_iforce_1atom_swizzle_v2r8(fix0,fiy0,fiz0,
328                                               f+i_coord_offset,fshift+i_shift_offset);
329
330         ggid                        = gid[iidx];
331         /* Update potential energies */
332         gmx_fjsp_update_1pot_v2r8(velecsum,kernel_data->energygrp_elec+ggid);
333         gmx_fjsp_update_1pot_v2r8(vgbsum,kernel_data->energygrp_polarization+ggid);
334         dvdasum = _fjsp_mul_v2r8(dvdasum, _fjsp_mul_v2r8(isai0,isai0));
335         gmx_fjsp_update_1pot_v2r8(dvdasum,dvda+inr);
336
337         /* Increment number of inner iterations */
338         inneriter                  += j_index_end - j_index_start;
339
340         /* Outer loop uses 9 flops */
341     }
342
343     /* Increment number of outer iterations */
344     outeriter        += nri;
345
346     /* Update outer/inner flops */
347
348     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_VF,outeriter*9 + inneriter*61);
349 }
350 /*
351  * Gromacs nonbonded kernel:   nb_kernel_ElecGB_VdwNone_GeomP1P1_F_sparc64_hpc_ace_double
352  * Electrostatics interaction: GeneralizedBorn
353  * VdW interaction:            None
354  * Geometry:                   Particle-Particle
355  * Calculate force/pot:        Force
356  */
357 void
358 nb_kernel_ElecGB_VdwNone_GeomP1P1_F_sparc64_hpc_ace_double
359                     (t_nblist                    * gmx_restrict       nlist,
360                      rvec                        * gmx_restrict          xx,
361                      rvec                        * gmx_restrict          ff,
362                      t_forcerec                  * gmx_restrict          fr,
363                      t_mdatoms                   * gmx_restrict     mdatoms,
364                      nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
365                      t_nrnb                      * gmx_restrict        nrnb)
366 {
367     /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
368      * just 0 for non-waters.
369      * Suffixes A,B refer to j loop unrolling done with double precision SIMD, e.g. for the two different
370      * jnr indices corresponding to data put in the four positions in the SIMD register.
371      */
372     int              i_shift_offset,i_coord_offset,outeriter,inneriter;
373     int              j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
374     int              jnrA,jnrB;
375     int              j_coord_offsetA,j_coord_offsetB;
376     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
377     real             rcutoff_scalar;
378     real             *shiftvec,*fshift,*x,*f;
379     _fjsp_v2r8       tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
380     int              vdwioffset0;
381     _fjsp_v2r8       ix0,iy0,iz0,fix0,fiy0,fiz0,iq0,isai0;
382     int              vdwjidx0A,vdwjidx0B;
383     _fjsp_v2r8       jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
384     _fjsp_v2r8       dx00,dy00,dz00,rsq00,rinv00,rinvsq00,r00,qq00,c6_00,c12_00;
385     _fjsp_v2r8       velec,felec,velecsum,facel,crf,krf,krf2;
386     real             *charge;
387     _fjsp_v2r8       vgb,fgb,vgbsum,dvdasum,gbscale,gbtabscale,isaprod,gbqqfactor,gbinvepsdiff,dvdaj,gbeps,twogbeps,dvdatmp;
388     _fjsp_v2r8       minushalf = gmx_fjsp_set1_v2r8(-0.5);
389     real             *invsqrta,*dvda,*gbtab;
390     _fjsp_v2r8       rt,vfeps,vftabscale,Y,F,G,H,Heps,Fp,VV,FF,twovfeps;
391     real             *vftab;
392     _fjsp_v2r8       itab_tmp;
393     _fjsp_v2r8       dummy_mask,cutoff_mask;
394     _fjsp_v2r8       one     = gmx_fjsp_set1_v2r8(1.0);
395     _fjsp_v2r8       two     = gmx_fjsp_set1_v2r8(2.0);
396     union { _fjsp_v2r8 simd; long long int i[2]; } vfconv,gbconv,ewconv;
397
398     x                = xx[0];
399     f                = ff[0];
400
401     nri              = nlist->nri;
402     iinr             = nlist->iinr;
403     jindex           = nlist->jindex;
404     jjnr             = nlist->jjnr;
405     shiftidx         = nlist->shift;
406     gid              = nlist->gid;
407     shiftvec         = fr->shift_vec[0];
408     fshift           = fr->fshift[0];
409     facel            = gmx_fjsp_set1_v2r8(fr->epsfac);
410     charge           = mdatoms->chargeA;
411
412     invsqrta         = fr->invsqrta;
413     dvda             = fr->dvda;
414     gbtabscale       = gmx_fjsp_set1_v2r8(fr->gbtab.scale);
415     gbtab            = fr->gbtab.data;
416     gbinvepsdiff     = gmx_fjsp_set1_v2r8((1.0/fr->epsilon_r) - (1.0/fr->gb_epsilon_solvent));
417
418     /* Avoid stupid compiler warnings */
419     jnrA = jnrB = 0;
420     j_coord_offsetA = 0;
421     j_coord_offsetB = 0;
422
423     outeriter        = 0;
424     inneriter        = 0;
425
426     /* Start outer loop over neighborlists */
427     for(iidx=0; iidx<nri; iidx++)
428     {
429         /* Load shift vector for this list */
430         i_shift_offset   = DIM*shiftidx[iidx];
431
432         /* Load limits for loop over neighbors */
433         j_index_start    = jindex[iidx];
434         j_index_end      = jindex[iidx+1];
435
436         /* Get outer coordinate index */
437         inr              = iinr[iidx];
438         i_coord_offset   = DIM*inr;
439
440         /* Load i particle coords and add shift vector */
441         gmx_fjsp_load_shift_and_1rvec_broadcast_v2r8(shiftvec+i_shift_offset,x+i_coord_offset,&ix0,&iy0,&iz0);
442
443         fix0             = _fjsp_setzero_v2r8();
444         fiy0             = _fjsp_setzero_v2r8();
445         fiz0             = _fjsp_setzero_v2r8();
446
447         /* Load parameters for i particles */
448         iq0              = _fjsp_mul_v2r8(facel,gmx_fjsp_load1_v2r8(charge+inr+0));
449         isai0            = gmx_fjsp_load1_v2r8(invsqrta+inr+0);
450
451         dvdasum          = _fjsp_setzero_v2r8();
452
453         /* Start inner kernel loop */
454         for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
455         {
456
457             /* Get j neighbor index, and coordinate index */
458             jnrA             = jjnr[jidx];
459             jnrB             = jjnr[jidx+1];
460             j_coord_offsetA  = DIM*jnrA;
461             j_coord_offsetB  = DIM*jnrB;
462
463             /* load j atom coordinates */
464             gmx_fjsp_load_1rvec_2ptr_swizzle_v2r8(x+j_coord_offsetA,x+j_coord_offsetB,
465                                               &jx0,&jy0,&jz0);
466
467             /* Calculate displacement vector */
468             dx00             = _fjsp_sub_v2r8(ix0,jx0);
469             dy00             = _fjsp_sub_v2r8(iy0,jy0);
470             dz00             = _fjsp_sub_v2r8(iz0,jz0);
471
472             /* Calculate squared distance and things based on it */
473             rsq00            = gmx_fjsp_calc_rsq_v2r8(dx00,dy00,dz00);
474
475             rinv00           = gmx_fjsp_invsqrt_v2r8(rsq00);
476
477             /* Load parameters for j particles */
478             jq0              = gmx_fjsp_load_2real_swizzle_v2r8(charge+jnrA+0,charge+jnrB+0);
479             isaj0            = gmx_fjsp_load_2real_swizzle_v2r8(invsqrta+jnrA+0,invsqrta+jnrB+0);
480
481             /**************************
482              * CALCULATE INTERACTIONS *
483              **************************/
484
485             r00              = _fjsp_mul_v2r8(rsq00,rinv00);
486
487             /* Compute parameters for interactions between i and j atoms */
488             qq00             = _fjsp_mul_v2r8(iq0,jq0);
489
490             /* GENERALIZED BORN AND COULOMB ELECTROSTATICS */
491             isaprod          = _fjsp_mul_v2r8(isai0,isaj0);
492             gbqqfactor       = _fjsp_neg_v2r8(_fjsp_mul_v2r8(qq00,_fjsp_mul_v2r8(isaprod,gbinvepsdiff)));
493             gbscale          = _fjsp_mul_v2r8(isaprod,gbtabscale);
494
495             /* Calculate generalized born table index - this is a separate table from the normal one,
496              * but we use the same procedure by multiplying r with scale and truncating to integer.
497              */
498             rt               = _fjsp_mul_v2r8(r00,gbscale);
499             itab_tmp         = _fjsp_dtox_v2r8(rt);
500             gbeps            = _fjsp_sub_v2r8(rt,_fjsp_xtod_v2r8(itab_tmp));
501             _fjsp_store_v2r8(&gbconv.simd,itab_tmp);
502
503             Y                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] );
504             F                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[1] );
505             GMX_FJSP_TRANSPOSE2_V2R8(Y,F);
506             G                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] +2);
507             H                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[1] +2);
508             GMX_FJSP_TRANSPOSE2_V2R8(G,H);
509             Fp               = _fjsp_madd_v2r8(gbeps,_fjsp_madd_v2r8(gbeps,H,G),F);
510             VV               = _fjsp_madd_v2r8(gbeps,Fp,Y);
511             vgb              = _fjsp_mul_v2r8(gbqqfactor,VV);
512
513             twogbeps         = _fjsp_add_v2r8(gbeps,gbeps);
514             FF               = _fjsp_madd_v2r8(_fjsp_madd_v2r8(twogbeps,H,G),gbeps,Fp);
515             fgb              = _fjsp_mul_v2r8(gbqqfactor,_fjsp_mul_v2r8(FF,gbscale));
516             dvdatmp          = _fjsp_mul_v2r8(minushalf,_fjsp_madd_v2r8(fgb,r00,vgb));
517             dvdasum          = _fjsp_add_v2r8(dvdasum,dvdatmp);
518             gmx_fjsp_increment_2real_swizzle_v2r8(dvda+jnrA,dvda+jnrB,_fjsp_mul_v2r8(dvdatmp,_fjsp_mul_v2r8(isaj0,isaj0)));
519             velec            = _fjsp_mul_v2r8(qq00,rinv00);
520             felec            = _fjsp_mul_v2r8(_fjsp_msub_v2r8(velec,rinv00,fgb),rinv00);
521
522             fscal            = felec;
523
524             /* Update vectorial force */
525             fix0             = _fjsp_madd_v2r8(dx00,fscal,fix0);
526             fiy0             = _fjsp_madd_v2r8(dy00,fscal,fiy0);
527             fiz0             = _fjsp_madd_v2r8(dz00,fscal,fiz0);
528             
529             gmx_fjsp_decrement_fma_1rvec_2ptr_swizzle_v2r8(f+j_coord_offsetA,f+j_coord_offsetB,fscal,dx00,dy00,dz00);
530
531             /* Inner loop uses 59 flops */
532         }
533
534         if(jidx<j_index_end)
535         {
536
537             jnrA             = jjnr[jidx];
538             j_coord_offsetA  = DIM*jnrA;
539
540             /* load j atom coordinates */
541             gmx_fjsp_load_1rvec_1ptr_swizzle_v2r8(x+j_coord_offsetA,
542                                               &jx0,&jy0,&jz0);
543
544             /* Calculate displacement vector */
545             dx00             = _fjsp_sub_v2r8(ix0,jx0);
546             dy00             = _fjsp_sub_v2r8(iy0,jy0);
547             dz00             = _fjsp_sub_v2r8(iz0,jz0);
548
549             /* Calculate squared distance and things based on it */
550             rsq00            = gmx_fjsp_calc_rsq_v2r8(dx00,dy00,dz00);
551
552             rinv00           = gmx_fjsp_invsqrt_v2r8(rsq00);
553
554             /* Load parameters for j particles */
555             jq0              = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(),charge+jnrA+0);
556             isaj0            = _fjsp_loadl_v2r8(_fjsp_setzero_v2r8(),invsqrta+jnrA+0);
557
558             /**************************
559              * CALCULATE INTERACTIONS *
560              **************************/
561
562             r00              = _fjsp_mul_v2r8(rsq00,rinv00);
563
564             /* Compute parameters for interactions between i and j atoms */
565             qq00             = _fjsp_mul_v2r8(iq0,jq0);
566
567             /* GENERALIZED BORN AND COULOMB ELECTROSTATICS */
568             isaprod          = _fjsp_mul_v2r8(isai0,isaj0);
569             gbqqfactor       = _fjsp_neg_v2r8(_fjsp_mul_v2r8(qq00,_fjsp_mul_v2r8(isaprod,gbinvepsdiff)));
570             gbscale          = _fjsp_mul_v2r8(isaprod,gbtabscale);
571
572             /* Calculate generalized born table index - this is a separate table from the normal one,
573              * but we use the same procedure by multiplying r with scale and truncating to integer.
574              */
575             rt               = _fjsp_mul_v2r8(r00,gbscale);
576             itab_tmp         = _fjsp_dtox_v2r8(rt);
577             gbeps            = _fjsp_sub_v2r8(rt,_fjsp_xtod_v2r8(itab_tmp));
578             _fjsp_store_v2r8(&gbconv.simd,itab_tmp);
579
580             Y                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] );
581             F                = _fjsp_setzero_v2r8();
582             GMX_FJSP_TRANSPOSE2_V2R8(Y,F);
583             G                = _fjsp_load_v2r8( gbtab + 4*gbconv.i[0] +2);
584             H                = _fjsp_setzero_v2r8();
585             GMX_FJSP_TRANSPOSE2_V2R8(G,H);
586             Fp               = _fjsp_madd_v2r8(gbeps,_fjsp_madd_v2r8(gbeps,H,G),F);
587             VV               = _fjsp_madd_v2r8(gbeps,Fp,Y);
588             vgb              = _fjsp_mul_v2r8(gbqqfactor,VV);
589
590             twogbeps         = _fjsp_add_v2r8(gbeps,gbeps);
591             FF               = _fjsp_madd_v2r8(_fjsp_madd_v2r8(twogbeps,H,G),gbeps,Fp);
592             fgb              = _fjsp_mul_v2r8(gbqqfactor,_fjsp_mul_v2r8(FF,gbscale));
593             dvdatmp          = _fjsp_mul_v2r8(minushalf,_fjsp_madd_v2r8(fgb,r00,vgb));
594             dvdasum          = _fjsp_add_v2r8(dvdasum,dvdatmp);
595             gmx_fjsp_increment_1real_v2r8(dvda+jnrA,_fjsp_mul_v2r8(dvdatmp,_fjsp_mul_v2r8(isaj0,isaj0)));
596             velec            = _fjsp_mul_v2r8(qq00,rinv00);
597             felec            = _fjsp_mul_v2r8(_fjsp_msub_v2r8(velec,rinv00,fgb),rinv00);
598
599             fscal            = felec;
600
601             fscal            = _fjsp_unpacklo_v2r8(fscal,_fjsp_setzero_v2r8());
602
603             /* Update vectorial force */
604             fix0             = _fjsp_madd_v2r8(dx00,fscal,fix0);
605             fiy0             = _fjsp_madd_v2r8(dy00,fscal,fiy0);
606             fiz0             = _fjsp_madd_v2r8(dz00,fscal,fiz0);
607             
608             gmx_fjsp_decrement_fma_1rvec_1ptr_swizzle_v2r8(f+j_coord_offsetA,fscal,dx00,dy00,dz00);
609
610             /* Inner loop uses 59 flops */
611         }
612
613         /* End of innermost loop */
614
615         gmx_fjsp_update_iforce_1atom_swizzle_v2r8(fix0,fiy0,fiz0,
616                                               f+i_coord_offset,fshift+i_shift_offset);
617
618         dvdasum = _fjsp_mul_v2r8(dvdasum, _fjsp_mul_v2r8(isai0,isai0));
619         gmx_fjsp_update_1pot_v2r8(dvdasum,dvda+inr);
620
621         /* Increment number of inner iterations */
622         inneriter                  += j_index_end - j_index_start;
623
624         /* Outer loop uses 7 flops */
625     }
626
627     /* Increment number of outer iterations */
628     outeriter        += nri;
629
630     /* Update outer/inner flops */
631
632     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_F,outeriter*7 + inneriter*59);
633 }