2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014,2015, 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.
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.
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.
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.
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.
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.
36 * Note: this file was generated by the GROMACS sse2_double kernel generator.
44 #include "../nb_kernel.h"
45 #include "gromacs/math/vec.h"
46 #include "gromacs/legacyheaders/nrnb.h"
48 #include "gromacs/simd/math_x86_sse2_double.h"
49 #include "kernelutil_x86_sse2_double.h"
52 * Gromacs nonbonded kernel: nb_kernel_ElecCoul_VdwLJ_GeomW3P1_VF_sse2_double
53 * Electrostatics interaction: Coulomb
54 * VdW interaction: LennardJones
55 * Geometry: Water3-Particle
56 * Calculate force/pot: PotentialAndForce
59 nb_kernel_ElecCoul_VdwLJ_GeomW3P1_VF_sse2_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
)
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 SSE double precision, e.g. for the two different
71 * jnr indices corresponding to data put in the four positions in the SIMD register.
73 int i_shift_offset
,i_coord_offset
,outeriter
,inneriter
;
74 int j_index_start
,j_index_end
,jidx
,nri
,inr
,ggid
,iidx
;
76 int j_coord_offsetA
,j_coord_offsetB
;
77 int *iinr
,*jindex
,*jjnr
,*shiftidx
,*gid
;
79 real
*shiftvec
,*fshift
,*x
,*f
;
80 __m128d tx
,ty
,tz
,fscal
,rcutoff
,rcutoff2
,jidxall
;
82 __m128d ix0
,iy0
,iz0
,fix0
,fiy0
,fiz0
,iq0
,isai0
;
84 __m128d ix1
,iy1
,iz1
,fix1
,fiy1
,fiz1
,iq1
,isai1
;
86 __m128d ix2
,iy2
,iz2
,fix2
,fiy2
,fiz2
,iq2
,isai2
;
87 int vdwjidx0A
,vdwjidx0B
;
88 __m128d jx0
,jy0
,jz0
,fjx0
,fjy0
,fjz0
,jq0
,isaj0
;
89 __m128d dx00
,dy00
,dz00
,rsq00
,rinv00
,rinvsq00
,r00
,qq00
,c6_00
,c12_00
;
90 __m128d dx10
,dy10
,dz10
,rsq10
,rinv10
,rinvsq10
,r10
,qq10
,c6_10
,c12_10
;
91 __m128d dx20
,dy20
,dz20
,rsq20
,rinv20
,rinvsq20
,r20
,qq20
,c6_20
,c12_20
;
92 __m128d velec
,felec
,velecsum
,facel
,crf
,krf
,krf2
;
95 __m128d rinvsix
,rvdw
,vvdw
,vvdw6
,vvdw12
,fvdw
,fvdw6
,fvdw12
,vvdwsum
,sh_vdw_invrcut6
;
98 __m128d one_sixth
= _mm_set1_pd(1.0/6.0);
99 __m128d one_twelfth
= _mm_set1_pd(1.0/12.0);
100 __m128d dummy_mask
,cutoff_mask
;
101 __m128d signbit
= gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
102 __m128d one
= _mm_set1_pd(1.0);
103 __m128d two
= _mm_set1_pd(2.0);
109 jindex
= nlist
->jindex
;
111 shiftidx
= nlist
->shift
;
113 shiftvec
= fr
->shift_vec
[0];
114 fshift
= fr
->fshift
[0];
115 facel
= _mm_set1_pd(fr
->epsfac
);
116 charge
= mdatoms
->chargeA
;
117 nvdwtype
= fr
->ntype
;
119 vdwtype
= mdatoms
->typeA
;
121 /* Setup water-specific parameters */
122 inr
= nlist
->iinr
[0];
123 iq0
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+0]));
124 iq1
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+1]));
125 iq2
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+2]));
126 vdwioffset0
= 2*nvdwtype
*vdwtype
[inr
+0];
128 /* Avoid stupid compiler warnings */
136 /* Start outer loop over neighborlists */
137 for(iidx
=0; iidx
<nri
; iidx
++)
139 /* Load shift vector for this list */
140 i_shift_offset
= DIM
*shiftidx
[iidx
];
142 /* Load limits for loop over neighbors */
143 j_index_start
= jindex
[iidx
];
144 j_index_end
= jindex
[iidx
+1];
146 /* Get outer coordinate index */
148 i_coord_offset
= DIM
*inr
;
150 /* Load i particle coords and add shift vector */
151 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec
+i_shift_offset
,x
+i_coord_offset
,
152 &ix0
,&iy0
,&iz0
,&ix1
,&iy1
,&iz1
,&ix2
,&iy2
,&iz2
);
154 fix0
= _mm_setzero_pd();
155 fiy0
= _mm_setzero_pd();
156 fiz0
= _mm_setzero_pd();
157 fix1
= _mm_setzero_pd();
158 fiy1
= _mm_setzero_pd();
159 fiz1
= _mm_setzero_pd();
160 fix2
= _mm_setzero_pd();
161 fiy2
= _mm_setzero_pd();
162 fiz2
= _mm_setzero_pd();
164 /* Reset potential sums */
165 velecsum
= _mm_setzero_pd();
166 vvdwsum
= _mm_setzero_pd();
168 /* Start inner kernel loop */
169 for(jidx
=j_index_start
; jidx
<j_index_end
-1; jidx
+=2)
172 /* Get j neighbor index, and coordinate index */
175 j_coord_offsetA
= DIM
*jnrA
;
176 j_coord_offsetB
= DIM
*jnrB
;
178 /* load j atom coordinates */
179 gmx_mm_load_1rvec_2ptr_swizzle_pd(x
+j_coord_offsetA
,x
+j_coord_offsetB
,
182 /* Calculate displacement vector */
183 dx00
= _mm_sub_pd(ix0
,jx0
);
184 dy00
= _mm_sub_pd(iy0
,jy0
);
185 dz00
= _mm_sub_pd(iz0
,jz0
);
186 dx10
= _mm_sub_pd(ix1
,jx0
);
187 dy10
= _mm_sub_pd(iy1
,jy0
);
188 dz10
= _mm_sub_pd(iz1
,jz0
);
189 dx20
= _mm_sub_pd(ix2
,jx0
);
190 dy20
= _mm_sub_pd(iy2
,jy0
);
191 dz20
= _mm_sub_pd(iz2
,jz0
);
193 /* Calculate squared distance and things based on it */
194 rsq00
= gmx_mm_calc_rsq_pd(dx00
,dy00
,dz00
);
195 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
196 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
198 rinv00
= gmx_mm_invsqrt_pd(rsq00
);
199 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
200 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
202 rinvsq00
= _mm_mul_pd(rinv00
,rinv00
);
203 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
204 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
206 /* Load parameters for j particles */
207 jq0
= gmx_mm_load_2real_swizzle_pd(charge
+jnrA
+0,charge
+jnrB
+0);
208 vdwjidx0A
= 2*vdwtype
[jnrA
+0];
209 vdwjidx0B
= 2*vdwtype
[jnrB
+0];
211 fjx0
= _mm_setzero_pd();
212 fjy0
= _mm_setzero_pd();
213 fjz0
= _mm_setzero_pd();
215 /**************************
216 * CALCULATE INTERACTIONS *
217 **************************/
219 /* Compute parameters for interactions between i and j atoms */
220 qq00
= _mm_mul_pd(iq0
,jq0
);
221 gmx_mm_load_2pair_swizzle_pd(vdwparam
+vdwioffset0
+vdwjidx0A
,
222 vdwparam
+vdwioffset0
+vdwjidx0B
,&c6_00
,&c12_00
);
224 /* COULOMB ELECTROSTATICS */
225 velec
= _mm_mul_pd(qq00
,rinv00
);
226 felec
= _mm_mul_pd(velec
,rinvsq00
);
228 /* LENNARD-JONES DISPERSION/REPULSION */
230 rinvsix
= _mm_mul_pd(_mm_mul_pd(rinvsq00
,rinvsq00
),rinvsq00
);
231 vvdw6
= _mm_mul_pd(c6_00
,rinvsix
);
232 vvdw12
= _mm_mul_pd(c12_00
,_mm_mul_pd(rinvsix
,rinvsix
));
233 vvdw
= _mm_sub_pd( _mm_mul_pd(vvdw12
,one_twelfth
) , _mm_mul_pd(vvdw6
,one_sixth
) );
234 fvdw
= _mm_mul_pd(_mm_sub_pd(vvdw12
,vvdw6
),rinvsq00
);
236 /* Update potential sum for this i atom from the interaction with this j atom. */
237 velecsum
= _mm_add_pd(velecsum
,velec
);
238 vvdwsum
= _mm_add_pd(vvdwsum
,vvdw
);
240 fscal
= _mm_add_pd(felec
,fvdw
);
242 /* Calculate temporary vectorial force */
243 tx
= _mm_mul_pd(fscal
,dx00
);
244 ty
= _mm_mul_pd(fscal
,dy00
);
245 tz
= _mm_mul_pd(fscal
,dz00
);
247 /* Update vectorial force */
248 fix0
= _mm_add_pd(fix0
,tx
);
249 fiy0
= _mm_add_pd(fiy0
,ty
);
250 fiz0
= _mm_add_pd(fiz0
,tz
);
252 fjx0
= _mm_add_pd(fjx0
,tx
);
253 fjy0
= _mm_add_pd(fjy0
,ty
);
254 fjz0
= _mm_add_pd(fjz0
,tz
);
256 /**************************
257 * CALCULATE INTERACTIONS *
258 **************************/
260 /* Compute parameters for interactions between i and j atoms */
261 qq10
= _mm_mul_pd(iq1
,jq0
);
263 /* COULOMB ELECTROSTATICS */
264 velec
= _mm_mul_pd(qq10
,rinv10
);
265 felec
= _mm_mul_pd(velec
,rinvsq10
);
267 /* Update potential sum for this i atom from the interaction with this j atom. */
268 velecsum
= _mm_add_pd(velecsum
,velec
);
272 /* Calculate temporary vectorial force */
273 tx
= _mm_mul_pd(fscal
,dx10
);
274 ty
= _mm_mul_pd(fscal
,dy10
);
275 tz
= _mm_mul_pd(fscal
,dz10
);
277 /* Update vectorial force */
278 fix1
= _mm_add_pd(fix1
,tx
);
279 fiy1
= _mm_add_pd(fiy1
,ty
);
280 fiz1
= _mm_add_pd(fiz1
,tz
);
282 fjx0
= _mm_add_pd(fjx0
,tx
);
283 fjy0
= _mm_add_pd(fjy0
,ty
);
284 fjz0
= _mm_add_pd(fjz0
,tz
);
286 /**************************
287 * CALCULATE INTERACTIONS *
288 **************************/
290 /* Compute parameters for interactions between i and j atoms */
291 qq20
= _mm_mul_pd(iq2
,jq0
);
293 /* COULOMB ELECTROSTATICS */
294 velec
= _mm_mul_pd(qq20
,rinv20
);
295 felec
= _mm_mul_pd(velec
,rinvsq20
);
297 /* Update potential sum for this i atom from the interaction with this j atom. */
298 velecsum
= _mm_add_pd(velecsum
,velec
);
302 /* Calculate temporary vectorial force */
303 tx
= _mm_mul_pd(fscal
,dx20
);
304 ty
= _mm_mul_pd(fscal
,dy20
);
305 tz
= _mm_mul_pd(fscal
,dz20
);
307 /* Update vectorial force */
308 fix2
= _mm_add_pd(fix2
,tx
);
309 fiy2
= _mm_add_pd(fiy2
,ty
);
310 fiz2
= _mm_add_pd(fiz2
,tz
);
312 fjx0
= _mm_add_pd(fjx0
,tx
);
313 fjy0
= _mm_add_pd(fjy0
,ty
);
314 fjz0
= _mm_add_pd(fjz0
,tz
);
316 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f
+j_coord_offsetA
,f
+j_coord_offsetB
,fjx0
,fjy0
,fjz0
);
318 /* Inner loop uses 99 flops */
325 j_coord_offsetA
= DIM
*jnrA
;
327 /* load j atom coordinates */
328 gmx_mm_load_1rvec_1ptr_swizzle_pd(x
+j_coord_offsetA
,
331 /* Calculate displacement vector */
332 dx00
= _mm_sub_pd(ix0
,jx0
);
333 dy00
= _mm_sub_pd(iy0
,jy0
);
334 dz00
= _mm_sub_pd(iz0
,jz0
);
335 dx10
= _mm_sub_pd(ix1
,jx0
);
336 dy10
= _mm_sub_pd(iy1
,jy0
);
337 dz10
= _mm_sub_pd(iz1
,jz0
);
338 dx20
= _mm_sub_pd(ix2
,jx0
);
339 dy20
= _mm_sub_pd(iy2
,jy0
);
340 dz20
= _mm_sub_pd(iz2
,jz0
);
342 /* Calculate squared distance and things based on it */
343 rsq00
= gmx_mm_calc_rsq_pd(dx00
,dy00
,dz00
);
344 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
345 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
347 rinv00
= gmx_mm_invsqrt_pd(rsq00
);
348 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
349 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
351 rinvsq00
= _mm_mul_pd(rinv00
,rinv00
);
352 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
353 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
355 /* Load parameters for j particles */
356 jq0
= _mm_load_sd(charge
+jnrA
+0);
357 vdwjidx0A
= 2*vdwtype
[jnrA
+0];
359 fjx0
= _mm_setzero_pd();
360 fjy0
= _mm_setzero_pd();
361 fjz0
= _mm_setzero_pd();
363 /**************************
364 * CALCULATE INTERACTIONS *
365 **************************/
367 /* Compute parameters for interactions between i and j atoms */
368 qq00
= _mm_mul_pd(iq0
,jq0
);
369 gmx_mm_load_1pair_swizzle_pd(vdwparam
+vdwioffset0
+vdwjidx0A
,&c6_00
,&c12_00
);
371 /* COULOMB ELECTROSTATICS */
372 velec
= _mm_mul_pd(qq00
,rinv00
);
373 felec
= _mm_mul_pd(velec
,rinvsq00
);
375 /* LENNARD-JONES DISPERSION/REPULSION */
377 rinvsix
= _mm_mul_pd(_mm_mul_pd(rinvsq00
,rinvsq00
),rinvsq00
);
378 vvdw6
= _mm_mul_pd(c6_00
,rinvsix
);
379 vvdw12
= _mm_mul_pd(c12_00
,_mm_mul_pd(rinvsix
,rinvsix
));
380 vvdw
= _mm_sub_pd( _mm_mul_pd(vvdw12
,one_twelfth
) , _mm_mul_pd(vvdw6
,one_sixth
) );
381 fvdw
= _mm_mul_pd(_mm_sub_pd(vvdw12
,vvdw6
),rinvsq00
);
383 /* Update potential sum for this i atom from the interaction with this j atom. */
384 velec
= _mm_unpacklo_pd(velec
,_mm_setzero_pd());
385 velecsum
= _mm_add_pd(velecsum
,velec
);
386 vvdw
= _mm_unpacklo_pd(vvdw
,_mm_setzero_pd());
387 vvdwsum
= _mm_add_pd(vvdwsum
,vvdw
);
389 fscal
= _mm_add_pd(felec
,fvdw
);
391 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
393 /* Calculate temporary vectorial force */
394 tx
= _mm_mul_pd(fscal
,dx00
);
395 ty
= _mm_mul_pd(fscal
,dy00
);
396 tz
= _mm_mul_pd(fscal
,dz00
);
398 /* Update vectorial force */
399 fix0
= _mm_add_pd(fix0
,tx
);
400 fiy0
= _mm_add_pd(fiy0
,ty
);
401 fiz0
= _mm_add_pd(fiz0
,tz
);
403 fjx0
= _mm_add_pd(fjx0
,tx
);
404 fjy0
= _mm_add_pd(fjy0
,ty
);
405 fjz0
= _mm_add_pd(fjz0
,tz
);
407 /**************************
408 * CALCULATE INTERACTIONS *
409 **************************/
411 /* Compute parameters for interactions between i and j atoms */
412 qq10
= _mm_mul_pd(iq1
,jq0
);
414 /* COULOMB ELECTROSTATICS */
415 velec
= _mm_mul_pd(qq10
,rinv10
);
416 felec
= _mm_mul_pd(velec
,rinvsq10
);
418 /* Update potential sum for this i atom from the interaction with this j atom. */
419 velec
= _mm_unpacklo_pd(velec
,_mm_setzero_pd());
420 velecsum
= _mm_add_pd(velecsum
,velec
);
424 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
426 /* Calculate temporary vectorial force */
427 tx
= _mm_mul_pd(fscal
,dx10
);
428 ty
= _mm_mul_pd(fscal
,dy10
);
429 tz
= _mm_mul_pd(fscal
,dz10
);
431 /* Update vectorial force */
432 fix1
= _mm_add_pd(fix1
,tx
);
433 fiy1
= _mm_add_pd(fiy1
,ty
);
434 fiz1
= _mm_add_pd(fiz1
,tz
);
436 fjx0
= _mm_add_pd(fjx0
,tx
);
437 fjy0
= _mm_add_pd(fjy0
,ty
);
438 fjz0
= _mm_add_pd(fjz0
,tz
);
440 /**************************
441 * CALCULATE INTERACTIONS *
442 **************************/
444 /* Compute parameters for interactions between i and j atoms */
445 qq20
= _mm_mul_pd(iq2
,jq0
);
447 /* COULOMB ELECTROSTATICS */
448 velec
= _mm_mul_pd(qq20
,rinv20
);
449 felec
= _mm_mul_pd(velec
,rinvsq20
);
451 /* Update potential sum for this i atom from the interaction with this j atom. */
452 velec
= _mm_unpacklo_pd(velec
,_mm_setzero_pd());
453 velecsum
= _mm_add_pd(velecsum
,velec
);
457 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
459 /* Calculate temporary vectorial force */
460 tx
= _mm_mul_pd(fscal
,dx20
);
461 ty
= _mm_mul_pd(fscal
,dy20
);
462 tz
= _mm_mul_pd(fscal
,dz20
);
464 /* Update vectorial force */
465 fix2
= _mm_add_pd(fix2
,tx
);
466 fiy2
= _mm_add_pd(fiy2
,ty
);
467 fiz2
= _mm_add_pd(fiz2
,tz
);
469 fjx0
= _mm_add_pd(fjx0
,tx
);
470 fjy0
= _mm_add_pd(fjy0
,ty
);
471 fjz0
= _mm_add_pd(fjz0
,tz
);
473 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f
+j_coord_offsetA
,fjx0
,fjy0
,fjz0
);
475 /* Inner loop uses 99 flops */
478 /* End of innermost loop */
480 gmx_mm_update_iforce_3atom_swizzle_pd(fix0
,fiy0
,fiz0
,fix1
,fiy1
,fiz1
,fix2
,fiy2
,fiz2
,
481 f
+i_coord_offset
,fshift
+i_shift_offset
);
484 /* Update potential energies */
485 gmx_mm_update_1pot_pd(velecsum
,kernel_data
->energygrp_elec
+ggid
);
486 gmx_mm_update_1pot_pd(vvdwsum
,kernel_data
->energygrp_vdw
+ggid
);
488 /* Increment number of inner iterations */
489 inneriter
+= j_index_end
- j_index_start
;
491 /* Outer loop uses 20 flops */
494 /* Increment number of outer iterations */
497 /* Update outer/inner flops */
499 inc_nrnb(nrnb
,eNR_NBKERNEL_ELEC_VDW_W3_VF
,outeriter
*20 + inneriter
*99);
502 * Gromacs nonbonded kernel: nb_kernel_ElecCoul_VdwLJ_GeomW3P1_F_sse2_double
503 * Electrostatics interaction: Coulomb
504 * VdW interaction: LennardJones
505 * Geometry: Water3-Particle
506 * Calculate force/pot: Force
509 nb_kernel_ElecCoul_VdwLJ_GeomW3P1_F_sse2_double
510 (t_nblist
* gmx_restrict nlist
,
511 rvec
* gmx_restrict xx
,
512 rvec
* gmx_restrict ff
,
513 t_forcerec
* gmx_restrict fr
,
514 t_mdatoms
* gmx_restrict mdatoms
,
515 nb_kernel_data_t gmx_unused
* gmx_restrict kernel_data
,
516 t_nrnb
* gmx_restrict nrnb
)
518 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
519 * just 0 for non-waters.
520 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
521 * jnr indices corresponding to data put in the four positions in the SIMD register.
523 int i_shift_offset
,i_coord_offset
,outeriter
,inneriter
;
524 int j_index_start
,j_index_end
,jidx
,nri
,inr
,ggid
,iidx
;
526 int j_coord_offsetA
,j_coord_offsetB
;
527 int *iinr
,*jindex
,*jjnr
,*shiftidx
,*gid
;
529 real
*shiftvec
,*fshift
,*x
,*f
;
530 __m128d tx
,ty
,tz
,fscal
,rcutoff
,rcutoff2
,jidxall
;
532 __m128d ix0
,iy0
,iz0
,fix0
,fiy0
,fiz0
,iq0
,isai0
;
534 __m128d ix1
,iy1
,iz1
,fix1
,fiy1
,fiz1
,iq1
,isai1
;
536 __m128d ix2
,iy2
,iz2
,fix2
,fiy2
,fiz2
,iq2
,isai2
;
537 int vdwjidx0A
,vdwjidx0B
;
538 __m128d jx0
,jy0
,jz0
,fjx0
,fjy0
,fjz0
,jq0
,isaj0
;
539 __m128d dx00
,dy00
,dz00
,rsq00
,rinv00
,rinvsq00
,r00
,qq00
,c6_00
,c12_00
;
540 __m128d dx10
,dy10
,dz10
,rsq10
,rinv10
,rinvsq10
,r10
,qq10
,c6_10
,c12_10
;
541 __m128d dx20
,dy20
,dz20
,rsq20
,rinv20
,rinvsq20
,r20
,qq20
,c6_20
,c12_20
;
542 __m128d velec
,felec
,velecsum
,facel
,crf
,krf
,krf2
;
545 __m128d rinvsix
,rvdw
,vvdw
,vvdw6
,vvdw12
,fvdw
,fvdw6
,fvdw12
,vvdwsum
,sh_vdw_invrcut6
;
548 __m128d one_sixth
= _mm_set1_pd(1.0/6.0);
549 __m128d one_twelfth
= _mm_set1_pd(1.0/12.0);
550 __m128d dummy_mask
,cutoff_mask
;
551 __m128d signbit
= gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
552 __m128d one
= _mm_set1_pd(1.0);
553 __m128d two
= _mm_set1_pd(2.0);
559 jindex
= nlist
->jindex
;
561 shiftidx
= nlist
->shift
;
563 shiftvec
= fr
->shift_vec
[0];
564 fshift
= fr
->fshift
[0];
565 facel
= _mm_set1_pd(fr
->epsfac
);
566 charge
= mdatoms
->chargeA
;
567 nvdwtype
= fr
->ntype
;
569 vdwtype
= mdatoms
->typeA
;
571 /* Setup water-specific parameters */
572 inr
= nlist
->iinr
[0];
573 iq0
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+0]));
574 iq1
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+1]));
575 iq2
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+2]));
576 vdwioffset0
= 2*nvdwtype
*vdwtype
[inr
+0];
578 /* Avoid stupid compiler warnings */
586 /* Start outer loop over neighborlists */
587 for(iidx
=0; iidx
<nri
; iidx
++)
589 /* Load shift vector for this list */
590 i_shift_offset
= DIM
*shiftidx
[iidx
];
592 /* Load limits for loop over neighbors */
593 j_index_start
= jindex
[iidx
];
594 j_index_end
= jindex
[iidx
+1];
596 /* Get outer coordinate index */
598 i_coord_offset
= DIM
*inr
;
600 /* Load i particle coords and add shift vector */
601 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec
+i_shift_offset
,x
+i_coord_offset
,
602 &ix0
,&iy0
,&iz0
,&ix1
,&iy1
,&iz1
,&ix2
,&iy2
,&iz2
);
604 fix0
= _mm_setzero_pd();
605 fiy0
= _mm_setzero_pd();
606 fiz0
= _mm_setzero_pd();
607 fix1
= _mm_setzero_pd();
608 fiy1
= _mm_setzero_pd();
609 fiz1
= _mm_setzero_pd();
610 fix2
= _mm_setzero_pd();
611 fiy2
= _mm_setzero_pd();
612 fiz2
= _mm_setzero_pd();
614 /* Start inner kernel loop */
615 for(jidx
=j_index_start
; jidx
<j_index_end
-1; jidx
+=2)
618 /* Get j neighbor index, and coordinate index */
621 j_coord_offsetA
= DIM
*jnrA
;
622 j_coord_offsetB
= DIM
*jnrB
;
624 /* load j atom coordinates */
625 gmx_mm_load_1rvec_2ptr_swizzle_pd(x
+j_coord_offsetA
,x
+j_coord_offsetB
,
628 /* Calculate displacement vector */
629 dx00
= _mm_sub_pd(ix0
,jx0
);
630 dy00
= _mm_sub_pd(iy0
,jy0
);
631 dz00
= _mm_sub_pd(iz0
,jz0
);
632 dx10
= _mm_sub_pd(ix1
,jx0
);
633 dy10
= _mm_sub_pd(iy1
,jy0
);
634 dz10
= _mm_sub_pd(iz1
,jz0
);
635 dx20
= _mm_sub_pd(ix2
,jx0
);
636 dy20
= _mm_sub_pd(iy2
,jy0
);
637 dz20
= _mm_sub_pd(iz2
,jz0
);
639 /* Calculate squared distance and things based on it */
640 rsq00
= gmx_mm_calc_rsq_pd(dx00
,dy00
,dz00
);
641 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
642 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
644 rinv00
= gmx_mm_invsqrt_pd(rsq00
);
645 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
646 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
648 rinvsq00
= _mm_mul_pd(rinv00
,rinv00
);
649 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
650 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
652 /* Load parameters for j particles */
653 jq0
= gmx_mm_load_2real_swizzle_pd(charge
+jnrA
+0,charge
+jnrB
+0);
654 vdwjidx0A
= 2*vdwtype
[jnrA
+0];
655 vdwjidx0B
= 2*vdwtype
[jnrB
+0];
657 fjx0
= _mm_setzero_pd();
658 fjy0
= _mm_setzero_pd();
659 fjz0
= _mm_setzero_pd();
661 /**************************
662 * CALCULATE INTERACTIONS *
663 **************************/
665 /* Compute parameters for interactions between i and j atoms */
666 qq00
= _mm_mul_pd(iq0
,jq0
);
667 gmx_mm_load_2pair_swizzle_pd(vdwparam
+vdwioffset0
+vdwjidx0A
,
668 vdwparam
+vdwioffset0
+vdwjidx0B
,&c6_00
,&c12_00
);
670 /* COULOMB ELECTROSTATICS */
671 velec
= _mm_mul_pd(qq00
,rinv00
);
672 felec
= _mm_mul_pd(velec
,rinvsq00
);
674 /* LENNARD-JONES DISPERSION/REPULSION */
676 rinvsix
= _mm_mul_pd(_mm_mul_pd(rinvsq00
,rinvsq00
),rinvsq00
);
677 fvdw
= _mm_mul_pd(_mm_sub_pd(_mm_mul_pd(c12_00
,rinvsix
),c6_00
),_mm_mul_pd(rinvsix
,rinvsq00
));
679 fscal
= _mm_add_pd(felec
,fvdw
);
681 /* Calculate temporary vectorial force */
682 tx
= _mm_mul_pd(fscal
,dx00
);
683 ty
= _mm_mul_pd(fscal
,dy00
);
684 tz
= _mm_mul_pd(fscal
,dz00
);
686 /* Update vectorial force */
687 fix0
= _mm_add_pd(fix0
,tx
);
688 fiy0
= _mm_add_pd(fiy0
,ty
);
689 fiz0
= _mm_add_pd(fiz0
,tz
);
691 fjx0
= _mm_add_pd(fjx0
,tx
);
692 fjy0
= _mm_add_pd(fjy0
,ty
);
693 fjz0
= _mm_add_pd(fjz0
,tz
);
695 /**************************
696 * CALCULATE INTERACTIONS *
697 **************************/
699 /* Compute parameters for interactions between i and j atoms */
700 qq10
= _mm_mul_pd(iq1
,jq0
);
702 /* COULOMB ELECTROSTATICS */
703 velec
= _mm_mul_pd(qq10
,rinv10
);
704 felec
= _mm_mul_pd(velec
,rinvsq10
);
708 /* Calculate temporary vectorial force */
709 tx
= _mm_mul_pd(fscal
,dx10
);
710 ty
= _mm_mul_pd(fscal
,dy10
);
711 tz
= _mm_mul_pd(fscal
,dz10
);
713 /* Update vectorial force */
714 fix1
= _mm_add_pd(fix1
,tx
);
715 fiy1
= _mm_add_pd(fiy1
,ty
);
716 fiz1
= _mm_add_pd(fiz1
,tz
);
718 fjx0
= _mm_add_pd(fjx0
,tx
);
719 fjy0
= _mm_add_pd(fjy0
,ty
);
720 fjz0
= _mm_add_pd(fjz0
,tz
);
722 /**************************
723 * CALCULATE INTERACTIONS *
724 **************************/
726 /* Compute parameters for interactions between i and j atoms */
727 qq20
= _mm_mul_pd(iq2
,jq0
);
729 /* COULOMB ELECTROSTATICS */
730 velec
= _mm_mul_pd(qq20
,rinv20
);
731 felec
= _mm_mul_pd(velec
,rinvsq20
);
735 /* Calculate temporary vectorial force */
736 tx
= _mm_mul_pd(fscal
,dx20
);
737 ty
= _mm_mul_pd(fscal
,dy20
);
738 tz
= _mm_mul_pd(fscal
,dz20
);
740 /* Update vectorial force */
741 fix2
= _mm_add_pd(fix2
,tx
);
742 fiy2
= _mm_add_pd(fiy2
,ty
);
743 fiz2
= _mm_add_pd(fiz2
,tz
);
745 fjx0
= _mm_add_pd(fjx0
,tx
);
746 fjy0
= _mm_add_pd(fjy0
,ty
);
747 fjz0
= _mm_add_pd(fjz0
,tz
);
749 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f
+j_coord_offsetA
,f
+j_coord_offsetB
,fjx0
,fjy0
,fjz0
);
751 /* Inner loop uses 91 flops */
758 j_coord_offsetA
= DIM
*jnrA
;
760 /* load j atom coordinates */
761 gmx_mm_load_1rvec_1ptr_swizzle_pd(x
+j_coord_offsetA
,
764 /* Calculate displacement vector */
765 dx00
= _mm_sub_pd(ix0
,jx0
);
766 dy00
= _mm_sub_pd(iy0
,jy0
);
767 dz00
= _mm_sub_pd(iz0
,jz0
);
768 dx10
= _mm_sub_pd(ix1
,jx0
);
769 dy10
= _mm_sub_pd(iy1
,jy0
);
770 dz10
= _mm_sub_pd(iz1
,jz0
);
771 dx20
= _mm_sub_pd(ix2
,jx0
);
772 dy20
= _mm_sub_pd(iy2
,jy0
);
773 dz20
= _mm_sub_pd(iz2
,jz0
);
775 /* Calculate squared distance and things based on it */
776 rsq00
= gmx_mm_calc_rsq_pd(dx00
,dy00
,dz00
);
777 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
778 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
780 rinv00
= gmx_mm_invsqrt_pd(rsq00
);
781 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
782 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
784 rinvsq00
= _mm_mul_pd(rinv00
,rinv00
);
785 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
786 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
788 /* Load parameters for j particles */
789 jq0
= _mm_load_sd(charge
+jnrA
+0);
790 vdwjidx0A
= 2*vdwtype
[jnrA
+0];
792 fjx0
= _mm_setzero_pd();
793 fjy0
= _mm_setzero_pd();
794 fjz0
= _mm_setzero_pd();
796 /**************************
797 * CALCULATE INTERACTIONS *
798 **************************/
800 /* Compute parameters for interactions between i and j atoms */
801 qq00
= _mm_mul_pd(iq0
,jq0
);
802 gmx_mm_load_1pair_swizzle_pd(vdwparam
+vdwioffset0
+vdwjidx0A
,&c6_00
,&c12_00
);
804 /* COULOMB ELECTROSTATICS */
805 velec
= _mm_mul_pd(qq00
,rinv00
);
806 felec
= _mm_mul_pd(velec
,rinvsq00
);
808 /* LENNARD-JONES DISPERSION/REPULSION */
810 rinvsix
= _mm_mul_pd(_mm_mul_pd(rinvsq00
,rinvsq00
),rinvsq00
);
811 fvdw
= _mm_mul_pd(_mm_sub_pd(_mm_mul_pd(c12_00
,rinvsix
),c6_00
),_mm_mul_pd(rinvsix
,rinvsq00
));
813 fscal
= _mm_add_pd(felec
,fvdw
);
815 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
817 /* Calculate temporary vectorial force */
818 tx
= _mm_mul_pd(fscal
,dx00
);
819 ty
= _mm_mul_pd(fscal
,dy00
);
820 tz
= _mm_mul_pd(fscal
,dz00
);
822 /* Update vectorial force */
823 fix0
= _mm_add_pd(fix0
,tx
);
824 fiy0
= _mm_add_pd(fiy0
,ty
);
825 fiz0
= _mm_add_pd(fiz0
,tz
);
827 fjx0
= _mm_add_pd(fjx0
,tx
);
828 fjy0
= _mm_add_pd(fjy0
,ty
);
829 fjz0
= _mm_add_pd(fjz0
,tz
);
831 /**************************
832 * CALCULATE INTERACTIONS *
833 **************************/
835 /* Compute parameters for interactions between i and j atoms */
836 qq10
= _mm_mul_pd(iq1
,jq0
);
838 /* COULOMB ELECTROSTATICS */
839 velec
= _mm_mul_pd(qq10
,rinv10
);
840 felec
= _mm_mul_pd(velec
,rinvsq10
);
844 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
846 /* Calculate temporary vectorial force */
847 tx
= _mm_mul_pd(fscal
,dx10
);
848 ty
= _mm_mul_pd(fscal
,dy10
);
849 tz
= _mm_mul_pd(fscal
,dz10
);
851 /* Update vectorial force */
852 fix1
= _mm_add_pd(fix1
,tx
);
853 fiy1
= _mm_add_pd(fiy1
,ty
);
854 fiz1
= _mm_add_pd(fiz1
,tz
);
856 fjx0
= _mm_add_pd(fjx0
,tx
);
857 fjy0
= _mm_add_pd(fjy0
,ty
);
858 fjz0
= _mm_add_pd(fjz0
,tz
);
860 /**************************
861 * CALCULATE INTERACTIONS *
862 **************************/
864 /* Compute parameters for interactions between i and j atoms */
865 qq20
= _mm_mul_pd(iq2
,jq0
);
867 /* COULOMB ELECTROSTATICS */
868 velec
= _mm_mul_pd(qq20
,rinv20
);
869 felec
= _mm_mul_pd(velec
,rinvsq20
);
873 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
875 /* Calculate temporary vectorial force */
876 tx
= _mm_mul_pd(fscal
,dx20
);
877 ty
= _mm_mul_pd(fscal
,dy20
);
878 tz
= _mm_mul_pd(fscal
,dz20
);
880 /* Update vectorial force */
881 fix2
= _mm_add_pd(fix2
,tx
);
882 fiy2
= _mm_add_pd(fiy2
,ty
);
883 fiz2
= _mm_add_pd(fiz2
,tz
);
885 fjx0
= _mm_add_pd(fjx0
,tx
);
886 fjy0
= _mm_add_pd(fjy0
,ty
);
887 fjz0
= _mm_add_pd(fjz0
,tz
);
889 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f
+j_coord_offsetA
,fjx0
,fjy0
,fjz0
);
891 /* Inner loop uses 91 flops */
894 /* End of innermost loop */
896 gmx_mm_update_iforce_3atom_swizzle_pd(fix0
,fiy0
,fiz0
,fix1
,fiy1
,fiz1
,fix2
,fiy2
,fiz2
,
897 f
+i_coord_offset
,fshift
+i_shift_offset
);
899 /* Increment number of inner iterations */
900 inneriter
+= j_index_end
- j_index_start
;
902 /* Outer loop uses 18 flops */
905 /* Increment number of outer iterations */
908 /* Update outer/inner flops */
910 inc_nrnb(nrnb
,eNR_NBKERNEL_ELEC_VDW_W3_F
,outeriter
*18 + inneriter
*91);