Merge release-4-6 into master
[alexxy/gromacs.git] / src / gromacs / gmxlib / nrnb.c
1 /*
2  * 
3  *                This source code is part of
4  * 
5  *                 G   R   O   M   A   C   S
6  * 
7  *          GROningen MAchine for Chemical Simulations
8  * 
9  *                        VERSION 3.2.0
10  * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12  * Copyright (c) 2001-2004, The GROMACS development team,
13  * check out http://www.gromacs.org for more information.
14
15  * This program is free software; you can redistribute it and/or
16  * modify it under the terms of the GNU General Public License
17  * as published by the Free Software Foundation; either version 2
18  * of the License, or (at your option) any later version.
19  * 
20  * If you want to redistribute modifications, please consider that
21  * scientific software is very special. Version control is crucial -
22  * bugs must be traceable. We will be happy to consider code for
23  * inclusion in the official distribution, but derived work must not
24  * be called official GROMACS. Details are found in the README & COPYING
25  * files - if they are missing, get the official version at www.gromacs.org.
26  * 
27  * To help us fund GROMACS development, we humbly ask that you cite
28  * the papers on the package - you can find them in the top README file.
29  * 
30  * For more info, check our website at http://www.gromacs.org
31  * 
32  * And Hey:
33  * GROningen Mixture of Alchemy and Childrens' Stories
34  */
35 #ifdef HAVE_CONFIG_H
36 #include <config.h>
37 #endif
38
39 #include <string.h>
40 #include "types/commrec.h"
41 #include "sysstuff.h"
42 #include "gmx_fatal.h"
43 #include "names.h"
44 #include "macros.h"
45 #include "nrnb.h"
46 #include "main.h"
47 #include "smalloc.h"
48 #include "copyrite.h"
49
50
51
52
53
54 typedef struct {
55   const char *name;
56   int  flop;
57 } t_nrnb_data;
58
59
60 static const t_nrnb_data nbdata[eNRNB] = {
61     /* These are re-used for different NB kernels, since there are so many.
62      * The actual number of flops is set dynamically.
63      */
64     { "NB VdW [V&F]",                    1 },
65     { "NB VdW [F]",                      1 },
66     { "NB Elec. [V&F]",                  1 },
67     { "NB Elec. [F]",                    1 },
68     { "NB Elec. [W3,V&F]",               1 },
69     { "NB Elec. [W3,F]",                 1 },
70     { "NB Elec. [W3-W3,V&F]",            1 },
71     { "NB Elec. [W3-W3,F]",              1 },
72     { "NB Elec. [W4,V&F]",               1 },
73     { "NB Elec. [W4,F]",                 1 },
74     { "NB Elec. [W4-W4,V&F]",            1 },
75     { "NB Elec. [W4-W4,F]",              1 },
76     { "NB VdW & Elec. [V&F]",            1 },
77     { "NB VdW & Elec. [F]",              1 },
78     { "NB VdW & Elec. [W3,V&F]",         1 },
79     { "NB VdW & Elec. [W3,F]",           1 },
80     { "NB VdW & Elec. [W3-W3,V&F]",      1 },
81     { "NB VdW & Elec. [W3-W3,F]",        1 },
82     { "NB VdW & Elec. [W4,V&F]",         1 },
83     { "NB VdW & Elec. [W4,F]",           1 },
84     { "NB VdW & Elec. [W4-W4,V&F]",      1 },
85     { "NB VdW & Elec. [W4-W4,F]",        1 },
86     
87     { "NB Generic kernel",               1 },
88     { "NB Free energy kernel",           1 },
89     { "NB All-vs-all",                   1 },
90     { "NB All-vs-all, GB",               1 },
91
92     { "Pair Search distance check",      9 }, /* nbnxn pair dist. check */
93     /* nbnxn kernel flops are based on inner-loops without exclusion checks.
94      * Plain Coulomb runs through the RF kernels, except with CUDA.
95      * invsqrt is counted as 6 flops: 1 for _mm_rsqt_ps + 5 for iteration.
96      * The flops are equal for plain-C, x86 SIMD and CUDA, except for:
97      * - plain-C kernel uses one flop more for Coulomb-only (F) than listed
98      * - x86 SIMD LJ geom-comb.rule kernels (fastest) use 2 more flops
99      * - x86 SIMD LJ LB-comb.rule kernels (fast) use 3 (8 for F+E) more flops
100      * - GPU always does exclusions, which requires 2-4 flops, but as invsqrt
101      *   is always counted as 6 flops, this roughly compensates.
102      */
103     { "NxN RF Elec. + VdW [F]",         38 }, /* nbnxn kernel LJ+RF, no ener */
104     { "NxN RF Elec. + VdW [V&F]",       54 },
105     { "NxN QSTab Elec. + VdW [F]",      41 }, /* nbnxn kernel LJ+tab, no en */
106     { "NxN QSTab Elec. + VdW [V&F]",    59 },
107     { "NxN Ewald Elec. + VdW [F]",      66 }, /* nbnxn kernel LJ+Ewald, no en */
108     { "NxN Ewald Elec. + VdW [V&F]",   107 },
109     { "NxN VdW [F]",                    33 }, /* nbnxn kernel LJ, no ener */
110     { "NxN VdW [V&F]",                  43 },
111     { "NxN RF Electrostatics [F]",      31 }, /* nbnxn kernel RF, no ener */
112     { "NxN RF Electrostatics [V&F]",    36 },
113     { "NxN QSTab Elec. [F]",            34 }, /* nbnxn kernel tab, no ener */
114     { "NxN QSTab Elec. [V&F]",          41 },
115     { "NxN Ewald Elec. [F]",            61 }, /* nbnxn kernel Ewald, no ener */
116     { "NxN Ewald Elec. [V&F]",          84 },
117     { "1,4 nonbonded interactions",     90 },
118     { "Born radii (Still)",             47 },
119     { "Born radii (HCT/OBC)",          183 },
120     { "Born force chain rule",          15 },
121     { "All-vs-All Still radii",          1 },
122     { "All-vs-All HCT/OBC radii",        1 },
123     { "All-vs-All Born chain rule",      1 },
124     { "Calc Weights",                   36 },
125     { "Spread Q",                        6 },
126     { "Spread Q Bspline",                2 }, 
127     { "Gather F",                      23  },
128     { "Gather F Bspline",              6   }, 
129     { "3D-FFT",                        8   },
130     { "Convolution",                   4   },
131     { "Solve PME",                     64  },
132     { "NS-Pairs",                      21  },
133     { "Reset In Box",                  3   },
134     { "Shift-X",                       6   },
135     { "CG-CoM",                        3   },
136     { "Sum Forces",                    1   },
137     { "Bonds",                         59  },
138     { "G96Bonds",                      44  },
139     { "FENE Bonds",                    58  },
140     { "Tab. Bonds",                    62  },
141     { "Restraint Potential",           86  },
142     { "Linear Angles",                 57  },
143     { "Angles",                        168 },
144     { "G96Angles",                     150 },
145     { "Quartic Angles",                160 },
146     { "Tab. Angles",                   169 },
147     { "Propers",                       229 },
148     { "Impropers",                     208 },
149     { "RB-Dihedrals",                  247 },
150     { "Four. Dihedrals",               247 },
151     { "Tab. Dihedrals",                227 },
152     { "Dist. Restr.",                  200 },
153     { "Orient. Restr.",                200 },
154     { "Dihedral Restr.",               200 },
155     { "Pos. Restr.",                   50  },
156     { "Flat-bottom posres",            50  },
157     { "Angle Restr.",                  191 },
158     { "Angle Restr. Z",                164 },
159     { "Morse Potent.",                 83  },
160     { "Cubic Bonds",                   54  },
161     { "Walls",                         31  },
162     { "Polarization",                  59  },
163     { "Anharmonic Polarization",       72  },
164     { "Water Pol.",                    62  },
165     { "Thole Pol.",                    296 },
166     { "Virial",                        18  },
167     { "Update",                        31  },
168     { "Ext.ens. Update",               54  },
169     { "Stop-CM",                       10  },
170     { "P-Coupling",                    6   },
171     { "Calc-Ekin",                     27  },
172     { "Lincs",                         60  },
173     { "Lincs-Mat",                     4   },
174     { "Shake",                         30  },
175     { "Constraint-V",                   8  },
176     { "Shake-Init",                    10  },
177     { "Constraint-Vir",                24  },
178     { "Settle",                        323 },
179     { "Virtual Site 2",                23  },
180     { "Virtual Site 3",                37  },
181     { "Virtual Site 3fd",              95  },
182     { "Virtual Site 3fad",             176 },
183     { "Virtual Site 3out",             87  },
184     { "Virtual Site 4fd",              110 }, 
185     { "Virtual Site 4fdn",             254 }, 
186     { "Virtual Site N",                 15 },
187     { "Mixed Generalized Born stuff",   10 } 
188 };
189
190
191 void init_nrnb(t_nrnb *nrnb)
192 {
193   int i;
194
195   for(i=0; (i<eNRNB); i++)
196     nrnb->n[i]=0.0;
197 }
198
199 void cp_nrnb(t_nrnb *dest, t_nrnb *src)
200 {
201   int i;
202
203   for(i=0; (i<eNRNB); i++)
204     dest->n[i]=src->n[i];
205 }
206
207 void add_nrnb(t_nrnb *dest, t_nrnb *s1, t_nrnb *s2)
208 {
209   int i;
210
211   for(i=0; (i<eNRNB); i++)
212     dest->n[i]=s1->n[i]+s2->n[i];
213 }
214
215 void print_nrnb(FILE *out, t_nrnb *nrnb)
216 {
217   int i;
218
219   for(i=0; (i<eNRNB); i++)
220     if (nrnb->n[i] > 0)
221       fprintf(out," %-26s %10.0f.\n",nbdata[i].name,nrnb->n[i]);
222 }
223
224 void _inc_nrnb(t_nrnb *nrnb,int enr,int inc,char *file,int line)
225 {
226   nrnb->n[enr]+=inc;
227 #ifdef DEBUG_NRNB
228   printf("nrnb %15s(%2d) incremented with %8d from file %s line %d\n",
229           nbdata[enr].name,enr,inc,file,line);
230 #endif
231 }
232
233 void print_flop(FILE *out,t_nrnb *nrnb,double *nbfs,double *mflop)
234 {
235   int    i;
236   double mni,frac,tfrac,tflop;
237   const char   *myline = "-----------------------------------------------------------------------------";
238   
239   *nbfs = 0.0;
240   for(i=0; (i<eNR_NBKERNEL_ALLVSALLGB); i++) {
241     if (strstr(nbdata[i].name,"W3-W3") != NULL)
242       *nbfs += 9e-6*nrnb->n[i];
243     else if (strstr(nbdata[i].name,"W3") != NULL)
244       *nbfs += 3e-6*nrnb->n[i];
245     else if (strstr(nbdata[i].name,"W4-W4") != NULL)
246       *nbfs += 10e-6*nrnb->n[i];
247     else if (strstr(nbdata[i].name,"W4") != NULL)
248       *nbfs += 4e-6*nrnb->n[i];
249     else
250       *nbfs += 1e-6*nrnb->n[i];
251   }
252   tflop=0;
253   for(i=0; (i<eNRNB); i++) 
254     tflop+=1e-6*nrnb->n[i]*nbdata[i].flop;
255   
256   if (tflop == 0) {
257     fprintf(out,"No MEGA Flopsen this time\n");
258     return;
259   }
260   if (out) {
261     fprintf(out,"\n\tM E G A - F L O P S   A C C O U N T I N G\n\n");
262   }
263
264   if (out)
265   {
266       fprintf(out," NB=Group-cutoff nonbonded kernels    NxN=N-by-N cluster Verlet kernels\n");
267       fprintf(out," RF=Reaction-Field  VdW=Van der Waals  QSTab=quadratic-spline table\n");
268       fprintf(out," W3=SPC/TIP3p  W4=TIP4p (single or pairs)\n");
269       fprintf(out," V&F=Potential and force  V=Potential only  F=Force only\n\n");
270
271       fprintf(out," %-32s %16s %15s  %7s\n",
272               "Computing:","M-Number","M-Flops","% Flops");
273       fprintf(out,"%s\n",myline);
274   }
275   *mflop=0.0;
276   tfrac=0.0;
277   for(i=0; (i<eNRNB); i++) {
278     mni     = 1e-6*nrnb->n[i];
279     *mflop += mni*nbdata[i].flop;
280     frac    = 100.0*mni*nbdata[i].flop/tflop;
281     tfrac  += frac;
282     if (out && mni != 0)
283       fprintf(out," %-32s %16.6f %15.3f  %6.1f\n",
284               nbdata[i].name,mni,mni*nbdata[i].flop,frac);
285   }
286   if (out) {
287     fprintf(out,"%s\n",myline);
288     fprintf(out," %-32s %16s %15.3f  %6.1f\n",
289             "Total","",*mflop,tfrac);
290     fprintf(out,"%s\n\n",myline);
291   }
292 }
293
294 void print_perf(FILE *out,double nodetime,double realtime,int nprocs,
295                 gmx_large_int_t nsteps,real delta_t,
296                 double nbfs,double mflop,
297                 int omp_nth_pp)
298 {
299   real runtime;
300
301   fprintf(out,"\n");
302
303   if (realtime > 0) 
304   {
305     fprintf(out,"%12s %12s %12s %10s\n","","Core t (s)","Wall t (s)","(%)");
306     fprintf(out,"%12s %12.3f %12.3f %10.1f\n","Time:",
307             nodetime, realtime, 100.0*nodetime/realtime);
308     /* only print day-hour-sec format if realtime is more than 30 min */
309     if (realtime > 30*60)
310     {
311       fprintf(out,"%12s %12s","","");
312       pr_difftime(out,realtime);
313     }
314     if (delta_t > 0) 
315     {
316       mflop = mflop/realtime;
317       runtime = nsteps*delta_t;
318
319       if (getenv("GMX_DETAILED_PERF_STATS") == NULL)
320       {
321           fprintf(out,"%12s %12s %12s\n",
322                   "","(ns/day)","(hour/ns)");
323           fprintf(out,"%12s %12.3f %12.3f\n","Performance:",
324                   runtime*24*3.6/realtime,1000*realtime/(3600*runtime));
325       }
326       else
327       {
328         fprintf(out,"%12s %12s %12s %12s %12s\n",
329                 "","(Mnbf/s)",(mflop > 1000) ? "(GFlops)" : "(MFlops)",
330                 "(ns/day)","(hour/ns)");
331         fprintf(out,"%12s %12.3f %12.3f %12.3f %12.3f\n","Performance:",
332                 nbfs/realtime,(mflop > 1000) ? (mflop/1000) : mflop,
333                 runtime*24*3.6/realtime,1000*realtime/(3600*runtime));
334       }
335     } 
336     else 
337     {
338       if (getenv("GMX_DETAILED_PERF_STATS") == NULL)
339       {
340           fprintf(out,"%12s %14s\n",
341                   "","(steps/hour)");
342           fprintf(out,"%12s %14.1f\n","Performance:",
343                   nsteps*3600.0/realtime);
344       }
345       else
346       {
347           fprintf(out,"%12s %12s %12s %14s\n",
348                   "","(Mnbf/s)",(mflop > 1000) ? "(GFlops)" : "(MFlops)",
349                   "(steps/hour)");
350           fprintf(out,"%12s %12.3f %12.3f %14.1f\n","Performance:",
351               nbfs/realtime,(mflop > 1000) ? (mflop/1000) : mflop,
352               nsteps*3600.0/realtime);
353       }
354     }
355   }
356 }
357
358 int cost_nrnb(int enr)
359 {
360   return nbdata[enr].flop;
361 }
362
363 const char *nrnb_str(int enr)
364 {
365   return nbdata[enr].name;
366 }
367
368 static const int    force_index[]={ 
369   eNR_BONDS,  eNR_ANGLES,  eNR_PROPER, eNR_IMPROPER, 
370   eNR_RB,     eNR_DISRES,  eNR_ORIRES, eNR_POSRES,
371   eNR_FBPOSRES, eNR_NS,
372 };
373 #define NFORCE_INDEX asize(force_index)
374
375 static const int    constr_index[]={ 
376   eNR_SHAKE,     eNR_SHAKE_RIJ, eNR_SETTLE,       eNR_UPDATE,       eNR_PCOUPL,
377   eNR_CONSTR_VIR,eNR_CONSTR_V
378 };
379 #define NCONSTR_INDEX asize(constr_index)
380
381 static double pr_av(FILE *log,t_commrec *cr,
382                     double fav,double ftot[],const char *title)
383 {
384   int    i,perc;
385   double dperc,unb;
386   
387   unb=0;
388   if (fav > 0) {
389     fav /= cr->nnodes - cr->npmenodes;
390     fprintf(log,"\n %-26s",title);
391     for(i=0; (i<cr->nnodes); i++) {
392         dperc=(100.0*ftot[i])/fav;
393         unb=max(unb,dperc);
394         perc=dperc;
395         fprintf(log,"%3d ",perc);
396     }
397     if (unb > 0) {
398       perc=10000.0/unb;
399       fprintf(log,"%6d%%\n\n",perc);
400     }
401     else
402       fprintf(log,"\n\n");
403   }
404   return unb;
405 }
406
407 void pr_load(FILE *log,t_commrec *cr,t_nrnb nrnb[])
408 {
409   int    i,j,perc;
410   double dperc,unb,uf,us;
411   double *ftot,fav;
412   double *stot,sav;
413   t_nrnb *av;
414
415   snew(av,1);
416   snew(ftot,cr->nnodes);
417   snew(stot,cr->nnodes);
418   init_nrnb(av);
419   for(i=0; (i<cr->nnodes); i++) {
420       add_nrnb(av,av,&(nrnb[i]));
421       /* Cost due to forces */
422       for(j=0; (j<eNR_NBKERNEL_ALLVSALLGB); j++)
423         ftot[i]+=nrnb[i].n[j]*cost_nrnb(j);
424       for(j=0; (j<NFORCE_INDEX); j++) 
425         ftot[i]+=nrnb[i].n[force_index[j]]*cost_nrnb(force_index[j]);
426       /* Due to shake */
427       for(j=0; (j<NCONSTR_INDEX); j++) {
428         stot[i]+=nrnb[i].n[constr_index[j]]*cost_nrnb(constr_index[j]);
429       }
430   }   
431   for(j=0; (j<eNRNB); j++)
432     av->n[j]=av->n[j]/(double)(cr->nnodes - cr->npmenodes);
433     
434     fprintf(log,"\nDetailed load balancing info in percentage of average\n");
435   
436   fprintf(log," Type                 NODE:");
437   for(i=0; (i<cr->nnodes); i++)
438       fprintf(log,"%3d ",i);
439   fprintf(log,"Scaling\n");
440   fprintf(log,"---------------------------");
441   for(i=0; (i<cr->nnodes); i++)
442       fprintf(log,"----");
443   fprintf(log,"-------\n");
444   
445   for(j=0; (j<eNRNB); j++) {
446     unb=100.0;
447     if (av->n[j] > 0) {
448       fprintf(log," %-26s",nrnb_str(j));
449       for(i=0; (i<cr->nnodes); i++) {
450           dperc=(100.0*nrnb[i].n[j])/av->n[j];
451           unb=max(unb,dperc);
452           perc=dperc;
453           fprintf(log,"%3d ",perc);
454       }
455       if (unb > 0) {
456         perc=10000.0/unb;
457         fprintf(log,"%6d%%\n",perc);
458       }
459       else
460         fprintf(log,"\n");
461     }   
462   }
463   fav=sav=0;
464   for(i=0; (i<cr->nnodes); i++) {
465     fav+=ftot[i];
466     sav+=stot[i];
467   }
468   uf=pr_av(log,cr,fav,ftot,"Total Force");
469   us=pr_av(log,cr,sav,stot,"Total Constr.");
470   
471   unb=(uf*fav+us*sav)/(fav+sav);
472   if (unb > 0) {
473     unb=10000.0/unb;
474     fprintf(log,"\nTotal Scaling: %.0f%% of max performance\n\n",unb);
475   }
476 }
477