?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Removed simple.h from nb_kernel_sse4_1_XX
[gromacs.git] / src / gromacs / gmxlib / nonbonded / nb_kernel_sse4_1_double / nb_kernel_ElecRFCut_VdwNone_GeomW4P1_sse4_1_double.c
blob9b59b049dbf61b332657429b51a1e298aa7591d8
1 /*
2 * This file is part of the GROMACS molecular simulation package.
3 *
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.
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 sse4_1_double kernel generator.
37 */
38 #include "gmxpre.h"
39
40 #include "config.h"
41
42 #include <math.h>
43
44 #include "../nb_kernel.h"
45 #include "gromacs/math/vec.h"
46 #include "gromacs/legacyheaders/nrnb.h"
47
48 #include "gromacs/simd/math_x86_sse4_1_double.h"
49 #include "kernelutil_x86_sse4_1_double.h"
50
51 /*
52 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwNone_GeomW4P1_VF_sse4_1_double
53 * Electrostatics interaction: ReactionField
54 * VdW interaction: None
55 * Geometry: Water4-Particle
56 * Calculate force/pot: PotentialAndForce
57 */
58 void
59 nb_kernel_ElecRFCut_VdwNone_GeomW4P1_VF_sse4_1_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 SSE double precision, 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 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
81 int vdwioffset1;
82 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
83 int vdwioffset2;
84 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
85 int vdwioffset3;
86 __m128d ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
87 int vdwjidx0A,vdwjidx0B;
88 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
89 __m128d dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
90 __m128d dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
91 __m128d dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
92 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
93 real *charge;
94 __m128d dummy_mask,cutoff_mask;
95 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
96 __m128d one = _mm_set1_pd(1.0);
97 __m128d two = _mm_set1_pd(2.0);
98 x = xx[0];
99 f = ff[0];
100
101 nri = nlist->nri;
102 iinr = nlist->iinr;
103 jindex = nlist->jindex;
104 jjnr = nlist->jjnr;
105 shiftidx = nlist->shift;
106 gid = nlist->gid;
107 shiftvec = fr->shift_vec[0];
108 fshift = fr->fshift[0];
109 facel = _mm_set1_pd(fr->epsfac);
110 charge = mdatoms->chargeA;
111 krf = _mm_set1_pd(fr->ic->k_rf);
112 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
113 crf = _mm_set1_pd(fr->ic->c_rf);
114
115 /* Setup water-specific parameters */
116 inr = nlist->iinr[0];
117 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
118 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
119 iq3 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+3]));
120
121 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
122 rcutoff_scalar = fr->rcoulomb;
123 rcutoff = _mm_set1_pd(rcutoff_scalar);
124 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
125
126 /* Avoid stupid compiler warnings */
127 jnrA = jnrB = 0;
128 j_coord_offsetA = 0;
129 j_coord_offsetB = 0;
130
131 outeriter = 0;
132 inneriter = 0;
133
134 /* Start outer loop over neighborlists */
135 for(iidx=0; iidx<nri; iidx++)
136 {
137 /* Load shift vector for this list */
138 i_shift_offset = DIM*shiftidx[iidx];
139
140 /* Load limits for loop over neighbors */
141 j_index_start = jindex[iidx];
142 j_index_end = jindex[iidx+1];
143
144 /* Get outer coordinate index */
145 inr = iinr[iidx];
146 i_coord_offset = DIM*inr;
147
148 /* Load i particle coords and add shift vector */
149 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset+DIM,
150 &ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
151
152 fix1 = _mm_setzero_pd();
153 fiy1 = _mm_setzero_pd();
154 fiz1 = _mm_setzero_pd();
155 fix2 = _mm_setzero_pd();
156 fiy2 = _mm_setzero_pd();
157 fiz2 = _mm_setzero_pd();
158 fix3 = _mm_setzero_pd();
159 fiy3 = _mm_setzero_pd();
160 fiz3 = _mm_setzero_pd();
161
162 /* Reset potential sums */
163 velecsum = _mm_setzero_pd();
164
165 /* Start inner kernel loop */
166 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
167 {
168
169 /* Get j neighbor index, and coordinate index */
170 jnrA = jjnr[jidx];
171 jnrB = jjnr[jidx+1];
172 j_coord_offsetA = DIM*jnrA;
173 j_coord_offsetB = DIM*jnrB;
174
175 /* load j atom coordinates */
176 gmx_mm_load_1rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
177 &jx0,&jy0,&jz0);
178
179 /* Calculate displacement vector */
180 dx10 = _mm_sub_pd(ix1,jx0);
181 dy10 = _mm_sub_pd(iy1,jy0);
182 dz10 = _mm_sub_pd(iz1,jz0);
183 dx20 = _mm_sub_pd(ix2,jx0);
184 dy20 = _mm_sub_pd(iy2,jy0);
185 dz20 = _mm_sub_pd(iz2,jz0);
186 dx30 = _mm_sub_pd(ix3,jx0);
187 dy30 = _mm_sub_pd(iy3,jy0);
188 dz30 = _mm_sub_pd(iz3,jz0);
189
190 /* Calculate squared distance and things based on it */
191 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
192 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
193 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
194
195 rinv10 = gmx_mm_invsqrt_pd(rsq10);
196 rinv20 = gmx_mm_invsqrt_pd(rsq20);
197 rinv30 = gmx_mm_invsqrt_pd(rsq30);
198
199 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
200 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
201 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
202
203 /* Load parameters for j particles */
204 jq0 = gmx_mm_load_2real_swizzle_pd(charge+jnrA+0,charge+jnrB+0);
205
206 fjx0 = _mm_setzero_pd();
207 fjy0 = _mm_setzero_pd();
208 fjz0 = _mm_setzero_pd();
209
210 /**************************
211 * CALCULATE INTERACTIONS *
212 **************************/
213
214 if (gmx_mm_any_lt(rsq10,rcutoff2))
215 {
216
217 /* Compute parameters for interactions between i and j atoms */
218 qq10 = _mm_mul_pd(iq1,jq0);
219
220 /* REACTION-FIELD ELECTROSTATICS */
221 velec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_add_pd(rinv10,_mm_mul_pd(krf,rsq10)),crf));
222 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
223
224 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
225
226 /* Update potential sum for this i atom from the interaction with this j atom. */
227 velec = _mm_and_pd(velec,cutoff_mask);
228 velecsum = _mm_add_pd(velecsum,velec);
229
230 fscal = felec;
231
232 fscal = _mm_and_pd(fscal,cutoff_mask);
233
234 /* Calculate temporary vectorial force */
235 tx = _mm_mul_pd(fscal,dx10);
236 ty = _mm_mul_pd(fscal,dy10);
237 tz = _mm_mul_pd(fscal,dz10);
238
239 /* Update vectorial force */
240 fix1 = _mm_add_pd(fix1,tx);
241 fiy1 = _mm_add_pd(fiy1,ty);
242 fiz1 = _mm_add_pd(fiz1,tz);
243
244 fjx0 = _mm_add_pd(fjx0,tx);
245 fjy0 = _mm_add_pd(fjy0,ty);
246 fjz0 = _mm_add_pd(fjz0,tz);
247
248 }
249
250 /**************************
251 * CALCULATE INTERACTIONS *
252 **************************/
253
254 if (gmx_mm_any_lt(rsq20,rcutoff2))
255 {
256
257 /* Compute parameters for interactions between i and j atoms */
258 qq20 = _mm_mul_pd(iq2,jq0);
259
260 /* REACTION-FIELD ELECTROSTATICS */
261 velec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_add_pd(rinv20,_mm_mul_pd(krf,rsq20)),crf));
262 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
263
264 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
265
266 /* Update potential sum for this i atom from the interaction with this j atom. */
267 velec = _mm_and_pd(velec,cutoff_mask);
268 velecsum = _mm_add_pd(velecsum,velec);
269
270 fscal = felec;
271
272 fscal = _mm_and_pd(fscal,cutoff_mask);
273
274 /* Calculate temporary vectorial force */
275 tx = _mm_mul_pd(fscal,dx20);
276 ty = _mm_mul_pd(fscal,dy20);
277 tz = _mm_mul_pd(fscal,dz20);
278
279 /* Update vectorial force */
280 fix2 = _mm_add_pd(fix2,tx);
281 fiy2 = _mm_add_pd(fiy2,ty);
282 fiz2 = _mm_add_pd(fiz2,tz);
283
284 fjx0 = _mm_add_pd(fjx0,tx);
285 fjy0 = _mm_add_pd(fjy0,ty);
286 fjz0 = _mm_add_pd(fjz0,tz);
287
288 }
289
290 /**************************
291 * CALCULATE INTERACTIONS *
292 **************************/
293
294 if (gmx_mm_any_lt(rsq30,rcutoff2))
295 {
296
297 /* Compute parameters for interactions between i and j atoms */
298 qq30 = _mm_mul_pd(iq3,jq0);
299
300 /* REACTION-FIELD ELECTROSTATICS */
301 velec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_add_pd(rinv30,_mm_mul_pd(krf,rsq30)),crf));
302 felec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_mul_pd(rinv30,rinvsq30),krf2));
303
304 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
305
306 /* Update potential sum for this i atom from the interaction with this j atom. */
307 velec = _mm_and_pd(velec,cutoff_mask);
308 velecsum = _mm_add_pd(velecsum,velec);
309
310 fscal = felec;
311
312 fscal = _mm_and_pd(fscal,cutoff_mask);
313
314 /* Calculate temporary vectorial force */
315 tx = _mm_mul_pd(fscal,dx30);
316 ty = _mm_mul_pd(fscal,dy30);
317 tz = _mm_mul_pd(fscal,dz30);
318
319 /* Update vectorial force */
320 fix3 = _mm_add_pd(fix3,tx);
321 fiy3 = _mm_add_pd(fiy3,ty);
322 fiz3 = _mm_add_pd(fiz3,tz);
323
324 fjx0 = _mm_add_pd(fjx0,tx);
325 fjy0 = _mm_add_pd(fjy0,ty);
326 fjz0 = _mm_add_pd(fjz0,tz);
327
328 }
329
330 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0);
331
332 /* Inner loop uses 111 flops */
333 }
334
335 if(jidx<j_index_end)
336 {
337
338 jnrA = jjnr[jidx];
339 j_coord_offsetA = DIM*jnrA;
340
341 /* load j atom coordinates */
342 gmx_mm_load_1rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
343 &jx0,&jy0,&jz0);
344
345 /* Calculate displacement vector */
346 dx10 = _mm_sub_pd(ix1,jx0);
347 dy10 = _mm_sub_pd(iy1,jy0);
348 dz10 = _mm_sub_pd(iz1,jz0);
349 dx20 = _mm_sub_pd(ix2,jx0);
350 dy20 = _mm_sub_pd(iy2,jy0);
351 dz20 = _mm_sub_pd(iz2,jz0);
352 dx30 = _mm_sub_pd(ix3,jx0);
353 dy30 = _mm_sub_pd(iy3,jy0);
354 dz30 = _mm_sub_pd(iz3,jz0);
355
356 /* Calculate squared distance and things based on it */
357 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
358 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
359 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
360
361 rinv10 = gmx_mm_invsqrt_pd(rsq10);
362 rinv20 = gmx_mm_invsqrt_pd(rsq20);
363 rinv30 = gmx_mm_invsqrt_pd(rsq30);
364
365 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
366 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
367 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
368
369 /* Load parameters for j particles */
370 jq0 = _mm_load_sd(charge+jnrA+0);
371
372 fjx0 = _mm_setzero_pd();
373 fjy0 = _mm_setzero_pd();
374 fjz0 = _mm_setzero_pd();
375
376 /**************************
377 * CALCULATE INTERACTIONS *
378 **************************/
379
380 if (gmx_mm_any_lt(rsq10,rcutoff2))
381 {
382
383 /* Compute parameters for interactions between i and j atoms */
384 qq10 = _mm_mul_pd(iq1,jq0);
385
386 /* REACTION-FIELD ELECTROSTATICS */
387 velec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_add_pd(rinv10,_mm_mul_pd(krf,rsq10)),crf));
388 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
389
390 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
391
392 /* Update potential sum for this i atom from the interaction with this j atom. */
393 velec = _mm_and_pd(velec,cutoff_mask);
394 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
395 velecsum = _mm_add_pd(velecsum,velec);
396
397 fscal = felec;
398
399 fscal = _mm_and_pd(fscal,cutoff_mask);
400
401 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
402
403 /* Calculate temporary vectorial force */
404 tx = _mm_mul_pd(fscal,dx10);
405 ty = _mm_mul_pd(fscal,dy10);
406 tz = _mm_mul_pd(fscal,dz10);
407
408 /* Update vectorial force */
409 fix1 = _mm_add_pd(fix1,tx);
410 fiy1 = _mm_add_pd(fiy1,ty);
411 fiz1 = _mm_add_pd(fiz1,tz);
412
413 fjx0 = _mm_add_pd(fjx0,tx);
414 fjy0 = _mm_add_pd(fjy0,ty);
415 fjz0 = _mm_add_pd(fjz0,tz);
416
417 }
418
419 /**************************
420 * CALCULATE INTERACTIONS *
421 **************************/
422
423 if (gmx_mm_any_lt(rsq20,rcutoff2))
424 {
425
426 /* Compute parameters for interactions between i and j atoms */
427 qq20 = _mm_mul_pd(iq2,jq0);
428
429 /* REACTION-FIELD ELECTROSTATICS */
430 velec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_add_pd(rinv20,_mm_mul_pd(krf,rsq20)),crf));
431 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
432
433 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
434
435 /* Update potential sum for this i atom from the interaction with this j atom. */
436 velec = _mm_and_pd(velec,cutoff_mask);
437 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
438 velecsum = _mm_add_pd(velecsum,velec);
439
440 fscal = felec;
441
442 fscal = _mm_and_pd(fscal,cutoff_mask);
443
444 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
445
446 /* Calculate temporary vectorial force */
447 tx = _mm_mul_pd(fscal,dx20);
448 ty = _mm_mul_pd(fscal,dy20);
449 tz = _mm_mul_pd(fscal,dz20);
450
451 /* Update vectorial force */
452 fix2 = _mm_add_pd(fix2,tx);
453 fiy2 = _mm_add_pd(fiy2,ty);
454 fiz2 = _mm_add_pd(fiz2,tz);
455
456 fjx0 = _mm_add_pd(fjx0,tx);
457 fjy0 = _mm_add_pd(fjy0,ty);
458 fjz0 = _mm_add_pd(fjz0,tz);
459
460 }
461
462 /**************************
463 * CALCULATE INTERACTIONS *
464 **************************/
465
466 if (gmx_mm_any_lt(rsq30,rcutoff2))
467 {
468
469 /* Compute parameters for interactions between i and j atoms */
470 qq30 = _mm_mul_pd(iq3,jq0);
471
472 /* REACTION-FIELD ELECTROSTATICS */
473 velec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_add_pd(rinv30,_mm_mul_pd(krf,rsq30)),crf));
474 felec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_mul_pd(rinv30,rinvsq30),krf2));
475
476 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
477
478 /* Update potential sum for this i atom from the interaction with this j atom. */
479 velec = _mm_and_pd(velec,cutoff_mask);
480 velec = _mm_unpacklo_pd(velec,_mm_setzero_pd());
481 velecsum = _mm_add_pd(velecsum,velec);
482
483 fscal = felec;
484
485 fscal = _mm_and_pd(fscal,cutoff_mask);
486
487 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
488
489 /* Calculate temporary vectorial force */
490 tx = _mm_mul_pd(fscal,dx30);
491 ty = _mm_mul_pd(fscal,dy30);
492 tz = _mm_mul_pd(fscal,dz30);
493
494 /* Update vectorial force */
495 fix3 = _mm_add_pd(fix3,tx);
496 fiy3 = _mm_add_pd(fiy3,ty);
497 fiz3 = _mm_add_pd(fiz3,tz);
498
499 fjx0 = _mm_add_pd(fjx0,tx);
500 fjy0 = _mm_add_pd(fjy0,ty);
501 fjz0 = _mm_add_pd(fjz0,tz);
502
503 }
504
505 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0);
506
507 /* Inner loop uses 111 flops */
508 }
509
510 /* End of innermost loop */
511
512 gmx_mm_update_iforce_3atom_swizzle_pd(fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
513 f+i_coord_offset+DIM,fshift+i_shift_offset);
514
515 ggid = gid[iidx];
516 /* Update potential energies */
517 gmx_mm_update_1pot_pd(velecsum,kernel_data->energygrp_elec+ggid);
518
519 /* Increment number of inner iterations */
520 inneriter += j_index_end - j_index_start;
521
522 /* Outer loop uses 19 flops */
523 }
524
525 /* Increment number of outer iterations */
526 outeriter += nri;
527
528 /* Update outer/inner flops */
529
530 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W4_VF,outeriter*19 + inneriter*111);
531 }
532 /*
533 * Gromacs nonbonded kernel: nb_kernel_ElecRFCut_VdwNone_GeomW4P1_F_sse4_1_double
534 * Electrostatics interaction: ReactionField
535 * VdW interaction: None
536 * Geometry: Water4-Particle
537 * Calculate force/pot: Force
538 */
539 void
540 nb_kernel_ElecRFCut_VdwNone_GeomW4P1_F_sse4_1_double
541 (t_nblist * gmx_restrict nlist,
542 rvec * gmx_restrict xx,
543 rvec * gmx_restrict ff,
544 t_forcerec * gmx_restrict fr,
545 t_mdatoms * gmx_restrict mdatoms,
546 nb_kernel_data_t gmx_unused * gmx_restrict kernel_data,
547 t_nrnb * gmx_restrict nrnb)
548 {
549 /* Suffixes 0,1,2,3 refer to particle indices for waters in the inner or outer loop, or
550 * just 0 for non-waters.
551 * Suffixes A,B refer to j loop unrolling done with SSE double precision, e.g. for the two different
552 * jnr indices corresponding to data put in the four positions in the SIMD register.
553 */
554 int i_shift_offset,i_coord_offset,outeriter,inneriter;
555 int j_index_start,j_index_end,jidx,nri,inr,ggid,iidx;
556 int jnrA,jnrB;
557 int j_coord_offsetA,j_coord_offsetB;
558 int *iinr,*jindex,*jjnr,*shiftidx,*gid;
559 real rcutoff_scalar;
560 real *shiftvec,*fshift,*x,*f;
561 __m128d tx,ty,tz,fscal,rcutoff,rcutoff2,jidxall;
562 int vdwioffset1;
563 __m128d ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
564 int vdwioffset2;
565 __m128d ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
566 int vdwioffset3;
567 __m128d ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
568 int vdwjidx0A,vdwjidx0B;
569 __m128d jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
570 __m128d dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10;
571 __m128d dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20;
572 __m128d dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30;
573 __m128d velec,felec,velecsum,facel,crf,krf,krf2;
574 real *charge;
575 __m128d dummy_mask,cutoff_mask;
576 __m128d signbit = gmx_mm_castsi128_pd( _mm_set_epi32(0x80000000,0x00000000,0x80000000,0x00000000) );
577 __m128d one = _mm_set1_pd(1.0);
578 __m128d two = _mm_set1_pd(2.0);
579 x = xx[0];
580 f = ff[0];
581
582 nri = nlist->nri;
583 iinr = nlist->iinr;
584 jindex = nlist->jindex;
585 jjnr = nlist->jjnr;
586 shiftidx = nlist->shift;
587 gid = nlist->gid;
588 shiftvec = fr->shift_vec[0];
589 fshift = fr->fshift[0];
590 facel = _mm_set1_pd(fr->epsfac);
591 charge = mdatoms->chargeA;
592 krf = _mm_set1_pd(fr->ic->k_rf);
593 krf2 = _mm_set1_pd(fr->ic->k_rf*2.0);
594 crf = _mm_set1_pd(fr->ic->c_rf);
595
596 /* Setup water-specific parameters */
597 inr = nlist->iinr[0];
598 iq1 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+1]));
599 iq2 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+2]));
600 iq3 = _mm_mul_pd(facel,_mm_set1_pd(charge[inr+3]));
601
602 /* When we use explicit cutoffs the value must be identical for elec and VdW, so use elec as an arbitrary choice */
603 rcutoff_scalar = fr->rcoulomb;
604 rcutoff = _mm_set1_pd(rcutoff_scalar);
605 rcutoff2 = _mm_mul_pd(rcutoff,rcutoff);
606
607 /* Avoid stupid compiler warnings */
608 jnrA = jnrB = 0;
609 j_coord_offsetA = 0;
610 j_coord_offsetB = 0;
611
612 outeriter = 0;
613 inneriter = 0;
614
615 /* Start outer loop over neighborlists */
616 for(iidx=0; iidx<nri; iidx++)
617 {
618 /* Load shift vector for this list */
619 i_shift_offset = DIM*shiftidx[iidx];
620
621 /* Load limits for loop over neighbors */
622 j_index_start = jindex[iidx];
623 j_index_end = jindex[iidx+1];
624
625 /* Get outer coordinate index */
626 inr = iinr[iidx];
627 i_coord_offset = DIM*inr;
628
629 /* Load i particle coords and add shift vector */
630 gmx_mm_load_shift_and_3rvec_broadcast_pd(shiftvec+i_shift_offset,x+i_coord_offset+DIM,
631 &ix1,&iy1,&iz1,&ix2,&iy2,&iz2,&ix3,&iy3,&iz3);
632
633 fix1 = _mm_setzero_pd();
634 fiy1 = _mm_setzero_pd();
635 fiz1 = _mm_setzero_pd();
636 fix2 = _mm_setzero_pd();
637 fiy2 = _mm_setzero_pd();
638 fiz2 = _mm_setzero_pd();
639 fix3 = _mm_setzero_pd();
640 fiy3 = _mm_setzero_pd();
641 fiz3 = _mm_setzero_pd();
642
643 /* Start inner kernel loop */
644 for(jidx=j_index_start; jidx<j_index_end-1; jidx+=2)
645 {
646
647 /* Get j neighbor index, and coordinate index */
648 jnrA = jjnr[jidx];
649 jnrB = jjnr[jidx+1];
650 j_coord_offsetA = DIM*jnrA;
651 j_coord_offsetB = DIM*jnrB;
652
653 /* load j atom coordinates */
654 gmx_mm_load_1rvec_2ptr_swizzle_pd(x+j_coord_offsetA,x+j_coord_offsetB,
655 &jx0,&jy0,&jz0);
656
657 /* Calculate displacement vector */
658 dx10 = _mm_sub_pd(ix1,jx0);
659 dy10 = _mm_sub_pd(iy1,jy0);
660 dz10 = _mm_sub_pd(iz1,jz0);
661 dx20 = _mm_sub_pd(ix2,jx0);
662 dy20 = _mm_sub_pd(iy2,jy0);
663 dz20 = _mm_sub_pd(iz2,jz0);
664 dx30 = _mm_sub_pd(ix3,jx0);
665 dy30 = _mm_sub_pd(iy3,jy0);
666 dz30 = _mm_sub_pd(iz3,jz0);
667
668 /* Calculate squared distance and things based on it */
669 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
670 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
671 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
672
673 rinv10 = gmx_mm_invsqrt_pd(rsq10);
674 rinv20 = gmx_mm_invsqrt_pd(rsq20);
675 rinv30 = gmx_mm_invsqrt_pd(rsq30);
676
677 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
678 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
679 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
680
681 /* Load parameters for j particles */
682 jq0 = gmx_mm_load_2real_swizzle_pd(charge+jnrA+0,charge+jnrB+0);
683
684 fjx0 = _mm_setzero_pd();
685 fjy0 = _mm_setzero_pd();
686 fjz0 = _mm_setzero_pd();
687
688 /**************************
689 * CALCULATE INTERACTIONS *
690 **************************/
691
692 if (gmx_mm_any_lt(rsq10,rcutoff2))
693 {
694
695 /* Compute parameters for interactions between i and j atoms */
696 qq10 = _mm_mul_pd(iq1,jq0);
697
698 /* REACTION-FIELD ELECTROSTATICS */
699 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
700
701 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
702
703 fscal = felec;
704
705 fscal = _mm_and_pd(fscal,cutoff_mask);
706
707 /* Calculate temporary vectorial force */
708 tx = _mm_mul_pd(fscal,dx10);
709 ty = _mm_mul_pd(fscal,dy10);
710 tz = _mm_mul_pd(fscal,dz10);
711
712 /* Update vectorial force */
713 fix1 = _mm_add_pd(fix1,tx);
714 fiy1 = _mm_add_pd(fiy1,ty);
715 fiz1 = _mm_add_pd(fiz1,tz);
716
717 fjx0 = _mm_add_pd(fjx0,tx);
718 fjy0 = _mm_add_pd(fjy0,ty);
719 fjz0 = _mm_add_pd(fjz0,tz);
720
721 }
722
723 /**************************
724 * CALCULATE INTERACTIONS *
725 **************************/
726
727 if (gmx_mm_any_lt(rsq20,rcutoff2))
728 {
729
730 /* Compute parameters for interactions between i and j atoms */
731 qq20 = _mm_mul_pd(iq2,jq0);
732
733 /* REACTION-FIELD ELECTROSTATICS */
734 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
735
736 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
737
738 fscal = felec;
739
740 fscal = _mm_and_pd(fscal,cutoff_mask);
741
742 /* Calculate temporary vectorial force */
743 tx = _mm_mul_pd(fscal,dx20);
744 ty = _mm_mul_pd(fscal,dy20);
745 tz = _mm_mul_pd(fscal,dz20);
746
747 /* Update vectorial force */
748 fix2 = _mm_add_pd(fix2,tx);
749 fiy2 = _mm_add_pd(fiy2,ty);
750 fiz2 = _mm_add_pd(fiz2,tz);
751
752 fjx0 = _mm_add_pd(fjx0,tx);
753 fjy0 = _mm_add_pd(fjy0,ty);
754 fjz0 = _mm_add_pd(fjz0,tz);
755
756 }
757
758 /**************************
759 * CALCULATE INTERACTIONS *
760 **************************/
761
762 if (gmx_mm_any_lt(rsq30,rcutoff2))
763 {
764
765 /* Compute parameters for interactions between i and j atoms */
766 qq30 = _mm_mul_pd(iq3,jq0);
767
768 /* REACTION-FIELD ELECTROSTATICS */
769 felec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_mul_pd(rinv30,rinvsq30),krf2));
770
771 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
772
773 fscal = felec;
774
775 fscal = _mm_and_pd(fscal,cutoff_mask);
776
777 /* Calculate temporary vectorial force */
778 tx = _mm_mul_pd(fscal,dx30);
779 ty = _mm_mul_pd(fscal,dy30);
780 tz = _mm_mul_pd(fscal,dz30);
781
782 /* Update vectorial force */
783 fix3 = _mm_add_pd(fix3,tx);
784 fiy3 = _mm_add_pd(fiy3,ty);
785 fiz3 = _mm_add_pd(fiz3,tz);
786
787 fjx0 = _mm_add_pd(fjx0,tx);
788 fjy0 = _mm_add_pd(fjy0,ty);
789 fjz0 = _mm_add_pd(fjz0,tz);
790
791 }
792
793 gmx_mm_decrement_1rvec_2ptr_swizzle_pd(f+j_coord_offsetA,f+j_coord_offsetB,fjx0,fjy0,fjz0);
794
795 /* Inner loop uses 93 flops */
796 }
797
798 if(jidx<j_index_end)
799 {
800
801 jnrA = jjnr[jidx];
802 j_coord_offsetA = DIM*jnrA;
803
804 /* load j atom coordinates */
805 gmx_mm_load_1rvec_1ptr_swizzle_pd(x+j_coord_offsetA,
806 &jx0,&jy0,&jz0);
807
808 /* Calculate displacement vector */
809 dx10 = _mm_sub_pd(ix1,jx0);
810 dy10 = _mm_sub_pd(iy1,jy0);
811 dz10 = _mm_sub_pd(iz1,jz0);
812 dx20 = _mm_sub_pd(ix2,jx0);
813 dy20 = _mm_sub_pd(iy2,jy0);
814 dz20 = _mm_sub_pd(iz2,jz0);
815 dx30 = _mm_sub_pd(ix3,jx0);
816 dy30 = _mm_sub_pd(iy3,jy0);
817 dz30 = _mm_sub_pd(iz3,jz0);
818
819 /* Calculate squared distance and things based on it */
820 rsq10 = gmx_mm_calc_rsq_pd(dx10,dy10,dz10);
821 rsq20 = gmx_mm_calc_rsq_pd(dx20,dy20,dz20);
822 rsq30 = gmx_mm_calc_rsq_pd(dx30,dy30,dz30);
823
824 rinv10 = gmx_mm_invsqrt_pd(rsq10);
825 rinv20 = gmx_mm_invsqrt_pd(rsq20);
826 rinv30 = gmx_mm_invsqrt_pd(rsq30);
827
828 rinvsq10 = _mm_mul_pd(rinv10,rinv10);
829 rinvsq20 = _mm_mul_pd(rinv20,rinv20);
830 rinvsq30 = _mm_mul_pd(rinv30,rinv30);
831
832 /* Load parameters for j particles */
833 jq0 = _mm_load_sd(charge+jnrA+0);
834
835 fjx0 = _mm_setzero_pd();
836 fjy0 = _mm_setzero_pd();
837 fjz0 = _mm_setzero_pd();
838
839 /**************************
840 * CALCULATE INTERACTIONS *
841 **************************/
842
843 if (gmx_mm_any_lt(rsq10,rcutoff2))
844 {
845
846 /* Compute parameters for interactions between i and j atoms */
847 qq10 = _mm_mul_pd(iq1,jq0);
848
849 /* REACTION-FIELD ELECTROSTATICS */
850 felec = _mm_mul_pd(qq10,_mm_sub_pd(_mm_mul_pd(rinv10,rinvsq10),krf2));
851
852 cutoff_mask = _mm_cmplt_pd(rsq10,rcutoff2);
853
854 fscal = felec;
855
856 fscal = _mm_and_pd(fscal,cutoff_mask);
857
858 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
859
860 /* Calculate temporary vectorial force */
861 tx = _mm_mul_pd(fscal,dx10);
862 ty = _mm_mul_pd(fscal,dy10);
863 tz = _mm_mul_pd(fscal,dz10);
864
865 /* Update vectorial force */
866 fix1 = _mm_add_pd(fix1,tx);
867 fiy1 = _mm_add_pd(fiy1,ty);
868 fiz1 = _mm_add_pd(fiz1,tz);
869
870 fjx0 = _mm_add_pd(fjx0,tx);
871 fjy0 = _mm_add_pd(fjy0,ty);
872 fjz0 = _mm_add_pd(fjz0,tz);
873
874 }
875
876 /**************************
877 * CALCULATE INTERACTIONS *
878 **************************/
879
880 if (gmx_mm_any_lt(rsq20,rcutoff2))
881 {
882
883 /* Compute parameters for interactions between i and j atoms */
884 qq20 = _mm_mul_pd(iq2,jq0);
885
886 /* REACTION-FIELD ELECTROSTATICS */
887 felec = _mm_mul_pd(qq20,_mm_sub_pd(_mm_mul_pd(rinv20,rinvsq20),krf2));
888
889 cutoff_mask = _mm_cmplt_pd(rsq20,rcutoff2);
890
891 fscal = felec;
892
893 fscal = _mm_and_pd(fscal,cutoff_mask);
894
895 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
896
897 /* Calculate temporary vectorial force */
898 tx = _mm_mul_pd(fscal,dx20);
899 ty = _mm_mul_pd(fscal,dy20);
900 tz = _mm_mul_pd(fscal,dz20);
901
902 /* Update vectorial force */
903 fix2 = _mm_add_pd(fix2,tx);
904 fiy2 = _mm_add_pd(fiy2,ty);
905 fiz2 = _mm_add_pd(fiz2,tz);
906
907 fjx0 = _mm_add_pd(fjx0,tx);
908 fjy0 = _mm_add_pd(fjy0,ty);
909 fjz0 = _mm_add_pd(fjz0,tz);
910
911 }
912
913 /**************************
914 * CALCULATE INTERACTIONS *
915 **************************/
916
917 if (gmx_mm_any_lt(rsq30,rcutoff2))
918 {
919
920 /* Compute parameters for interactions between i and j atoms */
921 qq30 = _mm_mul_pd(iq3,jq0);
922
923 /* REACTION-FIELD ELECTROSTATICS */
924 felec = _mm_mul_pd(qq30,_mm_sub_pd(_mm_mul_pd(rinv30,rinvsq30),krf2));
925
926 cutoff_mask = _mm_cmplt_pd(rsq30,rcutoff2);
927
928 fscal = felec;
929
930 fscal = _mm_and_pd(fscal,cutoff_mask);
931
932 fscal = _mm_unpacklo_pd(fscal,_mm_setzero_pd());
933
934 /* Calculate temporary vectorial force */
935 tx = _mm_mul_pd(fscal,dx30);
936 ty = _mm_mul_pd(fscal,dy30);
937 tz = _mm_mul_pd(fscal,dz30);
938
939 /* Update vectorial force */
940 fix3 = _mm_add_pd(fix3,tx);
941 fiy3 = _mm_add_pd(fiy3,ty);
942 fiz3 = _mm_add_pd(fiz3,tz);
943
944 fjx0 = _mm_add_pd(fjx0,tx);
945 fjy0 = _mm_add_pd(fjy0,ty);
946 fjz0 = _mm_add_pd(fjz0,tz);
947
948 }
949
950 gmx_mm_decrement_1rvec_1ptr_swizzle_pd(f+j_coord_offsetA,fjx0,fjy0,fjz0);
951
952 /* Inner loop uses 93 flops */
953 }
954
955 /* End of innermost loop */
956
957 gmx_mm_update_iforce_3atom_swizzle_pd(fix1,fiy1,fiz1,fix2,fiy2,fiz2,fix3,fiy3,fiz3,
958 f+i_coord_offset+DIM,fshift+i_shift_offset);
959
960 /* Increment number of inner iterations */
961 inneriter += j_index_end - j_index_start;
962
963 /* Outer loop uses 18 flops */
964 }
965
966 /* Increment number of outer iterations */
967 outeriter += nri;
968
969 /* Update outer/inner flops */
970
971 inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W4_F,outeriter*18 + inneriter*93);
972 }
The ultimate molecular dynamics simulation package