2 * This file is part of the GROMACS molecular simulation package.
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.
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 c kernel generator.
42 #include "../nb_kernel.h"
43 #include "types/simple.h"
44 #include "gromacs/math/vec.h"
48 * Gromacs nonbonded kernel: nb_kernel_ElecRF_VdwLJ_GeomW4W4_VF_c
49 * Electrostatics interaction: ReactionField
50 * VdW interaction: LennardJones
51 * Geometry: Water4-Water4
52 * Calculate force/pot: PotentialAndForce
55 nb_kernel_ElecRF_VdwLJ_GeomW4W4_VF_c
56 (t_nblist
* gmx_restrict nlist
,
57 rvec
* gmx_restrict xx
,
58 rvec
* gmx_restrict ff
,
59 t_forcerec
* gmx_restrict fr
,
60 t_mdatoms
* gmx_restrict mdatoms
,
61 nb_kernel_data_t gmx_unused
* gmx_restrict kernel_data
,
62 t_nrnb
* gmx_restrict nrnb
)
64 int i_shift_offset
,i_coord_offset
,j_coord_offset
;
65 int j_index_start
,j_index_end
;
66 int nri
,inr
,ggid
,iidx
,jidx
,jnr
,outeriter
,inneriter
;
67 real shX
,shY
,shZ
,tx
,ty
,tz
,fscal
,rcutoff
,rcutoff2
;
68 int *iinr
,*jindex
,*jjnr
,*shiftidx
,*gid
;
69 real
*shiftvec
,*fshift
,*x
,*f
;
71 real ix0
,iy0
,iz0
,fix0
,fiy0
,fiz0
,iq0
,isai0
;
73 real ix1
,iy1
,iz1
,fix1
,fiy1
,fiz1
,iq1
,isai1
;
75 real ix2
,iy2
,iz2
,fix2
,fiy2
,fiz2
,iq2
,isai2
;
77 real ix3
,iy3
,iz3
,fix3
,fiy3
,fiz3
,iq3
,isai3
;
79 real jx0
,jy0
,jz0
,fjx0
,fjy0
,fjz0
,jq0
,isaj0
;
81 real jx1
,jy1
,jz1
,fjx1
,fjy1
,fjz1
,jq1
,isaj1
;
83 real jx2
,jy2
,jz2
,fjx2
,fjy2
,fjz2
,jq2
,isaj2
;
85 real jx3
,jy3
,jz3
,fjx3
,fjy3
,fjz3
,jq3
,isaj3
;
86 real dx00
,dy00
,dz00
,rsq00
,rinv00
,rinvsq00
,r00
,qq00
,c6_00
,c12_00
,cexp1_00
,cexp2_00
;
87 real dx11
,dy11
,dz11
,rsq11
,rinv11
,rinvsq11
,r11
,qq11
,c6_11
,c12_11
,cexp1_11
,cexp2_11
;
88 real dx12
,dy12
,dz12
,rsq12
,rinv12
,rinvsq12
,r12
,qq12
,c6_12
,c12_12
,cexp1_12
,cexp2_12
;
89 real dx13
,dy13
,dz13
,rsq13
,rinv13
,rinvsq13
,r13
,qq13
,c6_13
,c12_13
,cexp1_13
,cexp2_13
;
90 real dx21
,dy21
,dz21
,rsq21
,rinv21
,rinvsq21
,r21
,qq21
,c6_21
,c12_21
,cexp1_21
,cexp2_21
;
91 real dx22
,dy22
,dz22
,rsq22
,rinv22
,rinvsq22
,r22
,qq22
,c6_22
,c12_22
,cexp1_22
,cexp2_22
;
92 real dx23
,dy23
,dz23
,rsq23
,rinv23
,rinvsq23
,r23
,qq23
,c6_23
,c12_23
,cexp1_23
,cexp2_23
;
93 real dx31
,dy31
,dz31
,rsq31
,rinv31
,rinvsq31
,r31
,qq31
,c6_31
,c12_31
,cexp1_31
,cexp2_31
;
94 real dx32
,dy32
,dz32
,rsq32
,rinv32
,rinvsq32
,r32
,qq32
,c6_32
,c12_32
,cexp1_32
,cexp2_32
;
95 real dx33
,dy33
,dz33
,rsq33
,rinv33
,rinvsq33
,r33
,qq33
,c6_33
,c12_33
,cexp1_33
,cexp2_33
;
96 real velec
,felec
,velecsum
,facel
,crf
,krf
,krf2
;
99 real rinvsix
,rvdw
,vvdw
,vvdw6
,vvdw12
,fvdw
,fvdw6
,fvdw12
,vvdwsum
,br
,vvdwexp
,sh_vdw_invrcut6
;
108 jindex
= nlist
->jindex
;
110 shiftidx
= nlist
->shift
;
112 shiftvec
= fr
->shift_vec
[0];
113 fshift
= fr
->fshift
[0];
115 charge
= mdatoms
->chargeA
;
119 nvdwtype
= fr
->ntype
;
121 vdwtype
= mdatoms
->typeA
;
123 /* Setup water-specific parameters */
124 inr
= nlist
->iinr
[0];
125 iq1
= facel
*charge
[inr
+1];
126 iq2
= facel
*charge
[inr
+2];
127 iq3
= facel
*charge
[inr
+3];
128 vdwioffset0
= 2*nvdwtype
*vdwtype
[inr
+0];
133 vdwjidx0
= 2*vdwtype
[inr
+0];
134 c6_00
= vdwparam
[vdwioffset0
+vdwjidx0
];
135 c12_00
= vdwparam
[vdwioffset0
+vdwjidx0
+1];
149 /* Start outer loop over neighborlists */
150 for(iidx
=0; iidx
<nri
; iidx
++)
152 /* Load shift vector for this list */
153 i_shift_offset
= DIM
*shiftidx
[iidx
];
154 shX
= shiftvec
[i_shift_offset
+XX
];
155 shY
= shiftvec
[i_shift_offset
+YY
];
156 shZ
= shiftvec
[i_shift_offset
+ZZ
];
158 /* Load limits for loop over neighbors */
159 j_index_start
= jindex
[iidx
];
160 j_index_end
= jindex
[iidx
+1];
162 /* Get outer coordinate index */
164 i_coord_offset
= DIM
*inr
;
166 /* Load i particle coords and add shift vector */
167 ix0
= shX
+ x
[i_coord_offset
+DIM
*0+XX
];
168 iy0
= shY
+ x
[i_coord_offset
+DIM
*0+YY
];
169 iz0
= shZ
+ x
[i_coord_offset
+DIM
*0+ZZ
];
170 ix1
= shX
+ x
[i_coord_offset
+DIM
*1+XX
];
171 iy1
= shY
+ x
[i_coord_offset
+DIM
*1+YY
];
172 iz1
= shZ
+ x
[i_coord_offset
+DIM
*1+ZZ
];
173 ix2
= shX
+ x
[i_coord_offset
+DIM
*2+XX
];
174 iy2
= shY
+ x
[i_coord_offset
+DIM
*2+YY
];
175 iz2
= shZ
+ x
[i_coord_offset
+DIM
*2+ZZ
];
176 ix3
= shX
+ x
[i_coord_offset
+DIM
*3+XX
];
177 iy3
= shY
+ x
[i_coord_offset
+DIM
*3+YY
];
178 iz3
= shZ
+ x
[i_coord_offset
+DIM
*3+ZZ
];
193 /* Reset potential sums */
197 /* Start inner kernel loop */
198 for(jidx
=j_index_start
; jidx
<j_index_end
; jidx
++)
200 /* Get j neighbor index, and coordinate index */
202 j_coord_offset
= DIM
*jnr
;
204 /* load j atom coordinates */
205 jx0
= x
[j_coord_offset
+DIM
*0+XX
];
206 jy0
= x
[j_coord_offset
+DIM
*0+YY
];
207 jz0
= x
[j_coord_offset
+DIM
*0+ZZ
];
208 jx1
= x
[j_coord_offset
+DIM
*1+XX
];
209 jy1
= x
[j_coord_offset
+DIM
*1+YY
];
210 jz1
= x
[j_coord_offset
+DIM
*1+ZZ
];
211 jx2
= x
[j_coord_offset
+DIM
*2+XX
];
212 jy2
= x
[j_coord_offset
+DIM
*2+YY
];
213 jz2
= x
[j_coord_offset
+DIM
*2+ZZ
];
214 jx3
= x
[j_coord_offset
+DIM
*3+XX
];
215 jy3
= x
[j_coord_offset
+DIM
*3+YY
];
216 jz3
= x
[j_coord_offset
+DIM
*3+ZZ
];
218 /* Calculate displacement vector */
250 /* Calculate squared distance and things based on it */
251 rsq00
= dx00
*dx00
+dy00
*dy00
+dz00
*dz00
;
252 rsq11
= dx11
*dx11
+dy11
*dy11
+dz11
*dz11
;
253 rsq12
= dx12
*dx12
+dy12
*dy12
+dz12
*dz12
;
254 rsq13
= dx13
*dx13
+dy13
*dy13
+dz13
*dz13
;
255 rsq21
= dx21
*dx21
+dy21
*dy21
+dz21
*dz21
;
256 rsq22
= dx22
*dx22
+dy22
*dy22
+dz22
*dz22
;
257 rsq23
= dx23
*dx23
+dy23
*dy23
+dz23
*dz23
;
258 rsq31
= dx31
*dx31
+dy31
*dy31
+dz31
*dz31
;
259 rsq32
= dx32
*dx32
+dy32
*dy32
+dz32
*dz32
;
260 rsq33
= dx33
*dx33
+dy33
*dy33
+dz33
*dz33
;
262 rinv11
= gmx_invsqrt(rsq11
);
263 rinv12
= gmx_invsqrt(rsq12
);
264 rinv13
= gmx_invsqrt(rsq13
);
265 rinv21
= gmx_invsqrt(rsq21
);
266 rinv22
= gmx_invsqrt(rsq22
);
267 rinv23
= gmx_invsqrt(rsq23
);
268 rinv31
= gmx_invsqrt(rsq31
);
269 rinv32
= gmx_invsqrt(rsq32
);
270 rinv33
= gmx_invsqrt(rsq33
);
272 rinvsq00
= 1.0/rsq00
;
273 rinvsq11
= rinv11
*rinv11
;
274 rinvsq12
= rinv12
*rinv12
;
275 rinvsq13
= rinv13
*rinv13
;
276 rinvsq21
= rinv21
*rinv21
;
277 rinvsq22
= rinv22
*rinv22
;
278 rinvsq23
= rinv23
*rinv23
;
279 rinvsq31
= rinv31
*rinv31
;
280 rinvsq32
= rinv32
*rinv32
;
281 rinvsq33
= rinv33
*rinv33
;
283 /**************************
284 * CALCULATE INTERACTIONS *
285 **************************/
287 /* LENNARD-JONES DISPERSION/REPULSION */
289 rinvsix
= rinvsq00
*rinvsq00
*rinvsq00
;
290 vvdw6
= c6_00
*rinvsix
;
291 vvdw12
= c12_00
*rinvsix
*rinvsix
;
292 vvdw
= vvdw12
*(1.0/12.0) - vvdw6
*(1.0/6.0);
293 fvdw
= (vvdw12
-vvdw6
)*rinvsq00
;
295 /* Update potential sums from outer loop */
300 /* Calculate temporary vectorial force */
305 /* Update vectorial force */
309 f
[j_coord_offset
+DIM
*0+XX
] -= tx
;
310 f
[j_coord_offset
+DIM
*0+YY
] -= ty
;
311 f
[j_coord_offset
+DIM
*0+ZZ
] -= tz
;
313 /**************************
314 * CALCULATE INTERACTIONS *
315 **************************/
317 /* REACTION-FIELD ELECTROSTATICS */
318 velec
= qq11
*(rinv11
+krf
*rsq11
-crf
);
319 felec
= qq11
*(rinv11
*rinvsq11
-krf2
);
321 /* Update potential sums from outer loop */
326 /* Calculate temporary vectorial force */
331 /* Update vectorial force */
335 f
[j_coord_offset
+DIM
*1+XX
] -= tx
;
336 f
[j_coord_offset
+DIM
*1+YY
] -= ty
;
337 f
[j_coord_offset
+DIM
*1+ZZ
] -= tz
;
339 /**************************
340 * CALCULATE INTERACTIONS *
341 **************************/
343 /* REACTION-FIELD ELECTROSTATICS */
344 velec
= qq12
*(rinv12
+krf
*rsq12
-crf
);
345 felec
= qq12
*(rinv12
*rinvsq12
-krf2
);
347 /* Update potential sums from outer loop */
352 /* Calculate temporary vectorial force */
357 /* Update vectorial force */
361 f
[j_coord_offset
+DIM
*2+XX
] -= tx
;
362 f
[j_coord_offset
+DIM
*2+YY
] -= ty
;
363 f
[j_coord_offset
+DIM
*2+ZZ
] -= tz
;
365 /**************************
366 * CALCULATE INTERACTIONS *
367 **************************/
369 /* REACTION-FIELD ELECTROSTATICS */
370 velec
= qq13
*(rinv13
+krf
*rsq13
-crf
);
371 felec
= qq13
*(rinv13
*rinvsq13
-krf2
);
373 /* Update potential sums from outer loop */
378 /* Calculate temporary vectorial force */
383 /* Update vectorial force */
387 f
[j_coord_offset
+DIM
*3+XX
] -= tx
;
388 f
[j_coord_offset
+DIM
*3+YY
] -= ty
;
389 f
[j_coord_offset
+DIM
*3+ZZ
] -= tz
;
391 /**************************
392 * CALCULATE INTERACTIONS *
393 **************************/
395 /* REACTION-FIELD ELECTROSTATICS */
396 velec
= qq21
*(rinv21
+krf
*rsq21
-crf
);
397 felec
= qq21
*(rinv21
*rinvsq21
-krf2
);
399 /* Update potential sums from outer loop */
404 /* Calculate temporary vectorial force */
409 /* Update vectorial force */
413 f
[j_coord_offset
+DIM
*1+XX
] -= tx
;
414 f
[j_coord_offset
+DIM
*1+YY
] -= ty
;
415 f
[j_coord_offset
+DIM
*1+ZZ
] -= tz
;
417 /**************************
418 * CALCULATE INTERACTIONS *
419 **************************/
421 /* REACTION-FIELD ELECTROSTATICS */
422 velec
= qq22
*(rinv22
+krf
*rsq22
-crf
);
423 felec
= qq22
*(rinv22
*rinvsq22
-krf2
);
425 /* Update potential sums from outer loop */
430 /* Calculate temporary vectorial force */
435 /* Update vectorial force */
439 f
[j_coord_offset
+DIM
*2+XX
] -= tx
;
440 f
[j_coord_offset
+DIM
*2+YY
] -= ty
;
441 f
[j_coord_offset
+DIM
*2+ZZ
] -= tz
;
443 /**************************
444 * CALCULATE INTERACTIONS *
445 **************************/
447 /* REACTION-FIELD ELECTROSTATICS */
448 velec
= qq23
*(rinv23
+krf
*rsq23
-crf
);
449 felec
= qq23
*(rinv23
*rinvsq23
-krf2
);
451 /* Update potential sums from outer loop */
456 /* Calculate temporary vectorial force */
461 /* Update vectorial force */
465 f
[j_coord_offset
+DIM
*3+XX
] -= tx
;
466 f
[j_coord_offset
+DIM
*3+YY
] -= ty
;
467 f
[j_coord_offset
+DIM
*3+ZZ
] -= tz
;
469 /**************************
470 * CALCULATE INTERACTIONS *
471 **************************/
473 /* REACTION-FIELD ELECTROSTATICS */
474 velec
= qq31
*(rinv31
+krf
*rsq31
-crf
);
475 felec
= qq31
*(rinv31
*rinvsq31
-krf2
);
477 /* Update potential sums from outer loop */
482 /* Calculate temporary vectorial force */
487 /* Update vectorial force */
491 f
[j_coord_offset
+DIM
*1+XX
] -= tx
;
492 f
[j_coord_offset
+DIM
*1+YY
] -= ty
;
493 f
[j_coord_offset
+DIM
*1+ZZ
] -= tz
;
495 /**************************
496 * CALCULATE INTERACTIONS *
497 **************************/
499 /* REACTION-FIELD ELECTROSTATICS */
500 velec
= qq32
*(rinv32
+krf
*rsq32
-crf
);
501 felec
= qq32
*(rinv32
*rinvsq32
-krf2
);
503 /* Update potential sums from outer loop */
508 /* Calculate temporary vectorial force */
513 /* Update vectorial force */
517 f
[j_coord_offset
+DIM
*2+XX
] -= tx
;
518 f
[j_coord_offset
+DIM
*2+YY
] -= ty
;
519 f
[j_coord_offset
+DIM
*2+ZZ
] -= tz
;
521 /**************************
522 * CALCULATE INTERACTIONS *
523 **************************/
525 /* REACTION-FIELD ELECTROSTATICS */
526 velec
= qq33
*(rinv33
+krf
*rsq33
-crf
);
527 felec
= qq33
*(rinv33
*rinvsq33
-krf2
);
529 /* Update potential sums from outer loop */
534 /* Calculate temporary vectorial force */
539 /* Update vectorial force */
543 f
[j_coord_offset
+DIM
*3+XX
] -= tx
;
544 f
[j_coord_offset
+DIM
*3+YY
] -= ty
;
545 f
[j_coord_offset
+DIM
*3+ZZ
] -= tz
;
547 /* Inner loop uses 311 flops */
549 /* End of innermost loop */
552 f
[i_coord_offset
+DIM
*0+XX
] += fix0
;
553 f
[i_coord_offset
+DIM
*0+YY
] += fiy0
;
554 f
[i_coord_offset
+DIM
*0+ZZ
] += fiz0
;
558 f
[i_coord_offset
+DIM
*1+XX
] += fix1
;
559 f
[i_coord_offset
+DIM
*1+YY
] += fiy1
;
560 f
[i_coord_offset
+DIM
*1+ZZ
] += fiz1
;
564 f
[i_coord_offset
+DIM
*2+XX
] += fix2
;
565 f
[i_coord_offset
+DIM
*2+YY
] += fiy2
;
566 f
[i_coord_offset
+DIM
*2+ZZ
] += fiz2
;
570 f
[i_coord_offset
+DIM
*3+XX
] += fix3
;
571 f
[i_coord_offset
+DIM
*3+YY
] += fiy3
;
572 f
[i_coord_offset
+DIM
*3+ZZ
] += fiz3
;
576 fshift
[i_shift_offset
+XX
] += tx
;
577 fshift
[i_shift_offset
+YY
] += ty
;
578 fshift
[i_shift_offset
+ZZ
] += tz
;
581 /* Update potential energies */
582 kernel_data
->energygrp_elec
[ggid
] += velecsum
;
583 kernel_data
->energygrp_vdw
[ggid
] += vvdwsum
;
585 /* Increment number of inner iterations */
586 inneriter
+= j_index_end
- j_index_start
;
588 /* Outer loop uses 41 flops */
591 /* Increment number of outer iterations */
594 /* Update outer/inner flops */
596 inc_nrnb(nrnb
,eNR_NBKERNEL_ELEC_VDW_W4W4_VF
,outeriter
*41 + inneriter
*311);
599 * Gromacs nonbonded kernel: nb_kernel_ElecRF_VdwLJ_GeomW4W4_F_c
600 * Electrostatics interaction: ReactionField
601 * VdW interaction: LennardJones
602 * Geometry: Water4-Water4
603 * Calculate force/pot: Force
606 nb_kernel_ElecRF_VdwLJ_GeomW4W4_F_c
607 (t_nblist
* gmx_restrict nlist
,
608 rvec
* gmx_restrict xx
,
609 rvec
* gmx_restrict ff
,
610 t_forcerec
* gmx_restrict fr
,
611 t_mdatoms
* gmx_restrict mdatoms
,
612 nb_kernel_data_t gmx_unused
* gmx_restrict kernel_data
,
613 t_nrnb
* gmx_restrict nrnb
)
615 int i_shift_offset
,i_coord_offset
,j_coord_offset
;
616 int j_index_start
,j_index_end
;
617 int nri
,inr
,ggid
,iidx
,jidx
,jnr
,outeriter
,inneriter
;
618 real shX
,shY
,shZ
,tx
,ty
,tz
,fscal
,rcutoff
,rcutoff2
;
619 int *iinr
,*jindex
,*jjnr
,*shiftidx
,*gid
;
620 real
*shiftvec
,*fshift
,*x
,*f
;
622 real ix0
,iy0
,iz0
,fix0
,fiy0
,fiz0
,iq0
,isai0
;
624 real ix1
,iy1
,iz1
,fix1
,fiy1
,fiz1
,iq1
,isai1
;
626 real ix2
,iy2
,iz2
,fix2
,fiy2
,fiz2
,iq2
,isai2
;
628 real ix3
,iy3
,iz3
,fix3
,fiy3
,fiz3
,iq3
,isai3
;
630 real jx0
,jy0
,jz0
,fjx0
,fjy0
,fjz0
,jq0
,isaj0
;
632 real jx1
,jy1
,jz1
,fjx1
,fjy1
,fjz1
,jq1
,isaj1
;
634 real jx2
,jy2
,jz2
,fjx2
,fjy2
,fjz2
,jq2
,isaj2
;
636 real jx3
,jy3
,jz3
,fjx3
,fjy3
,fjz3
,jq3
,isaj3
;
637 real dx00
,dy00
,dz00
,rsq00
,rinv00
,rinvsq00
,r00
,qq00
,c6_00
,c12_00
,cexp1_00
,cexp2_00
;
638 real dx11
,dy11
,dz11
,rsq11
,rinv11
,rinvsq11
,r11
,qq11
,c6_11
,c12_11
,cexp1_11
,cexp2_11
;
639 real dx12
,dy12
,dz12
,rsq12
,rinv12
,rinvsq12
,r12
,qq12
,c6_12
,c12_12
,cexp1_12
,cexp2_12
;
640 real dx13
,dy13
,dz13
,rsq13
,rinv13
,rinvsq13
,r13
,qq13
,c6_13
,c12_13
,cexp1_13
,cexp2_13
;
641 real dx21
,dy21
,dz21
,rsq21
,rinv21
,rinvsq21
,r21
,qq21
,c6_21
,c12_21
,cexp1_21
,cexp2_21
;
642 real dx22
,dy22
,dz22
,rsq22
,rinv22
,rinvsq22
,r22
,qq22
,c6_22
,c12_22
,cexp1_22
,cexp2_22
;
643 real dx23
,dy23
,dz23
,rsq23
,rinv23
,rinvsq23
,r23
,qq23
,c6_23
,c12_23
,cexp1_23
,cexp2_23
;
644 real dx31
,dy31
,dz31
,rsq31
,rinv31
,rinvsq31
,r31
,qq31
,c6_31
,c12_31
,cexp1_31
,cexp2_31
;
645 real dx32
,dy32
,dz32
,rsq32
,rinv32
,rinvsq32
,r32
,qq32
,c6_32
,c12_32
,cexp1_32
,cexp2_32
;
646 real dx33
,dy33
,dz33
,rsq33
,rinv33
,rinvsq33
,r33
,qq33
,c6_33
,c12_33
,cexp1_33
,cexp2_33
;
647 real velec
,felec
,velecsum
,facel
,crf
,krf
,krf2
;
650 real rinvsix
,rvdw
,vvdw
,vvdw6
,vvdw12
,fvdw
,fvdw6
,fvdw12
,vvdwsum
,br
,vvdwexp
,sh_vdw_invrcut6
;
659 jindex
= nlist
->jindex
;
661 shiftidx
= nlist
->shift
;
663 shiftvec
= fr
->shift_vec
[0];
664 fshift
= fr
->fshift
[0];
666 charge
= mdatoms
->chargeA
;
670 nvdwtype
= fr
->ntype
;
672 vdwtype
= mdatoms
->typeA
;
674 /* Setup water-specific parameters */
675 inr
= nlist
->iinr
[0];
676 iq1
= facel
*charge
[inr
+1];
677 iq2
= facel
*charge
[inr
+2];
678 iq3
= facel
*charge
[inr
+3];
679 vdwioffset0
= 2*nvdwtype
*vdwtype
[inr
+0];
684 vdwjidx0
= 2*vdwtype
[inr
+0];
685 c6_00
= vdwparam
[vdwioffset0
+vdwjidx0
];
686 c12_00
= vdwparam
[vdwioffset0
+vdwjidx0
+1];
700 /* Start outer loop over neighborlists */
701 for(iidx
=0; iidx
<nri
; iidx
++)
703 /* Load shift vector for this list */
704 i_shift_offset
= DIM
*shiftidx
[iidx
];
705 shX
= shiftvec
[i_shift_offset
+XX
];
706 shY
= shiftvec
[i_shift_offset
+YY
];
707 shZ
= shiftvec
[i_shift_offset
+ZZ
];
709 /* Load limits for loop over neighbors */
710 j_index_start
= jindex
[iidx
];
711 j_index_end
= jindex
[iidx
+1];
713 /* Get outer coordinate index */
715 i_coord_offset
= DIM
*inr
;
717 /* Load i particle coords and add shift vector */
718 ix0
= shX
+ x
[i_coord_offset
+DIM
*0+XX
];
719 iy0
= shY
+ x
[i_coord_offset
+DIM
*0+YY
];
720 iz0
= shZ
+ x
[i_coord_offset
+DIM
*0+ZZ
];
721 ix1
= shX
+ x
[i_coord_offset
+DIM
*1+XX
];
722 iy1
= shY
+ x
[i_coord_offset
+DIM
*1+YY
];
723 iz1
= shZ
+ x
[i_coord_offset
+DIM
*1+ZZ
];
724 ix2
= shX
+ x
[i_coord_offset
+DIM
*2+XX
];
725 iy2
= shY
+ x
[i_coord_offset
+DIM
*2+YY
];
726 iz2
= shZ
+ x
[i_coord_offset
+DIM
*2+ZZ
];
727 ix3
= shX
+ x
[i_coord_offset
+DIM
*3+XX
];
728 iy3
= shY
+ x
[i_coord_offset
+DIM
*3+YY
];
729 iz3
= shZ
+ x
[i_coord_offset
+DIM
*3+ZZ
];
744 /* Start inner kernel loop */
745 for(jidx
=j_index_start
; jidx
<j_index_end
; jidx
++)
747 /* Get j neighbor index, and coordinate index */
749 j_coord_offset
= DIM
*jnr
;
751 /* load j atom coordinates */
752 jx0
= x
[j_coord_offset
+DIM
*0+XX
];
753 jy0
= x
[j_coord_offset
+DIM
*0+YY
];
754 jz0
= x
[j_coord_offset
+DIM
*0+ZZ
];
755 jx1
= x
[j_coord_offset
+DIM
*1+XX
];
756 jy1
= x
[j_coord_offset
+DIM
*1+YY
];
757 jz1
= x
[j_coord_offset
+DIM
*1+ZZ
];
758 jx2
= x
[j_coord_offset
+DIM
*2+XX
];
759 jy2
= x
[j_coord_offset
+DIM
*2+YY
];
760 jz2
= x
[j_coord_offset
+DIM
*2+ZZ
];
761 jx3
= x
[j_coord_offset
+DIM
*3+XX
];
762 jy3
= x
[j_coord_offset
+DIM
*3+YY
];
763 jz3
= x
[j_coord_offset
+DIM
*3+ZZ
];
765 /* Calculate displacement vector */
797 /* Calculate squared distance and things based on it */
798 rsq00
= dx00
*dx00
+dy00
*dy00
+dz00
*dz00
;
799 rsq11
= dx11
*dx11
+dy11
*dy11
+dz11
*dz11
;
800 rsq12
= dx12
*dx12
+dy12
*dy12
+dz12
*dz12
;
801 rsq13
= dx13
*dx13
+dy13
*dy13
+dz13
*dz13
;
802 rsq21
= dx21
*dx21
+dy21
*dy21
+dz21
*dz21
;
803 rsq22
= dx22
*dx22
+dy22
*dy22
+dz22
*dz22
;
804 rsq23
= dx23
*dx23
+dy23
*dy23
+dz23
*dz23
;
805 rsq31
= dx31
*dx31
+dy31
*dy31
+dz31
*dz31
;
806 rsq32
= dx32
*dx32
+dy32
*dy32
+dz32
*dz32
;
807 rsq33
= dx33
*dx33
+dy33
*dy33
+dz33
*dz33
;
809 rinv11
= gmx_invsqrt(rsq11
);
810 rinv12
= gmx_invsqrt(rsq12
);
811 rinv13
= gmx_invsqrt(rsq13
);
812 rinv21
= gmx_invsqrt(rsq21
);
813 rinv22
= gmx_invsqrt(rsq22
);
814 rinv23
= gmx_invsqrt(rsq23
);
815 rinv31
= gmx_invsqrt(rsq31
);
816 rinv32
= gmx_invsqrt(rsq32
);
817 rinv33
= gmx_invsqrt(rsq33
);
819 rinvsq00
= 1.0/rsq00
;
820 rinvsq11
= rinv11
*rinv11
;
821 rinvsq12
= rinv12
*rinv12
;
822 rinvsq13
= rinv13
*rinv13
;
823 rinvsq21
= rinv21
*rinv21
;
824 rinvsq22
= rinv22
*rinv22
;
825 rinvsq23
= rinv23
*rinv23
;
826 rinvsq31
= rinv31
*rinv31
;
827 rinvsq32
= rinv32
*rinv32
;
828 rinvsq33
= rinv33
*rinv33
;
830 /**************************
831 * CALCULATE INTERACTIONS *
832 **************************/
834 /* LENNARD-JONES DISPERSION/REPULSION */
836 rinvsix
= rinvsq00
*rinvsq00
*rinvsq00
;
837 fvdw
= (c12_00
*rinvsix
-c6_00
)*rinvsix
*rinvsq00
;
841 /* Calculate temporary vectorial force */
846 /* Update vectorial force */
850 f
[j_coord_offset
+DIM
*0+XX
] -= tx
;
851 f
[j_coord_offset
+DIM
*0+YY
] -= ty
;
852 f
[j_coord_offset
+DIM
*0+ZZ
] -= tz
;
854 /**************************
855 * CALCULATE INTERACTIONS *
856 **************************/
858 /* REACTION-FIELD ELECTROSTATICS */
859 felec
= qq11
*(rinv11
*rinvsq11
-krf2
);
863 /* Calculate temporary vectorial force */
868 /* Update vectorial force */
872 f
[j_coord_offset
+DIM
*1+XX
] -= tx
;
873 f
[j_coord_offset
+DIM
*1+YY
] -= ty
;
874 f
[j_coord_offset
+DIM
*1+ZZ
] -= tz
;
876 /**************************
877 * CALCULATE INTERACTIONS *
878 **************************/
880 /* REACTION-FIELD ELECTROSTATICS */
881 felec
= qq12
*(rinv12
*rinvsq12
-krf2
);
885 /* Calculate temporary vectorial force */
890 /* Update vectorial force */
894 f
[j_coord_offset
+DIM
*2+XX
] -= tx
;
895 f
[j_coord_offset
+DIM
*2+YY
] -= ty
;
896 f
[j_coord_offset
+DIM
*2+ZZ
] -= tz
;
898 /**************************
899 * CALCULATE INTERACTIONS *
900 **************************/
902 /* REACTION-FIELD ELECTROSTATICS */
903 felec
= qq13
*(rinv13
*rinvsq13
-krf2
);
907 /* Calculate temporary vectorial force */
912 /* Update vectorial force */
916 f
[j_coord_offset
+DIM
*3+XX
] -= tx
;
917 f
[j_coord_offset
+DIM
*3+YY
] -= ty
;
918 f
[j_coord_offset
+DIM
*3+ZZ
] -= tz
;
920 /**************************
921 * CALCULATE INTERACTIONS *
922 **************************/
924 /* REACTION-FIELD ELECTROSTATICS */
925 felec
= qq21
*(rinv21
*rinvsq21
-krf2
);
929 /* Calculate temporary vectorial force */
934 /* Update vectorial force */
938 f
[j_coord_offset
+DIM
*1+XX
] -= tx
;
939 f
[j_coord_offset
+DIM
*1+YY
] -= ty
;
940 f
[j_coord_offset
+DIM
*1+ZZ
] -= tz
;
942 /**************************
943 * CALCULATE INTERACTIONS *
944 **************************/
946 /* REACTION-FIELD ELECTROSTATICS */
947 felec
= qq22
*(rinv22
*rinvsq22
-krf2
);
951 /* Calculate temporary vectorial force */
956 /* Update vectorial force */
960 f
[j_coord_offset
+DIM
*2+XX
] -= tx
;
961 f
[j_coord_offset
+DIM
*2+YY
] -= ty
;
962 f
[j_coord_offset
+DIM
*2+ZZ
] -= tz
;
964 /**************************
965 * CALCULATE INTERACTIONS *
966 **************************/
968 /* REACTION-FIELD ELECTROSTATICS */
969 felec
= qq23
*(rinv23
*rinvsq23
-krf2
);
973 /* Calculate temporary vectorial force */
978 /* Update vectorial force */
982 f
[j_coord_offset
+DIM
*3+XX
] -= tx
;
983 f
[j_coord_offset
+DIM
*3+YY
] -= ty
;
984 f
[j_coord_offset
+DIM
*3+ZZ
] -= tz
;
986 /**************************
987 * CALCULATE INTERACTIONS *
988 **************************/
990 /* REACTION-FIELD ELECTROSTATICS */
991 felec
= qq31
*(rinv31
*rinvsq31
-krf2
);
995 /* Calculate temporary vectorial force */
1000 /* Update vectorial force */
1004 f
[j_coord_offset
+DIM
*1+XX
] -= tx
;
1005 f
[j_coord_offset
+DIM
*1+YY
] -= ty
;
1006 f
[j_coord_offset
+DIM
*1+ZZ
] -= tz
;
1008 /**************************
1009 * CALCULATE INTERACTIONS *
1010 **************************/
1012 /* REACTION-FIELD ELECTROSTATICS */
1013 felec
= qq32
*(rinv32
*rinvsq32
-krf2
);
1017 /* Calculate temporary vectorial force */
1022 /* Update vectorial force */
1026 f
[j_coord_offset
+DIM
*2+XX
] -= tx
;
1027 f
[j_coord_offset
+DIM
*2+YY
] -= ty
;
1028 f
[j_coord_offset
+DIM
*2+ZZ
] -= tz
;
1030 /**************************
1031 * CALCULATE INTERACTIONS *
1032 **************************/
1034 /* REACTION-FIELD ELECTROSTATICS */
1035 felec
= qq33
*(rinv33
*rinvsq33
-krf2
);
1039 /* Calculate temporary vectorial force */
1044 /* Update vectorial force */
1048 f
[j_coord_offset
+DIM
*3+XX
] -= tx
;
1049 f
[j_coord_offset
+DIM
*3+YY
] -= ty
;
1050 f
[j_coord_offset
+DIM
*3+ZZ
] -= tz
;
1052 /* Inner loop uses 261 flops */
1054 /* End of innermost loop */
1057 f
[i_coord_offset
+DIM
*0+XX
] += fix0
;
1058 f
[i_coord_offset
+DIM
*0+YY
] += fiy0
;
1059 f
[i_coord_offset
+DIM
*0+ZZ
] += fiz0
;
1063 f
[i_coord_offset
+DIM
*1+XX
] += fix1
;
1064 f
[i_coord_offset
+DIM
*1+YY
] += fiy1
;
1065 f
[i_coord_offset
+DIM
*1+ZZ
] += fiz1
;
1069 f
[i_coord_offset
+DIM
*2+XX
] += fix2
;
1070 f
[i_coord_offset
+DIM
*2+YY
] += fiy2
;
1071 f
[i_coord_offset
+DIM
*2+ZZ
] += fiz2
;
1075 f
[i_coord_offset
+DIM
*3+XX
] += fix3
;
1076 f
[i_coord_offset
+DIM
*3+YY
] += fiy3
;
1077 f
[i_coord_offset
+DIM
*3+ZZ
] += fiz3
;
1081 fshift
[i_shift_offset
+XX
] += tx
;
1082 fshift
[i_shift_offset
+YY
] += ty
;
1083 fshift
[i_shift_offset
+ZZ
] += tz
;
1085 /* Increment number of inner iterations */
1086 inneriter
+= j_index_end
- j_index_start
;
1088 /* Outer loop uses 39 flops */
1091 /* Increment number of outer iterations */
1094 /* Update outer/inner flops */
1096 inc_nrnb(nrnb
,eNR_NBKERNEL_ELEC_VDW_W4W4_F
,outeriter
*39 + inneriter
*261);