ebb7ac0ebfb0ea8ff08efec78d09c726272975b8
[alexxy/gromacs.git] / src / gmxlib / nonbonded / nb_kernel_c / nb_kernel_ElecCoul_VdwNone_GeomW4P1_c.c
1 /*
2  * Note: this file was generated by the Gromacs c kernel generator.
3  *
4  *                This source code is part of
5  *
6  *                 G   R   O   M   A   C   S
7  *
8  * Copyright (c) 2001-2012, The GROMACS Development Team
9  *
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
13  *
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
17  * later version.
18  *
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.
21  */
22 #ifdef HAVE_CONFIG_H
23 #include <config.h>
24 #endif
25
26 #include <math.h>
27
28 #include "../nb_kernel.h"
29 #include "types/simple.h"
30 #include "vec.h"
31 #include "nrnb.h"
32
33 /*
34  * Gromacs nonbonded kernel:   nb_kernel_ElecCoul_VdwNone_GeomW4P1_VF_c
35  * Electrostatics interaction: Coulomb
36  * VdW interaction:            None
37  * Geometry:                   Water4-Particle
38  * Calculate force/pot:        PotentialAndForce
39  */
40 void
41 nb_kernel_ElecCoul_VdwNone_GeomW4P1_VF_c
42                     (t_nblist * gmx_restrict                nlist,
43                      rvec * gmx_restrict                    xx,
44                      rvec * gmx_restrict                    ff,
45                      t_forcerec * gmx_restrict              fr,
46                      t_mdatoms * gmx_restrict               mdatoms,
47                      nb_kernel_data_t * gmx_restrict        kernel_data,
48                      t_nrnb * gmx_restrict                  nrnb)
49 {
50     int              i_shift_offset,i_coord_offset,j_coord_offset;
51     int              j_index_start,j_index_end;
52     int              nri,inr,ggid,iidx,jidx,jnr,outeriter,inneriter;
53     real             shX,shY,shZ,tx,ty,tz,fscal,rcutoff,rcutoff2;
54     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
55     real             *shiftvec,*fshift,*x,*f;
56     int              vdwioffset1;
57     real             ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
58     int              vdwioffset2;
59     real             ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
60     int              vdwioffset3;
61     real             ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
62     int              vdwjidx0;
63     real             jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
64     real             dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10,cexp1_10,cexp2_10;
65     real             dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20,cexp1_20,cexp2_20;
66     real             dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30,cexp1_30,cexp2_30;
67     real             velec,felec,velecsum,facel,crf,krf,krf2;
68     real             *charge;
69
70     x                = xx[0];
71     f                = ff[0];
72
73     nri              = nlist->nri;
74     iinr             = nlist->iinr;
75     jindex           = nlist->jindex;
76     jjnr             = nlist->jjnr;
77     shiftidx         = nlist->shift;
78     gid              = nlist->gid;
79     shiftvec         = fr->shift_vec[0];
80     fshift           = fr->fshift[0];
81     facel            = fr->epsfac;
82     charge           = mdatoms->chargeA;
83
84     /* Setup water-specific parameters */
85     inr              = nlist->iinr[0];
86     iq1              = facel*charge[inr+1];
87     iq2              = facel*charge[inr+2];
88     iq3              = facel*charge[inr+3];
89
90     outeriter        = 0;
91     inneriter        = 0;
92
93     /* Start outer loop over neighborlists */
94     for(iidx=0; iidx<nri; iidx++)
95     {
96         /* Load shift vector for this list */
97         i_shift_offset   = DIM*shiftidx[iidx];
98         shX              = shiftvec[i_shift_offset+XX];
99         shY              = shiftvec[i_shift_offset+YY];
100         shZ              = shiftvec[i_shift_offset+ZZ];
101
102         /* Load limits for loop over neighbors */
103         j_index_start    = jindex[iidx];
104         j_index_end      = jindex[iidx+1];
105
106         /* Get outer coordinate index */
107         inr              = iinr[iidx];
108         i_coord_offset   = DIM*inr;
109
110         /* Load i particle coords and add shift vector */
111         ix1              = shX + x[i_coord_offset+DIM*1+XX];
112         iy1              = shY + x[i_coord_offset+DIM*1+YY];
113         iz1              = shZ + x[i_coord_offset+DIM*1+ZZ];
114         ix2              = shX + x[i_coord_offset+DIM*2+XX];
115         iy2              = shY + x[i_coord_offset+DIM*2+YY];
116         iz2              = shZ + x[i_coord_offset+DIM*2+ZZ];
117         ix3              = shX + x[i_coord_offset+DIM*3+XX];
118         iy3              = shY + x[i_coord_offset+DIM*3+YY];
119         iz3              = shZ + x[i_coord_offset+DIM*3+ZZ];
120
121         fix1             = 0.0;
122         fiy1             = 0.0;
123         fiz1             = 0.0;
124         fix2             = 0.0;
125         fiy2             = 0.0;
126         fiz2             = 0.0;
127         fix3             = 0.0;
128         fiy3             = 0.0;
129         fiz3             = 0.0;
130
131         /* Reset potential sums */
132         velecsum         = 0.0;
133
134         /* Start inner kernel loop */
135         for(jidx=j_index_start; jidx<j_index_end; jidx++)
136         {
137             /* Get j neighbor index, and coordinate index */
138             jnr              = jjnr[jidx];
139             j_coord_offset   = DIM*jnr;
140
141             /* load j atom coordinates */
142             jx0              = x[j_coord_offset+DIM*0+XX];
143             jy0              = x[j_coord_offset+DIM*0+YY];
144             jz0              = x[j_coord_offset+DIM*0+ZZ];
145
146             /* Calculate displacement vector */
147             dx10             = ix1 - jx0;
148             dy10             = iy1 - jy0;
149             dz10             = iz1 - jz0;
150             dx20             = ix2 - jx0;
151             dy20             = iy2 - jy0;
152             dz20             = iz2 - jz0;
153             dx30             = ix3 - jx0;
154             dy30             = iy3 - jy0;
155             dz30             = iz3 - jz0;
156
157             /* Calculate squared distance and things based on it */
158             rsq10            = dx10*dx10+dy10*dy10+dz10*dz10;
159             rsq20            = dx20*dx20+dy20*dy20+dz20*dz20;
160             rsq30            = dx30*dx30+dy30*dy30+dz30*dz30;
161
162             rinv10           = gmx_invsqrt(rsq10);
163             rinv20           = gmx_invsqrt(rsq20);
164             rinv30           = gmx_invsqrt(rsq30);
165
166             rinvsq10         = rinv10*rinv10;
167             rinvsq20         = rinv20*rinv20;
168             rinvsq30         = rinv30*rinv30;
169
170             /* Load parameters for j particles */
171             jq0              = charge[jnr+0];
172
173             /**************************
174              * CALCULATE INTERACTIONS *
175              **************************/
176
177             qq10             = iq1*jq0;
178
179             /* COULOMB ELECTROSTATICS */
180             velec            = qq10*rinv10;
181             felec            = velec*rinvsq10;
182
183             /* Update potential sums from outer loop */
184             velecsum        += velec;
185
186             fscal            = felec;
187
188             /* Calculate temporary vectorial force */
189             tx               = fscal*dx10;
190             ty               = fscal*dy10;
191             tz               = fscal*dz10;
192
193             /* Update vectorial force */
194             fix1            += tx;
195             fiy1            += ty;
196             fiz1            += tz;
197             f[j_coord_offset+DIM*0+XX] -= tx;
198             f[j_coord_offset+DIM*0+YY] -= ty;
199             f[j_coord_offset+DIM*0+ZZ] -= tz;
200
201             /**************************
202              * CALCULATE INTERACTIONS *
203              **************************/
204
205             qq20             = iq2*jq0;
206
207             /* COULOMB ELECTROSTATICS */
208             velec            = qq20*rinv20;
209             felec            = velec*rinvsq20;
210
211             /* Update potential sums from outer loop */
212             velecsum        += velec;
213
214             fscal            = felec;
215
216             /* Calculate temporary vectorial force */
217             tx               = fscal*dx20;
218             ty               = fscal*dy20;
219             tz               = fscal*dz20;
220
221             /* Update vectorial force */
222             fix2            += tx;
223             fiy2            += ty;
224             fiz2            += tz;
225             f[j_coord_offset+DIM*0+XX] -= tx;
226             f[j_coord_offset+DIM*0+YY] -= ty;
227             f[j_coord_offset+DIM*0+ZZ] -= tz;
228
229             /**************************
230              * CALCULATE INTERACTIONS *
231              **************************/
232
233             qq30             = iq3*jq0;
234
235             /* COULOMB ELECTROSTATICS */
236             velec            = qq30*rinv30;
237             felec            = velec*rinvsq30;
238
239             /* Update potential sums from outer loop */
240             velecsum        += velec;
241
242             fscal            = felec;
243
244             /* Calculate temporary vectorial force */
245             tx               = fscal*dx30;
246             ty               = fscal*dy30;
247             tz               = fscal*dz30;
248
249             /* Update vectorial force */
250             fix3            += tx;
251             fiy3            += ty;
252             fiz3            += tz;
253             f[j_coord_offset+DIM*0+XX] -= tx;
254             f[j_coord_offset+DIM*0+YY] -= ty;
255             f[j_coord_offset+DIM*0+ZZ] -= tz;
256
257             /* Inner loop uses 84 flops */
258         }
259         /* End of innermost loop */
260
261         tx = ty = tz = 0;
262         f[i_coord_offset+DIM*1+XX] += fix1;
263         f[i_coord_offset+DIM*1+YY] += fiy1;
264         f[i_coord_offset+DIM*1+ZZ] += fiz1;
265         tx                         += fix1;
266         ty                         += fiy1;
267         tz                         += fiz1;
268         f[i_coord_offset+DIM*2+XX] += fix2;
269         f[i_coord_offset+DIM*2+YY] += fiy2;
270         f[i_coord_offset+DIM*2+ZZ] += fiz2;
271         tx                         += fix2;
272         ty                         += fiy2;
273         tz                         += fiz2;
274         f[i_coord_offset+DIM*3+XX] += fix3;
275         f[i_coord_offset+DIM*3+YY] += fiy3;
276         f[i_coord_offset+DIM*3+ZZ] += fiz3;
277         tx                         += fix3;
278         ty                         += fiy3;
279         tz                         += fiz3;
280         fshift[i_shift_offset+XX]  += tx;
281         fshift[i_shift_offset+YY]  += ty;
282         fshift[i_shift_offset+ZZ]  += tz;
283
284         ggid                        = gid[iidx];
285         /* Update potential energies */
286         kernel_data->energygrp_elec[ggid] += velecsum;
287
288         /* Increment number of inner iterations */
289         inneriter                  += j_index_end - j_index_start;
290
291         /* Outer loop uses 31 flops */
292     }
293
294     /* Increment number of outer iterations */
295     outeriter        += nri;
296
297     /* Update outer/inner flops */
298
299     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W4_VF,outeriter*31 + inneriter*84);
300 }
301 /*
302  * Gromacs nonbonded kernel:   nb_kernel_ElecCoul_VdwNone_GeomW4P1_F_c
303  * Electrostatics interaction: Coulomb
304  * VdW interaction:            None
305  * Geometry:                   Water4-Particle
306  * Calculate force/pot:        Force
307  */
308 void
309 nb_kernel_ElecCoul_VdwNone_GeomW4P1_F_c
310                     (t_nblist * gmx_restrict                nlist,
311                      rvec * gmx_restrict                    xx,
312                      rvec * gmx_restrict                    ff,
313                      t_forcerec * gmx_restrict              fr,
314                      t_mdatoms * gmx_restrict               mdatoms,
315                      nb_kernel_data_t * gmx_restrict        kernel_data,
316                      t_nrnb * gmx_restrict                  nrnb)
317 {
318     int              i_shift_offset,i_coord_offset,j_coord_offset;
319     int              j_index_start,j_index_end;
320     int              nri,inr,ggid,iidx,jidx,jnr,outeriter,inneriter;
321     real             shX,shY,shZ,tx,ty,tz,fscal,rcutoff,rcutoff2;
322     int              *iinr,*jindex,*jjnr,*shiftidx,*gid;
323     real             *shiftvec,*fshift,*x,*f;
324     int              vdwioffset1;
325     real             ix1,iy1,iz1,fix1,fiy1,fiz1,iq1,isai1;
326     int              vdwioffset2;
327     real             ix2,iy2,iz2,fix2,fiy2,fiz2,iq2,isai2;
328     int              vdwioffset3;
329     real             ix3,iy3,iz3,fix3,fiy3,fiz3,iq3,isai3;
330     int              vdwjidx0;
331     real             jx0,jy0,jz0,fjx0,fjy0,fjz0,jq0,isaj0;
332     real             dx10,dy10,dz10,rsq10,rinv10,rinvsq10,r10,qq10,c6_10,c12_10,cexp1_10,cexp2_10;
333     real             dx20,dy20,dz20,rsq20,rinv20,rinvsq20,r20,qq20,c6_20,c12_20,cexp1_20,cexp2_20;
334     real             dx30,dy30,dz30,rsq30,rinv30,rinvsq30,r30,qq30,c6_30,c12_30,cexp1_30,cexp2_30;
335     real             velec,felec,velecsum,facel,crf,krf,krf2;
336     real             *charge;
337
338     x                = xx[0];
339     f                = ff[0];
340
341     nri              = nlist->nri;
342     iinr             = nlist->iinr;
343     jindex           = nlist->jindex;
344     jjnr             = nlist->jjnr;
345     shiftidx         = nlist->shift;
346     gid              = nlist->gid;
347     shiftvec         = fr->shift_vec[0];
348     fshift           = fr->fshift[0];
349     facel            = fr->epsfac;
350     charge           = mdatoms->chargeA;
351
352     /* Setup water-specific parameters */
353     inr              = nlist->iinr[0];
354     iq1              = facel*charge[inr+1];
355     iq2              = facel*charge[inr+2];
356     iq3              = facel*charge[inr+3];
357
358     outeriter        = 0;
359     inneriter        = 0;
360
361     /* Start outer loop over neighborlists */
362     for(iidx=0; iidx<nri; iidx++)
363     {
364         /* Load shift vector for this list */
365         i_shift_offset   = DIM*shiftidx[iidx];
366         shX              = shiftvec[i_shift_offset+XX];
367         shY              = shiftvec[i_shift_offset+YY];
368         shZ              = shiftvec[i_shift_offset+ZZ];
369
370         /* Load limits for loop over neighbors */
371         j_index_start    = jindex[iidx];
372         j_index_end      = jindex[iidx+1];
373
374         /* Get outer coordinate index */
375         inr              = iinr[iidx];
376         i_coord_offset   = DIM*inr;
377
378         /* Load i particle coords and add shift vector */
379         ix1              = shX + x[i_coord_offset+DIM*1+XX];
380         iy1              = shY + x[i_coord_offset+DIM*1+YY];
381         iz1              = shZ + x[i_coord_offset+DIM*1+ZZ];
382         ix2              = shX + x[i_coord_offset+DIM*2+XX];
383         iy2              = shY + x[i_coord_offset+DIM*2+YY];
384         iz2              = shZ + x[i_coord_offset+DIM*2+ZZ];
385         ix3              = shX + x[i_coord_offset+DIM*3+XX];
386         iy3              = shY + x[i_coord_offset+DIM*3+YY];
387         iz3              = shZ + x[i_coord_offset+DIM*3+ZZ];
388
389         fix1             = 0.0;
390         fiy1             = 0.0;
391         fiz1             = 0.0;
392         fix2             = 0.0;
393         fiy2             = 0.0;
394         fiz2             = 0.0;
395         fix3             = 0.0;
396         fiy3             = 0.0;
397         fiz3             = 0.0;
398
399         /* Start inner kernel loop */
400         for(jidx=j_index_start; jidx<j_index_end; jidx++)
401         {
402             /* Get j neighbor index, and coordinate index */
403             jnr              = jjnr[jidx];
404             j_coord_offset   = DIM*jnr;
405
406             /* load j atom coordinates */
407             jx0              = x[j_coord_offset+DIM*0+XX];
408             jy0              = x[j_coord_offset+DIM*0+YY];
409             jz0              = x[j_coord_offset+DIM*0+ZZ];
410
411             /* Calculate displacement vector */
412             dx10             = ix1 - jx0;
413             dy10             = iy1 - jy0;
414             dz10             = iz1 - jz0;
415             dx20             = ix2 - jx0;
416             dy20             = iy2 - jy0;
417             dz20             = iz2 - jz0;
418             dx30             = ix3 - jx0;
419             dy30             = iy3 - jy0;
420             dz30             = iz3 - jz0;
421
422             /* Calculate squared distance and things based on it */
423             rsq10            = dx10*dx10+dy10*dy10+dz10*dz10;
424             rsq20            = dx20*dx20+dy20*dy20+dz20*dz20;
425             rsq30            = dx30*dx30+dy30*dy30+dz30*dz30;
426
427             rinv10           = gmx_invsqrt(rsq10);
428             rinv20           = gmx_invsqrt(rsq20);
429             rinv30           = gmx_invsqrt(rsq30);
430
431             rinvsq10         = rinv10*rinv10;
432             rinvsq20         = rinv20*rinv20;
433             rinvsq30         = rinv30*rinv30;
434
435             /* Load parameters for j particles */
436             jq0              = charge[jnr+0];
437
438             /**************************
439              * CALCULATE INTERACTIONS *
440              **************************/
441
442             qq10             = iq1*jq0;
443
444             /* COULOMB ELECTROSTATICS */
445             velec            = qq10*rinv10;
446             felec            = velec*rinvsq10;
447
448             fscal            = felec;
449
450             /* Calculate temporary vectorial force */
451             tx               = fscal*dx10;
452             ty               = fscal*dy10;
453             tz               = fscal*dz10;
454
455             /* Update vectorial force */
456             fix1            += tx;
457             fiy1            += ty;
458             fiz1            += tz;
459             f[j_coord_offset+DIM*0+XX] -= tx;
460             f[j_coord_offset+DIM*0+YY] -= ty;
461             f[j_coord_offset+DIM*0+ZZ] -= tz;
462
463             /**************************
464              * CALCULATE INTERACTIONS *
465              **************************/
466
467             qq20             = iq2*jq0;
468
469             /* COULOMB ELECTROSTATICS */
470             velec            = qq20*rinv20;
471             felec            = velec*rinvsq20;
472
473             fscal            = felec;
474
475             /* Calculate temporary vectorial force */
476             tx               = fscal*dx20;
477             ty               = fscal*dy20;
478             tz               = fscal*dz20;
479
480             /* Update vectorial force */
481             fix2            += tx;
482             fiy2            += ty;
483             fiz2            += tz;
484             f[j_coord_offset+DIM*0+XX] -= tx;
485             f[j_coord_offset+DIM*0+YY] -= ty;
486             f[j_coord_offset+DIM*0+ZZ] -= tz;
487
488             /**************************
489              * CALCULATE INTERACTIONS *
490              **************************/
491
492             qq30             = iq3*jq0;
493
494             /* COULOMB ELECTROSTATICS */
495             velec            = qq30*rinv30;
496             felec            = velec*rinvsq30;
497
498             fscal            = felec;
499
500             /* Calculate temporary vectorial force */
501             tx               = fscal*dx30;
502             ty               = fscal*dy30;
503             tz               = fscal*dz30;
504
505             /* Update vectorial force */
506             fix3            += tx;
507             fiy3            += ty;
508             fiz3            += tz;
509             f[j_coord_offset+DIM*0+XX] -= tx;
510             f[j_coord_offset+DIM*0+YY] -= ty;
511             f[j_coord_offset+DIM*0+ZZ] -= tz;
512
513             /* Inner loop uses 81 flops */
514         }
515         /* End of innermost loop */
516
517         tx = ty = tz = 0;
518         f[i_coord_offset+DIM*1+XX] += fix1;
519         f[i_coord_offset+DIM*1+YY] += fiy1;
520         f[i_coord_offset+DIM*1+ZZ] += fiz1;
521         tx                         += fix1;
522         ty                         += fiy1;
523         tz                         += fiz1;
524         f[i_coord_offset+DIM*2+XX] += fix2;
525         f[i_coord_offset+DIM*2+YY] += fiy2;
526         f[i_coord_offset+DIM*2+ZZ] += fiz2;
527         tx                         += fix2;
528         ty                         += fiy2;
529         tz                         += fiz2;
530         f[i_coord_offset+DIM*3+XX] += fix3;
531         f[i_coord_offset+DIM*3+YY] += fiy3;
532         f[i_coord_offset+DIM*3+ZZ] += fiz3;
533         tx                         += fix3;
534         ty                         += fiy3;
535         tz                         += fiz3;
536         fshift[i_shift_offset+XX]  += tx;
537         fshift[i_shift_offset+YY]  += ty;
538         fshift[i_shift_offset+ZZ]  += tz;
539
540         /* Increment number of inner iterations */
541         inneriter                  += j_index_end - j_index_start;
542
543         /* Outer loop uses 30 flops */
544     }
545
546     /* Increment number of outer iterations */
547     outeriter        += nri;
548
549     /* Update outer/inner flops */
550
551     inc_nrnb(nrnb,eNR_NBKERNEL_ELEC_W4_F,outeriter*30 + inneriter*81);
552 }