2 * Note: this file was generated by the Gromacs sse2_double kernel generator.
4 * This source code is part of
8 * Copyright (c) 2001-2012, The GROMACS Development Team
10 * Gromacs is a library for molecular simulation and trajectory analysis,
11 * written by Erik Lindahl, David van der Spoel, Berk Hess, and others - for
12 * a full list of developers and information, check out http://www.gromacs.org
14 * This program is free software; you can redistribute it and/or modify it under
15 * the terms of the GNU Lesser General Public License as published by the Free
16 * Software Foundation; either version 2 of the License, or (at your option) any
19 * To help fund GROMACS development, we humbly ask that you cite
20 * the papers people have written on it - you can find them on the website.
28 #include "../nb_kernel.h"
29 #include "types/simple.h"
33 #include "gmx_math_x86_sse2_double.h"
34 #include "kernelutil_x86_sse2_double.h"
37 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwNone_GeomW4P1_VF_sse2_double
38 * Electrostatics interaction: ReactionField
39 * VdW interaction: None
40 * Geometry: Water4-Particle
41 * Calculate force/pot: PotentialAndForce
44 nb_kernel_ElecRFCut_VdwNone_GeomW4P1_VF_sse2_double
45 (t_nblist
* gmx_restrict nlist
,
46 rvec
* gmx_restrict xx
,
47 rvec
* gmx_restrict ff
,
48 t_forcerec
* gmx_restrict fr
,
49 t_mdatoms
* gmx_restrict mdatoms
,
50 nb_kernel_data_t
* gmx_restrict kernel_data
,
51 t_nrnb
* gmx_restrict nrnb
)
53 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
54 * just 0 for non-waters.
55 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
56 * jnr indices corresponding to data put in the four positions in the SIMD register.
58 int i_shift_offset
,i_coord_offset
,outeriter
,inneriter
;
59 int j_index_start
,j_index_end
,jidx
,nri
,inr
,ggid
,iidx
;
61 int j_coord_offsetA
,j_coord_offsetB
;
62 int *iinr
,*jindex
,*jjnr
,*shiftidx
,*gid
;
64 real
*shiftvec
,*fshift
,*x
,*f
;
65 __m128d tx
,ty
,tz
,fscal
,rcutoff
,rcutoff2
,jidxall
;
67 __m128d ix1
,iy1
,iz1
,fix1
,fiy1
,fiz1
,iq1
,isai1
;
69 __m128d ix2
,iy2
,iz2
,fix2
,fiy2
,fiz2
,iq2
,isai2
;
71 __m128d ix3
,iy3
,iz3
,fix3
,fiy3
,fiz3
,iq3
,isai3
;
72 int vdwjidx0A
,vdwjidx0B
;
73 __m128d jx0
,jy0
,jz0
,fjx0
,fjy0
,fjz0
,jq0
,isaj0
;
74 __m128d dx10
,dy10
,dz10
,rsq10
,rinv10
,rinvsq10
,r10
,qq10
,c6_10
,c12_10
;
75 __m128d dx20
,dy20
,dz20
,rsq20
,rinv20
,rinvsq20
,r20
,qq20
,c6_20
,c12_20
;
76 __m128d dx30
,dy30
,dz30
,rsq30
,rinv30
,rinvsq30
,r30
,qq30
,c6_30
,c12_30
;
77 __m128d velec
,felec
,velecsum
,facel
,crf
,krf
,krf2
;
79 __m128d dummy_mask
,cutoff_mask
;
80 __m128d signbit
= gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
81 __m128d one
= _mm_set1_pd(1.0);
82 __m128d two
= _mm_set1_pd(2.0);
88 jindex
= nlist
->jindex
;
90 shiftidx
= nlist
->shift
;
92 shiftvec
= fr
->shift_vec
[0];
93 fshift
= fr
->fshift
[0];
94 facel
= _mm_set1_pd(fr
->epsfac
);
95 charge
= mdatoms
->chargeA
;
96 krf
= _mm_set1_pd(fr
->ic
->k_rf
);
97 krf2
= _mm_set1_pd(fr
->ic
->k_rf
*2.0);
98 crf
= _mm_set1_pd(fr
->ic
->c_rf
);
100 /* Setup water-specific parameters */
101 inr
= nlist
->iinr
[0];
102 iq1
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+1]));
103 iq2
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+2]));
104 iq3
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+3]));
106 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
107 rcutoff_scalar
= fr
->rcoulomb
;
108 rcutoff
= _mm_set1_pd(rcutoff_scalar
);
109 rcutoff2
= _mm_mul_pd(rcutoff
,rcutoff
);
111 /* Avoid stupid compiler warnings */
119 /* Start outer loop over neighborlists */
120 for(iidx
=0; iidx
<nri
; iidx
++)
122 /* Load shift vector for this list */
123 i_shift_offset
= DIM
*shiftidx
[iidx
];
125 /* Load limits for loop over neighbors */
126 j_index_start
= jindex
[iidx
];
127 j_index_end
= jindex
[iidx
+1];
129 /* Get outer coordinate index */
131 i_coord_offset
= DIM
*inr
;
133 /* Load i particle coords and add shift vector */
134 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec
+i_shift_offset
,x
+i_coord_offset
+DIM
,
135 &ix1
,&iy1
,&iz1
,&ix2
,&iy2
,&iz2
,&ix3
,&iy3
,&iz3
);
137 fix1
= _mm_setzero_pd();
138 fiy1
= _mm_setzero_pd();
139 fiz1
= _mm_setzero_pd();
140 fix2
= _mm_setzero_pd();
141 fiy2
= _mm_setzero_pd();
142 fiz2
= _mm_setzero_pd();
143 fix3
= _mm_setzero_pd();
144 fiy3
= _mm_setzero_pd();
145 fiz3
= _mm_setzero_pd();
147 /* Reset potential sums */
148 velecsum
= _mm_setzero_pd();
150 /* Start inner kernel loop */
151 for(jidx
=j_index_start
; jidx
<j_index_end
-1; jidx
+=2)
154 /* Get j neighbor index, and coordinate index */
157 j_coord_offsetA
= DIM
*jnrA
;
158 j_coord_offsetB
= DIM
*jnrB
;
160 /* load j atom coordinates */
161 gmx_mm_load_1rvec_2ptr_swizzle_pd(x
+j_coord_offsetA
,x
+j_coord_offsetB
,
164 /* Calculate displacement vector */
165 dx10
= _mm_sub_pd(ix1
,jx0
);
166 dy10
= _mm_sub_pd(iy1
,jy0
);
167 dz10
= _mm_sub_pd(iz1
,jz0
);
168 dx20
= _mm_sub_pd(ix2
,jx0
);
169 dy20
= _mm_sub_pd(iy2
,jy0
);
170 dz20
= _mm_sub_pd(iz2
,jz0
);
171 dx30
= _mm_sub_pd(ix3
,jx0
);
172 dy30
= _mm_sub_pd(iy3
,jy0
);
173 dz30
= _mm_sub_pd(iz3
,jz0
);
175 /* Calculate squared distance and things based on it */
176 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
177 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
178 rsq30
= gmx_mm_calc_rsq_pd(dx30
,dy30
,dz30
);
180 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
181 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
182 rinv30
= gmx_mm_invsqrt_pd(rsq30
);
184 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
185 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
186 rinvsq30
= _mm_mul_pd(rinv30
,rinv30
);
188 /* Load parameters for j particles */
189 jq0
= gmx_mm_load_2real_swizzle_pd(charge
+jnrA
+0,charge
+jnrB
+0);
191 fjx0
= _mm_setzero_pd();
192 fjy0
= _mm_setzero_pd();
193 fjz0
= _mm_setzero_pd();
195 /**************************
196 * CALCULATE INTERACTIONS *
197 **************************/
199 if (gmx_mm_any_lt(rsq10
,rcutoff2
))
202 /* Compute parameters for interactions between i and j atoms */
203 qq10
= _mm_mul_pd(iq1
,jq0
);
205 /* REACTION-FIELD ELECTROSTATICS */
206 velec
= _mm_mul_pd(qq10
,_mm_sub_pd(_mm_add_pd(rinv10
,_mm_mul_pd(krf
,rsq10
)),crf
));
207 felec
= _mm_mul_pd(qq10
,_mm_sub_pd(_mm_mul_pd(rinv10
,rinvsq10
),krf2
));
209 cutoff_mask
= _mm_cmplt_pd(rsq10
,rcutoff2
);
211 /* Update potential sum for this i atom from the interaction with this j atom. */
212 velec
= _mm_and_pd(velec
,cutoff_mask
);
213 velecsum
= _mm_add_pd(velecsum
,velec
);
217 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
219 /* Calculate temporary vectorial force */
220 tx
= _mm_mul_pd(fscal
,dx10
);
221 ty
= _mm_mul_pd(fscal
,dy10
);
222 tz
= _mm_mul_pd(fscal
,dz10
);
224 /* Update vectorial force */
225 fix1
= _mm_add_pd(fix1
,tx
);
226 fiy1
= _mm_add_pd(fiy1
,ty
);
227 fiz1
= _mm_add_pd(fiz1
,tz
);
229 fjx0
= _mm_add_pd(fjx0
,tx
);
230 fjy0
= _mm_add_pd(fjy0
,ty
);
231 fjz0
= _mm_add_pd(fjz0
,tz
);
235 /**************************
236 * CALCULATE INTERACTIONS *
237 **************************/
239 if (gmx_mm_any_lt(rsq20
,rcutoff2
))
242 /* Compute parameters for interactions between i and j atoms */
243 qq20
= _mm_mul_pd(iq2
,jq0
);
245 /* REACTION-FIELD ELECTROSTATICS */
246 velec
= _mm_mul_pd(qq20
,_mm_sub_pd(_mm_add_pd(rinv20
,_mm_mul_pd(krf
,rsq20
)),crf
));
247 felec
= _mm_mul_pd(qq20
,_mm_sub_pd(_mm_mul_pd(rinv20
,rinvsq20
),krf2
));
249 cutoff_mask
= _mm_cmplt_pd(rsq20
,rcutoff2
);
251 /* Update potential sum for this i atom from the interaction with this j atom. */
252 velec
= _mm_and_pd(velec
,cutoff_mask
);
253 velecsum
= _mm_add_pd(velecsum
,velec
);
257 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
259 /* Calculate temporary vectorial force */
260 tx
= _mm_mul_pd(fscal
,dx20
);
261 ty
= _mm_mul_pd(fscal
,dy20
);
262 tz
= _mm_mul_pd(fscal
,dz20
);
264 /* Update vectorial force */
265 fix2
= _mm_add_pd(fix2
,tx
);
266 fiy2
= _mm_add_pd(fiy2
,ty
);
267 fiz2
= _mm_add_pd(fiz2
,tz
);
269 fjx0
= _mm_add_pd(fjx0
,tx
);
270 fjy0
= _mm_add_pd(fjy0
,ty
);
271 fjz0
= _mm_add_pd(fjz0
,tz
);
275 /**************************
276 * CALCULATE INTERACTIONS *
277 **************************/
279 if (gmx_mm_any_lt(rsq30
,rcutoff2
))
282 /* Compute parameters for interactions between i and j atoms */
283 qq30
= _mm_mul_pd(iq3
,jq0
);
285 /* REACTION-FIELD ELECTROSTATICS */
286 velec
= _mm_mul_pd(qq30
,_mm_sub_pd(_mm_add_pd(rinv30
,_mm_mul_pd(krf
,rsq30
)),crf
));
287 felec
= _mm_mul_pd(qq30
,_mm_sub_pd(_mm_mul_pd(rinv30
,rinvsq30
),krf2
));
289 cutoff_mask
= _mm_cmplt_pd(rsq30
,rcutoff2
);
291 /* Update potential sum for this i atom from the interaction with this j atom. */
292 velec
= _mm_and_pd(velec
,cutoff_mask
);
293 velecsum
= _mm_add_pd(velecsum
,velec
);
297 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
299 /* Calculate temporary vectorial force */
300 tx
= _mm_mul_pd(fscal
,dx30
);
301 ty
= _mm_mul_pd(fscal
,dy30
);
302 tz
= _mm_mul_pd(fscal
,dz30
);
304 /* Update vectorial force */
305 fix3
= _mm_add_pd(fix3
,tx
);
306 fiy3
= _mm_add_pd(fiy3
,ty
);
307 fiz3
= _mm_add_pd(fiz3
,tz
);
309 fjx0
= _mm_add_pd(fjx0
,tx
);
310 fjy0
= _mm_add_pd(fjy0
,ty
);
311 fjz0
= _mm_add_pd(fjz0
,tz
);
315 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f
+j_coord_offsetA
,f
+j_coord_offsetB
,fjx0
,fjy0
,fjz0
);
317 /* Inner loop uses 111 flops */
324 j_coord_offsetA
= DIM
*jnrA
;
326 /* load j atom coordinates */
327 gmx_mm_load_1rvec_1ptr_swizzle_pd(x
+j_coord_offsetA
,
330 /* Calculate displacement vector */
331 dx10
= _mm_sub_pd(ix1
,jx0
);
332 dy10
= _mm_sub_pd(iy1
,jy0
);
333 dz10
= _mm_sub_pd(iz1
,jz0
);
334 dx20
= _mm_sub_pd(ix2
,jx0
);
335 dy20
= _mm_sub_pd(iy2
,jy0
);
336 dz20
= _mm_sub_pd(iz2
,jz0
);
337 dx30
= _mm_sub_pd(ix3
,jx0
);
338 dy30
= _mm_sub_pd(iy3
,jy0
);
339 dz30
= _mm_sub_pd(iz3
,jz0
);
341 /* Calculate squared distance and things based on it */
342 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
343 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
344 rsq30
= gmx_mm_calc_rsq_pd(dx30
,dy30
,dz30
);
346 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
347 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
348 rinv30
= gmx_mm_invsqrt_pd(rsq30
);
350 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
351 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
352 rinvsq30
= _mm_mul_pd(rinv30
,rinv30
);
354 /* Load parameters for j particles */
355 jq0
= _mm_load_sd(charge
+jnrA
+0);
357 fjx0
= _mm_setzero_pd();
358 fjy0
= _mm_setzero_pd();
359 fjz0
= _mm_setzero_pd();
361 /**************************
362 * CALCULATE INTERACTIONS *
363 **************************/
365 if (gmx_mm_any_lt(rsq10
,rcutoff2
))
368 /* Compute parameters for interactions between i and j atoms */
369 qq10
= _mm_mul_pd(iq1
,jq0
);
371 /* REACTION-FIELD ELECTROSTATICS */
372 velec
= _mm_mul_pd(qq10
,_mm_sub_pd(_mm_add_pd(rinv10
,_mm_mul_pd(krf
,rsq10
)),crf
));
373 felec
= _mm_mul_pd(qq10
,_mm_sub_pd(_mm_mul_pd(rinv10
,rinvsq10
),krf2
));
375 cutoff_mask
= _mm_cmplt_pd(rsq10
,rcutoff2
);
377 /* Update potential sum for this i atom from the interaction with this j atom. */
378 velec
= _mm_and_pd(velec
,cutoff_mask
);
379 velec
= _mm_unpacklo_pd(velec
,_mm_setzero_pd());
380 velecsum
= _mm_add_pd(velecsum
,velec
);
384 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
386 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
388 /* Calculate temporary vectorial force */
389 tx
= _mm_mul_pd(fscal
,dx10
);
390 ty
= _mm_mul_pd(fscal
,dy10
);
391 tz
= _mm_mul_pd(fscal
,dz10
);
393 /* Update vectorial force */
394 fix1
= _mm_add_pd(fix1
,tx
);
395 fiy1
= _mm_add_pd(fiy1
,ty
);
396 fiz1
= _mm_add_pd(fiz1
,tz
);
398 fjx0
= _mm_add_pd(fjx0
,tx
);
399 fjy0
= _mm_add_pd(fjy0
,ty
);
400 fjz0
= _mm_add_pd(fjz0
,tz
);
404 /**************************
405 * CALCULATE INTERACTIONS *
406 **************************/
408 if (gmx_mm_any_lt(rsq20
,rcutoff2
))
411 /* Compute parameters for interactions between i and j atoms */
412 qq20
= _mm_mul_pd(iq2
,jq0
);
414 /* REACTION-FIELD ELECTROSTATICS */
415 velec
= _mm_mul_pd(qq20
,_mm_sub_pd(_mm_add_pd(rinv20
,_mm_mul_pd(krf
,rsq20
)),crf
));
416 felec
= _mm_mul_pd(qq20
,_mm_sub_pd(_mm_mul_pd(rinv20
,rinvsq20
),krf2
));
418 cutoff_mask
= _mm_cmplt_pd(rsq20
,rcutoff2
);
420 /* Update potential sum for this i atom from the interaction with this j atom. */
421 velec
= _mm_and_pd(velec
,cutoff_mask
);
422 velec
= _mm_unpacklo_pd(velec
,_mm_setzero_pd());
423 velecsum
= _mm_add_pd(velecsum
,velec
);
427 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
429 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
431 /* Calculate temporary vectorial force */
432 tx
= _mm_mul_pd(fscal
,dx20
);
433 ty
= _mm_mul_pd(fscal
,dy20
);
434 tz
= _mm_mul_pd(fscal
,dz20
);
436 /* Update vectorial force */
437 fix2
= _mm_add_pd(fix2
,tx
);
438 fiy2
= _mm_add_pd(fiy2
,ty
);
439 fiz2
= _mm_add_pd(fiz2
,tz
);
441 fjx0
= _mm_add_pd(fjx0
,tx
);
442 fjy0
= _mm_add_pd(fjy0
,ty
);
443 fjz0
= _mm_add_pd(fjz0
,tz
);
447 /**************************
448 * CALCULATE INTERACTIONS *
449 **************************/
451 if (gmx_mm_any_lt(rsq30
,rcutoff2
))
454 /* Compute parameters for interactions between i and j atoms */
455 qq30
= _mm_mul_pd(iq3
,jq0
);
457 /* REACTION-FIELD ELECTROSTATICS */
458 velec
= _mm_mul_pd(qq30
,_mm_sub_pd(_mm_add_pd(rinv30
,_mm_mul_pd(krf
,rsq30
)),crf
));
459 felec
= _mm_mul_pd(qq30
,_mm_sub_pd(_mm_mul_pd(rinv30
,rinvsq30
),krf2
));
461 cutoff_mask
= _mm_cmplt_pd(rsq30
,rcutoff2
);
463 /* Update potential sum for this i atom from the interaction with this j atom. */
464 velec
= _mm_and_pd(velec
,cutoff_mask
);
465 velec
= _mm_unpacklo_pd(velec
,_mm_setzero_pd());
466 velecsum
= _mm_add_pd(velecsum
,velec
);
470 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
472 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
474 /* Calculate temporary vectorial force */
475 tx
= _mm_mul_pd(fscal
,dx30
);
476 ty
= _mm_mul_pd(fscal
,dy30
);
477 tz
= _mm_mul_pd(fscal
,dz30
);
479 /* Update vectorial force */
480 fix3
= _mm_add_pd(fix3
,tx
);
481 fiy3
= _mm_add_pd(fiy3
,ty
);
482 fiz3
= _mm_add_pd(fiz3
,tz
);
484 fjx0
= _mm_add_pd(fjx0
,tx
);
485 fjy0
= _mm_add_pd(fjy0
,ty
);
486 fjz0
= _mm_add_pd(fjz0
,tz
);
490 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f
+j_coord_offsetA
,fjx0
,fjy0
,fjz0
);
492 /* Inner loop uses 111 flops */
495 /* End of innermost loop */
497 gmx_mm_update_iforce_3atom_swizzle_pd(fix1
,fiy1
,fiz1
,fix2
,fiy2
,fiz2
,fix3
,fiy3
,fiz3
,
498 f
+i_coord_offset
+DIM
,fshift
+i_shift_offset
);
501 /* Update potential energies */
502 gmx_mm_update_1pot_pd(velecsum
,kernel_data
->energygrp_elec
+ggid
);
504 /* Increment number of inner iterations */
505 inneriter
+= j_index_end
- j_index_start
;
507 /* Outer loop uses 19 flops */
510 /* Increment number of outer iterations */
513 /* Update outer/inner flops */
515 inc_nrnb(nrnb
,eNR_NBKERNEL_ELEC_W4_VF
,outeriter
*19 + inneriter
*111);
518 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwNone_GeomW4P1_F_sse2_double
519 * Electrostatics interaction: ReactionField
520 * VdW interaction: None
521 * Geometry: Water4-Particle
522 * Calculate force/pot: Force
525 nb_kernel_ElecRFCut_VdwNone_GeomW4P1_F_sse2_double
526 (t_nblist
* gmx_restrict nlist
,
527 rvec
* gmx_restrict xx
,
528 rvec
* gmx_restrict ff
,
529 t_forcerec
* gmx_restrict fr
,
530 t_mdatoms
* gmx_restrict mdatoms
,
531 nb_kernel_data_t
* gmx_restrict kernel_data
,
532 t_nrnb
* gmx_restrict nrnb
)
534 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
535 * just 0 for non-waters.
536 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
537 * jnr indices corresponding to data put in the four positions in the SIMD register.
539 int i_shift_offset
,i_coord_offset
,outeriter
,inneriter
;
540 int j_index_start
,j_index_end
,jidx
,nri
,inr
,ggid
,iidx
;
542 int j_coord_offsetA
,j_coord_offsetB
;
543 int *iinr
,*jindex
,*jjnr
,*shiftidx
,*gid
;
545 real
*shiftvec
,*fshift
,*x
,*f
;
546 __m128d tx
,ty
,tz
,fscal
,rcutoff
,rcutoff2
,jidxall
;
548 __m128d ix1
,iy1
,iz1
,fix1
,fiy1
,fiz1
,iq1
,isai1
;
550 __m128d ix2
,iy2
,iz2
,fix2
,fiy2
,fiz2
,iq2
,isai2
;
552 __m128d ix3
,iy3
,iz3
,fix3
,fiy3
,fiz3
,iq3
,isai3
;
553 int vdwjidx0A
,vdwjidx0B
;
554 __m128d jx0
,jy0
,jz0
,fjx0
,fjy0
,fjz0
,jq0
,isaj0
;
555 __m128d dx10
,dy10
,dz10
,rsq10
,rinv10
,rinvsq10
,r10
,qq10
,c6_10
,c12_10
;
556 __m128d dx20
,dy20
,dz20
,rsq20
,rinv20
,rinvsq20
,r20
,qq20
,c6_20
,c12_20
;
557 __m128d dx30
,dy30
,dz30
,rsq30
,rinv30
,rinvsq30
,r30
,qq30
,c6_30
,c12_30
;
558 __m128d velec
,felec
,velecsum
,facel
,crf
,krf
,krf2
;
560 __m128d dummy_mask
,cutoff_mask
;
561 __m128d signbit
= gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
562 __m128d one
= _mm_set1_pd(1.0);
563 __m128d two
= _mm_set1_pd(2.0);
569 jindex
= nlist
->jindex
;
571 shiftidx
= nlist
->shift
;
573 shiftvec
= fr
->shift_vec
[0];
574 fshift
= fr
->fshift
[0];
575 facel
= _mm_set1_pd(fr
->epsfac
);
576 charge
= mdatoms
->chargeA
;
577 krf
= _mm_set1_pd(fr
->ic
->k_rf
);
578 krf2
= _mm_set1_pd(fr
->ic
->k_rf
*2.0);
579 crf
= _mm_set1_pd(fr
->ic
->c_rf
);
581 /* Setup water-specific parameters */
582 inr
= nlist
->iinr
[0];
583 iq1
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+1]));
584 iq2
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+2]));
585 iq3
= _mm_mul_pd(facel
,_mm_set1_pd(charge
[inr
+3]));
587 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
588 rcutoff_scalar
= fr
->rcoulomb
;
589 rcutoff
= _mm_set1_pd(rcutoff_scalar
);
590 rcutoff2
= _mm_mul_pd(rcutoff
,rcutoff
);
592 /* Avoid stupid compiler warnings */
600 /* Start outer loop over neighborlists */
601 for(iidx
=0; iidx
<nri
; iidx
++)
603 /* Load shift vector for this list */
604 i_shift_offset
= DIM
*shiftidx
[iidx
];
606 /* Load limits for loop over neighbors */
607 j_index_start
= jindex
[iidx
];
608 j_index_end
= jindex
[iidx
+1];
610 /* Get outer coordinate index */
612 i_coord_offset
= DIM
*inr
;
614 /* Load i particle coords and add shift vector */
615 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec
+i_shift_offset
,x
+i_coord_offset
+DIM
,
616 &ix1
,&iy1
,&iz1
,&ix2
,&iy2
,&iz2
,&ix3
,&iy3
,&iz3
);
618 fix1
= _mm_setzero_pd();
619 fiy1
= _mm_setzero_pd();
620 fiz1
= _mm_setzero_pd();
621 fix2
= _mm_setzero_pd();
622 fiy2
= _mm_setzero_pd();
623 fiz2
= _mm_setzero_pd();
624 fix3
= _mm_setzero_pd();
625 fiy3
= _mm_setzero_pd();
626 fiz3
= _mm_setzero_pd();
628 /* Start inner kernel loop */
629 for(jidx
=j_index_start
; jidx
<j_index_end
-1; jidx
+=2)
632 /* Get j neighbor index, and coordinate index */
635 j_coord_offsetA
= DIM
*jnrA
;
636 j_coord_offsetB
= DIM
*jnrB
;
638 /* load j atom coordinates */
639 gmx_mm_load_1rvec_2ptr_swizzle_pd(x
+j_coord_offsetA
,x
+j_coord_offsetB
,
642 /* Calculate displacement vector */
643 dx10
= _mm_sub_pd(ix1
,jx0
);
644 dy10
= _mm_sub_pd(iy1
,jy0
);
645 dz10
= _mm_sub_pd(iz1
,jz0
);
646 dx20
= _mm_sub_pd(ix2
,jx0
);
647 dy20
= _mm_sub_pd(iy2
,jy0
);
648 dz20
= _mm_sub_pd(iz2
,jz0
);
649 dx30
= _mm_sub_pd(ix3
,jx0
);
650 dy30
= _mm_sub_pd(iy3
,jy0
);
651 dz30
= _mm_sub_pd(iz3
,jz0
);
653 /* Calculate squared distance and things based on it */
654 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
655 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
656 rsq30
= gmx_mm_calc_rsq_pd(dx30
,dy30
,dz30
);
658 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
659 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
660 rinv30
= gmx_mm_invsqrt_pd(rsq30
);
662 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
663 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
664 rinvsq30
= _mm_mul_pd(rinv30
,rinv30
);
666 /* Load parameters for j particles */
667 jq0
= gmx_mm_load_2real_swizzle_pd(charge
+jnrA
+0,charge
+jnrB
+0);
669 fjx0
= _mm_setzero_pd();
670 fjy0
= _mm_setzero_pd();
671 fjz0
= _mm_setzero_pd();
673 /**************************
674 * CALCULATE INTERACTIONS *
675 **************************/
677 if (gmx_mm_any_lt(rsq10
,rcutoff2
))
680 /* Compute parameters for interactions between i and j atoms */
681 qq10
= _mm_mul_pd(iq1
,jq0
);
683 /* REACTION-FIELD ELECTROSTATICS */
684 felec
= _mm_mul_pd(qq10
,_mm_sub_pd(_mm_mul_pd(rinv10
,rinvsq10
),krf2
));
686 cutoff_mask
= _mm_cmplt_pd(rsq10
,rcutoff2
);
690 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
692 /* Calculate temporary vectorial force */
693 tx
= _mm_mul_pd(fscal
,dx10
);
694 ty
= _mm_mul_pd(fscal
,dy10
);
695 tz
= _mm_mul_pd(fscal
,dz10
);
697 /* Update vectorial force */
698 fix1
= _mm_add_pd(fix1
,tx
);
699 fiy1
= _mm_add_pd(fiy1
,ty
);
700 fiz1
= _mm_add_pd(fiz1
,tz
);
702 fjx0
= _mm_add_pd(fjx0
,tx
);
703 fjy0
= _mm_add_pd(fjy0
,ty
);
704 fjz0
= _mm_add_pd(fjz0
,tz
);
708 /**************************
709 * CALCULATE INTERACTIONS *
710 **************************/
712 if (gmx_mm_any_lt(rsq20
,rcutoff2
))
715 /* Compute parameters for interactions between i and j atoms */
716 qq20
= _mm_mul_pd(iq2
,jq0
);
718 /* REACTION-FIELD ELECTROSTATICS */
719 felec
= _mm_mul_pd(qq20
,_mm_sub_pd(_mm_mul_pd(rinv20
,rinvsq20
),krf2
));
721 cutoff_mask
= _mm_cmplt_pd(rsq20
,rcutoff2
);
725 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
727 /* Calculate temporary vectorial force */
728 tx
= _mm_mul_pd(fscal
,dx20
);
729 ty
= _mm_mul_pd(fscal
,dy20
);
730 tz
= _mm_mul_pd(fscal
,dz20
);
732 /* Update vectorial force */
733 fix2
= _mm_add_pd(fix2
,tx
);
734 fiy2
= _mm_add_pd(fiy2
,ty
);
735 fiz2
= _mm_add_pd(fiz2
,tz
);
737 fjx0
= _mm_add_pd(fjx0
,tx
);
738 fjy0
= _mm_add_pd(fjy0
,ty
);
739 fjz0
= _mm_add_pd(fjz0
,tz
);
743 /**************************
744 * CALCULATE INTERACTIONS *
745 **************************/
747 if (gmx_mm_any_lt(rsq30
,rcutoff2
))
750 /* Compute parameters for interactions between i and j atoms */
751 qq30
= _mm_mul_pd(iq3
,jq0
);
753 /* REACTION-FIELD ELECTROSTATICS */
754 felec
= _mm_mul_pd(qq30
,_mm_sub_pd(_mm_mul_pd(rinv30
,rinvsq30
),krf2
));
756 cutoff_mask
= _mm_cmplt_pd(rsq30
,rcutoff2
);
760 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
762 /* Calculate temporary vectorial force */
763 tx
= _mm_mul_pd(fscal
,dx30
);
764 ty
= _mm_mul_pd(fscal
,dy30
);
765 tz
= _mm_mul_pd(fscal
,dz30
);
767 /* Update vectorial force */
768 fix3
= _mm_add_pd(fix3
,tx
);
769 fiy3
= _mm_add_pd(fiy3
,ty
);
770 fiz3
= _mm_add_pd(fiz3
,tz
);
772 fjx0
= _mm_add_pd(fjx0
,tx
);
773 fjy0
= _mm_add_pd(fjy0
,ty
);
774 fjz0
= _mm_add_pd(fjz0
,tz
);
778 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f
+j_coord_offsetA
,f
+j_coord_offsetB
,fjx0
,fjy0
,fjz0
);
780 /* Inner loop uses 93 flops */
787 j_coord_offsetA
= DIM
*jnrA
;
789 /* load j atom coordinates */
790 gmx_mm_load_1rvec_1ptr_swizzle_pd(x
+j_coord_offsetA
,
793 /* Calculate displacement vector */
794 dx10
= _mm_sub_pd(ix1
,jx0
);
795 dy10
= _mm_sub_pd(iy1
,jy0
);
796 dz10
= _mm_sub_pd(iz1
,jz0
);
797 dx20
= _mm_sub_pd(ix2
,jx0
);
798 dy20
= _mm_sub_pd(iy2
,jy0
);
799 dz20
= _mm_sub_pd(iz2
,jz0
);
800 dx30
= _mm_sub_pd(ix3
,jx0
);
801 dy30
= _mm_sub_pd(iy3
,jy0
);
802 dz30
= _mm_sub_pd(iz3
,jz0
);
804 /* Calculate squared distance and things based on it */
805 rsq10
= gmx_mm_calc_rsq_pd(dx10
,dy10
,dz10
);
806 rsq20
= gmx_mm_calc_rsq_pd(dx20
,dy20
,dz20
);
807 rsq30
= gmx_mm_calc_rsq_pd(dx30
,dy30
,dz30
);
809 rinv10
= gmx_mm_invsqrt_pd(rsq10
);
810 rinv20
= gmx_mm_invsqrt_pd(rsq20
);
811 rinv30
= gmx_mm_invsqrt_pd(rsq30
);
813 rinvsq10
= _mm_mul_pd(rinv10
,rinv10
);
814 rinvsq20
= _mm_mul_pd(rinv20
,rinv20
);
815 rinvsq30
= _mm_mul_pd(rinv30
,rinv30
);
817 /* Load parameters for j particles */
818 jq0
= _mm_load_sd(charge
+jnrA
+0);
820 fjx0
= _mm_setzero_pd();
821 fjy0
= _mm_setzero_pd();
822 fjz0
= _mm_setzero_pd();
824 /**************************
825 * CALCULATE INTERACTIONS *
826 **************************/
828 if (gmx_mm_any_lt(rsq10
,rcutoff2
))
831 /* Compute parameters for interactions between i and j atoms */
832 qq10
= _mm_mul_pd(iq1
,jq0
);
834 /* REACTION-FIELD ELECTROSTATICS */
835 felec
= _mm_mul_pd(qq10
,_mm_sub_pd(_mm_mul_pd(rinv10
,rinvsq10
),krf2
));
837 cutoff_mask
= _mm_cmplt_pd(rsq10
,rcutoff2
);
841 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
843 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
845 /* Calculate temporary vectorial force */
846 tx
= _mm_mul_pd(fscal
,dx10
);
847 ty
= _mm_mul_pd(fscal
,dy10
);
848 tz
= _mm_mul_pd(fscal
,dz10
);
850 /* Update vectorial force */
851 fix1
= _mm_add_pd(fix1
,tx
);
852 fiy1
= _mm_add_pd(fiy1
,ty
);
853 fiz1
= _mm_add_pd(fiz1
,tz
);
855 fjx0
= _mm_add_pd(fjx0
,tx
);
856 fjy0
= _mm_add_pd(fjy0
,ty
);
857 fjz0
= _mm_add_pd(fjz0
,tz
);
861 /**************************
862 * CALCULATE INTERACTIONS *
863 **************************/
865 if (gmx_mm_any_lt(rsq20
,rcutoff2
))
868 /* Compute parameters for interactions between i and j atoms */
869 qq20
= _mm_mul_pd(iq2
,jq0
);
871 /* REACTION-FIELD ELECTROSTATICS */
872 felec
= _mm_mul_pd(qq20
,_mm_sub_pd(_mm_mul_pd(rinv20
,rinvsq20
),krf2
));
874 cutoff_mask
= _mm_cmplt_pd(rsq20
,rcutoff2
);
878 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
880 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
882 /* Calculate temporary vectorial force */
883 tx
= _mm_mul_pd(fscal
,dx20
);
884 ty
= _mm_mul_pd(fscal
,dy20
);
885 tz
= _mm_mul_pd(fscal
,dz20
);
887 /* Update vectorial force */
888 fix2
= _mm_add_pd(fix2
,tx
);
889 fiy2
= _mm_add_pd(fiy2
,ty
);
890 fiz2
= _mm_add_pd(fiz2
,tz
);
892 fjx0
= _mm_add_pd(fjx0
,tx
);
893 fjy0
= _mm_add_pd(fjy0
,ty
);
894 fjz0
= _mm_add_pd(fjz0
,tz
);
898 /**************************
899 * CALCULATE INTERACTIONS *
900 **************************/
902 if (gmx_mm_any_lt(rsq30
,rcutoff2
))
905 /* Compute parameters for interactions between i and j atoms */
906 qq30
= _mm_mul_pd(iq3
,jq0
);
908 /* REACTION-FIELD ELECTROSTATICS */
909 felec
= _mm_mul_pd(qq30
,_mm_sub_pd(_mm_mul_pd(rinv30
,rinvsq30
),krf2
));
911 cutoff_mask
= _mm_cmplt_pd(rsq30
,rcutoff2
);
915 fscal
= _mm_and_pd(fscal
,cutoff_mask
);
917 fscal
= _mm_unpacklo_pd(fscal
,_mm_setzero_pd());
919 /* Calculate temporary vectorial force */
920 tx
= _mm_mul_pd(fscal
,dx30
);
921 ty
= _mm_mul_pd(fscal
,dy30
);
922 tz
= _mm_mul_pd(fscal
,dz30
);
924 /* Update vectorial force */
925 fix3
= _mm_add_pd(fix3
,tx
);
926 fiy3
= _mm_add_pd(fiy3
,ty
);
927 fiz3
= _mm_add_pd(fiz3
,tz
);
929 fjx0
= _mm_add_pd(fjx0
,tx
);
930 fjy0
= _mm_add_pd(fjy0
,ty
);
931 fjz0
= _mm_add_pd(fjz0
,tz
);
935 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f
+j_coord_offsetA
,fjx0
,fjy0
,fjz0
);
937 /* Inner loop uses 93 flops */
940 /* End of innermost loop */
942 gmx_mm_update_iforce_3atom_swizzle_pd(fix1
,fiy1
,fiz1
,fix2
,fiy2
,fiz2
,fix3
,fiy3
,fiz3
,
943 f
+i_coord_offset
+DIM
,fshift
+i_shift_offset
);
945 /* Increment number of inner iterations */
946 inneriter
+= j_index_end
- j_index_start
;
948 /* Outer loop uses 18 flops */
951 /* Increment number of outer iterations */
954 /* Update outer/inner flops */
956 inc_nrnb(nrnb
,eNR_NBKERNEL_ELEC_W4_F
,outeriter
*18 + inneriter
*93);