6fb67e57ccb953407ea4cfe4102d8c7ab8164336
[alexxy/gromacs.git] / src / gmxlib / network.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 "gmx_fatal.h"
41 #include "main.h"
42 #include "smalloc.h"
43 #include "network.h"
44 #include "copyrite.h"
45 #include "statutil.h"
46 #include "ctype.h"
47 #include "macros.h"
48
49 #ifdef GMX_LIB_MPI
50 #include <mpi.h>
51 #endif
52
53 #ifdef GMX_THREADS
54 #include "tmpi.h"
55 #endif
56
57 #include "mpelogging.h"
58
59 /* The source code in this file should be thread-safe. 
60       Please keep it that way. */
61
62 bool gmx_mpi_initialized(void)
63 {
64   int n;
65 #ifndef GMX_MPI
66   return 0;
67 #else
68   MPI_Initialized(&n);
69   
70   return n;
71 #endif
72 }
73
74 int gmx_setup(int *argc,char **argv,int *nnodes)
75 {
76 #ifndef GMX_MPI
77   gmx_call("gmx_setup");
78   return 0;
79 #else
80   char   buf[256];
81   int    resultlen;               /* actual length of node name      */
82   int    i,flag;
83   int  mpi_num_nodes;
84   int  mpi_my_rank;
85   char mpi_hostname[MPI_MAX_PROCESSOR_NAME];
86
87   /* Call the MPI routines */
88 #ifdef GMX_LIB_MPI
89 #ifdef GMX_FAHCORE
90   (void) fah_MPI_Init(argc,&argv);
91 #else
92   (void) MPI_Init(argc,&argv);
93 #endif
94 #endif
95   (void) MPI_Comm_size( MPI_COMM_WORLD, &mpi_num_nodes );
96   (void) MPI_Comm_rank( MPI_COMM_WORLD, &mpi_my_rank );
97   (void) MPI_Get_processor_name( mpi_hostname, &resultlen );
98
99
100 #ifdef USE_MPE
101   /* MPE logging routines. Get event IDs from MPE: */
102   /* General events */
103   ev_timestep1               = MPE_Log_get_event_number( );
104   ev_timestep2               = MPE_Log_get_event_number( );
105   ev_force_start             = MPE_Log_get_event_number( );
106   ev_force_finish            = MPE_Log_get_event_number( );
107   ev_do_fnbf_start           = MPE_Log_get_event_number( );
108   ev_do_fnbf_finish          = MPE_Log_get_event_number( );
109   ev_ns_start                = MPE_Log_get_event_number( );
110   ev_ns_finish               = MPE_Log_get_event_number( );
111   ev_calc_bonds_start        = MPE_Log_get_event_number( );
112   ev_calc_bonds_finish       = MPE_Log_get_event_number( );
113   ev_global_stat_start       = MPE_Log_get_event_number( );
114   ev_global_stat_finish      = MPE_Log_get_event_number( );
115   ev_virial_start            = MPE_Log_get_event_number( );
116   ev_virial_finish           = MPE_Log_get_event_number( );
117   
118   /* Shift related events */
119   ev_shift_start             = MPE_Log_get_event_number( );
120   ev_shift_finish            = MPE_Log_get_event_number( );
121   ev_unshift_start           = MPE_Log_get_event_number( );
122   ev_unshift_finish          = MPE_Log_get_event_number( );
123   ev_mk_mshift_start         = MPE_Log_get_event_number( );
124   ev_mk_mshift_finish        = MPE_Log_get_event_number( );
125   
126   /* PME related events */
127   ev_pme_start               = MPE_Log_get_event_number( );
128   ev_pme_finish              = MPE_Log_get_event_number( );
129   ev_spread_on_grid_start    = MPE_Log_get_event_number( );
130   ev_spread_on_grid_finish   = MPE_Log_get_event_number( );
131   ev_sum_qgrid_start         = MPE_Log_get_event_number( );
132   ev_sum_qgrid_finish        = MPE_Log_get_event_number( );
133   ev_gmxfft3d_start          = MPE_Log_get_event_number( );
134   ev_gmxfft3d_finish         = MPE_Log_get_event_number( );
135   ev_solve_pme_start         = MPE_Log_get_event_number( );
136   ev_solve_pme_finish        = MPE_Log_get_event_number( );
137   ev_gather_f_bsplines_start = MPE_Log_get_event_number( );
138   ev_gather_f_bsplines_finish= MPE_Log_get_event_number( );
139   ev_reduce_start            = MPE_Log_get_event_number( );
140   ev_reduce_finish           = MPE_Log_get_event_number( );
141   ev_rscatter_start          = MPE_Log_get_event_number( );
142   ev_rscatter_finish         = MPE_Log_get_event_number( );
143   ev_alltoall_start          = MPE_Log_get_event_number( );
144   ev_alltoall_finish         = MPE_Log_get_event_number( );
145   ev_pmeredist_start         = MPE_Log_get_event_number( );
146   ev_pmeredist_finish        = MPE_Log_get_event_number( );
147   ev_init_pme_start          = MPE_Log_get_event_number( );      
148   ev_init_pme_finish         = MPE_Log_get_event_number( );
149   ev_send_coordinates_start  = MPE_Log_get_event_number( );
150   ev_send_coordinates_finish = MPE_Log_get_event_number( );
151   ev_update_fr_start         = MPE_Log_get_event_number( );
152   ev_update_fr_finish        = MPE_Log_get_event_number( );
153   ev_clear_rvecs_start       = MPE_Log_get_event_number( );
154   ev_clear_rvecs_finish      = MPE_Log_get_event_number( ); 
155   ev_update_start            = MPE_Log_get_event_number( ); 
156   ev_update_finish           = MPE_Log_get_event_number( ); 
157   ev_output_start            = MPE_Log_get_event_number( ); 
158   ev_output_finish           = MPE_Log_get_event_number( ); 
159   ev_sum_lrforces_start      = MPE_Log_get_event_number( ); 
160   ev_sum_lrforces_finish     = MPE_Log_get_event_number( ); 
161   ev_sort_start              = MPE_Log_get_event_number( );
162   ev_sort_finish             = MPE_Log_get_event_number( );
163   ev_sum_qgrid_start         = MPE_Log_get_event_number( );
164   ev_sum_qgrid_finish        = MPE_Log_get_event_number( );
165   
166   /* Essential dynamics related events */
167   ev_edsam_start             = MPE_Log_get_event_number( );
168   ev_edsam_finish            = MPE_Log_get_event_number( );
169   ev_get_coords_start        = MPE_Log_get_event_number( );
170   ev_get_coords_finish       = MPE_Log_get_event_number( );
171   ev_ed_apply_cons_start     = MPE_Log_get_event_number( );
172   ev_ed_apply_cons_finish    = MPE_Log_get_event_number( );
173   ev_fit_to_reference_start  = MPE_Log_get_event_number( );
174   ev_fit_to_reference_finish = MPE_Log_get_event_number( );
175   
176   /* describe events: */
177   if ( mpi_my_rank == 0 ) 
178   {
179     /* General events */
180     MPE_Describe_state(ev_timestep1,               ev_timestep2,                "timestep START",  "magenta" );
181     MPE_Describe_state(ev_force_start,             ev_force_finish,             "force",           "cornflower blue" );
182     MPE_Describe_state(ev_do_fnbf_start,           ev_do_fnbf_finish,           "do_fnbf",         "navy" );
183     MPE_Describe_state(ev_ns_start,                ev_ns_finish,                "neighbor search", "tomato" );
184     MPE_Describe_state(ev_calc_bonds_start,        ev_calc_bonds_finish,        "bonded forces",   "slate blue" );
185     MPE_Describe_state(ev_global_stat_start,       ev_global_stat_finish,       "global stat",     "firebrick3");
186     MPE_Describe_state(ev_update_fr_start,         ev_update_fr_finish,         "update forcerec", "goldenrod");
187     MPE_Describe_state(ev_clear_rvecs_start,       ev_clear_rvecs_finish,       "clear rvecs",     "bisque");
188     MPE_Describe_state(ev_update_start,            ev_update_finish,            "update",          "cornsilk");
189     MPE_Describe_state(ev_output_start,            ev_output_finish,            "output",          "black");
190     MPE_Describe_state(ev_virial_start,            ev_virial_finish,            "calc_virial",     "thistle4");
191     
192     /* PME related events */
193     MPE_Describe_state(ev_pme_start,               ev_pme_finish,               "doing PME",       "grey" );
194     MPE_Describe_state(ev_spread_on_grid_start,    ev_spread_on_grid_finish,    "spread",          "dark orange" );   
195     MPE_Describe_state(ev_sum_qgrid_start,         ev_sum_qgrid_finish,         "sum qgrid",       "slate blue");
196     MPE_Describe_state(ev_gmxfft3d_start,          ev_gmxfft3d_finish,          "fft3d",           "snow2" );   
197     MPE_Describe_state(ev_solve_pme_start,         ev_solve_pme_finish,         "solve PME",       "indian red" );   
198     MPE_Describe_state(ev_gather_f_bsplines_start, ev_gather_f_bsplines_finish, "bsplines",        "light sea green" );   
199     MPE_Describe_state(ev_reduce_start,            ev_reduce_finish,            "reduce",          "cyan1" );
200     MPE_Describe_state(ev_rscatter_start,          ev_rscatter_finish,          "rscatter",        "cyan3" );
201     MPE_Describe_state(ev_alltoall_start,          ev_alltoall_finish,          "alltoall",        "LightCyan4" );
202     MPE_Describe_state(ev_pmeredist_start,         ev_pmeredist_finish,         "pmeredist",       "thistle" );
203     MPE_Describe_state(ev_init_pme_start,          ev_init_pme_finish,          "init PME",        "snow4");
204     MPE_Describe_state(ev_send_coordinates_start,  ev_send_coordinates_finish,  "send_coordinates","blue");
205     MPE_Describe_state(ev_sum_lrforces_start,      ev_sum_lrforces_finish,      "sum_LRforces",    "lime green");
206     MPE_Describe_state(ev_sort_start,              ev_sort_finish,              "sort pme atoms",  "brown");
207     MPE_Describe_state(ev_sum_qgrid_start,         ev_sum_qgrid_finish,         "sum charge grid", "medium orchid");
208     
209     /* Shift related events */
210     MPE_Describe_state(ev_shift_start,             ev_shift_finish,             "shift",           "orange");
211     MPE_Describe_state(ev_unshift_start,           ev_unshift_finish,           "unshift",         "dark orange");    
212     MPE_Describe_state(ev_mk_mshift_start,         ev_mk_mshift_finish,         "mk_mshift",       "maroon");
213         
214     /* Essential dynamics related events */
215     MPE_Describe_state(ev_edsam_start,             ev_edsam_finish,             "EDSAM",           "deep sky blue");
216     MPE_Describe_state(ev_get_coords_start,        ev_get_coords_finish,        "ED get coords",   "steel blue");
217     MPE_Describe_state(ev_ed_apply_cons_start,     ev_ed_apply_cons_finish,     "ED apply constr", "forest green");
218     MPE_Describe_state(ev_fit_to_reference_start,  ev_fit_to_reference_finish,  "ED fit to ref",   "lavender");
219        
220   }
221   MPE_Init_log();
222 #endif
223  
224 #ifdef GMX_LIB_MPI 
225   fprintf(stderr,"NNODES=%d, MYRANK=%d, HOSTNAME=%s\n",
226           mpi_num_nodes,mpi_my_rank,mpi_hostname);
227 #endif
228   
229   *nnodes=mpi_num_nodes;
230   
231   return mpi_my_rank;
232 #endif
233 }
234
235 int  gmx_node_num(void)
236 {
237 #ifndef GMX_MPI
238   return 1;
239 #else
240   int i;
241   (void) MPI_Comm_size(MPI_COMM_WORLD, &i);
242   return i;
243 #endif
244 }
245
246 int gmx_node_rank(void)
247 {
248 #ifndef GMX_MPI
249   return 0;
250 #else
251   int i;
252   (void) MPI_Comm_rank(MPI_COMM_WORLD, &i);
253   return i;
254 #endif
255 }
256
257 void gmx_setup_nodecomm(FILE *fplog,t_commrec *cr)
258 {
259   gmx_nodecomm_t *nc;
260   int  n,rank,resultlen,hostnum,i,j,ng,ni;
261 #ifdef GMX_MPI
262   char mpi_hostname[MPI_MAX_PROCESSOR_NAME],num[MPI_MAX_PROCESSOR_NAME];
263 #endif
264
265   /* Many MPI implementations do not optimize MPI_Allreduce
266    * (and probably also other global communication calls)
267    * for multi-core nodes connected by a network.
268    * We can optimize such communication by using one MPI call
269    * within each node and one between the nodes.
270    * For MVAPICH2 and Intel MPI this reduces the time for
271    * the global_stat communication by 25%
272    * for 2x2-core 3 GHz Woodcrest connected by mixed DDR/SDR Infiniband.
273    * B. Hess, November 2007
274    */
275
276   nc = &cr->nc;
277
278   nc->bUse = FALSE;
279 #ifndef GMX_THREADS
280   if (getenv("GMX_NO_NODECOMM") == NULL) {
281 #ifdef GMX_MPI
282     MPI_Comm_size(cr->mpi_comm_mygroup,&n);
283     MPI_Comm_rank(cr->mpi_comm_mygroup,&rank);
284     MPI_Get_processor_name(mpi_hostname,&resultlen);
285     /* This procedure can only differentiate nodes with host names
286      * that end on unique numbers.
287      */
288     i = 0;
289     j = 0;
290     /* Only parse the host name up to the first dot */
291     while(i < resultlen && mpi_hostname[i] != '.') {
292       if (isdigit(mpi_hostname[i])) {
293         num[j++] = mpi_hostname[i];
294       }
295       i++;
296     }
297     num[j] = '\0';
298     if (j == 0) {
299       hostnum = 0;
300     } else {
301       /* Use only the last 9 decimals, so we don't overflow an int */
302       hostnum = strtol(num + max(0,j-9), NULL, 10); 
303     }
304
305     if (debug) {
306       fprintf(debug,
307               "In gmx_setup_nodecomm: splitting communicator of size %d\n",
308               n);
309       fprintf(debug,"In gmx_setup_nodecomm: hostname '%s', hostnum %d\n",
310               mpi_hostname,hostnum);
311     }
312
313     /* The intra-node communicator, split on node number */
314     MPI_Comm_split(cr->mpi_comm_mygroup,hostnum,rank,&nc->comm_intra);
315     MPI_Comm_rank(nc->comm_intra,&nc->rank_intra);
316     if (debug) {
317       fprintf(debug,"In gmx_setup_nodecomm: node rank %d rank_intra %d\n",
318               rank,nc->rank_intra);
319     }
320     /* The inter-node communicator, split on rank_intra.
321      * We actually only need the one for rank=0,
322      * but it is easier to create them all.
323      */
324     MPI_Comm_split(cr->mpi_comm_mygroup,nc->rank_intra,rank,&nc->comm_inter);
325     /* Check if this really created two step communication */
326     MPI_Comm_size(nc->comm_inter,&ng);
327     MPI_Comm_size(nc->comm_intra,&ni);
328     if (debug) {
329       fprintf(debug,"In gmx_setup_nodecomm: groups %d, my group size %d\n",
330               ng,ni);
331     }
332     if ((ng > 1 && ng < n) || (ni > 1 && ni < n)) {
333       nc->bUse = TRUE;
334       if (fplog)
335         fprintf(fplog,"Using two step summing over %d groups of on average %.1f processes\n\n",ng,(real)n/(real)ng);
336       if (nc->rank_intra > 0)
337         MPI_Comm_free(&nc->comm_inter);
338     } else {
339       /* One group or all processes in a separate group, use normal summing */
340       MPI_Comm_free(&nc->comm_inter);
341       MPI_Comm_free(&nc->comm_intra);
342     }
343 #endif
344   }
345 #endif
346 }
347
348 void gmx_barrier(const t_commrec *cr)
349 {
350 #ifndef GMX_MPI
351   gmx_call("gmx_barrier");
352 #else
353   MPI_Barrier(cr->mpi_comm_mygroup);
354 #endif
355 }
356
357 void gmx_abort(int noderank,int nnodes,int errorno)
358 {
359 #ifndef GMX_MPI
360   gmx_call("gmx_abort");
361 #else
362 #ifdef GMX_THREADS
363   fprintf(stderr,"Halting program %s\n",ShortProgram());
364   thanx(stderr);
365   exit(1);
366 #else
367   if (nnodes > 1)
368   {
369       fprintf(stderr,"Halting parallel program %s on CPU %d out of %d\n",
370               ShortProgram(),noderank,nnodes);
371   }
372   else
373   {
374       fprintf(stderr,"Halting program %s\n",ShortProgram());
375   }
376
377   thanx(stderr);
378   MPI_Abort(MPI_COMM_WORLD,errorno);
379   exit(1);
380 #endif
381 #endif
382 }
383
384 void gmx_bcast(int nbytes,void *b,const t_commrec *cr)
385 {
386 #ifndef GMX_MPI
387   gmx_call("gmx_bast");
388 #else
389   MPI_Bcast(b,nbytes,MPI_BYTE,MASTERRANK(cr),cr->mpi_comm_mygroup);
390 #endif
391 }
392
393 void gmx_bcast_sim(int nbytes,void *b,const t_commrec *cr)
394 {
395 #ifndef GMX_MPI
396   gmx_call("gmx_bast");
397 #else
398   MPI_Bcast(b,nbytes,MPI_BYTE,MASTERRANK(cr),cr->mpi_comm_mysim);
399 #endif
400 }
401
402 void gmx_sumd(int nr,double r[],const t_commrec *cr)
403 {
404 #ifndef GMX_MPI
405     gmx_call("gmx_sumd");
406 #else
407 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREADS)
408     if (cr->nc.bUse) {
409         if (cr->nc.rank_intra == 0)
410         {
411             /* Use two step summing. */
412             MPI_Reduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM,0,
413                        cr->nc.comm_intra);
414             /* Sum the roots of the internal (intra) buffers. */
415             MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM,
416                           cr->nc.comm_inter);
417         }
418         else
419         {
420             /* This is here because of the silly MPI specification
421                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
422             MPI_Reduce(r,NULL,nr,MPI_DOUBLE,MPI_SUM,0,cr->nc.comm_intra);
423         }
424         MPI_Bcast(r,nr,MPI_DOUBLE,0,cr->nc.comm_intra);
425     } 
426     else 
427     {
428         MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM, 
429                       cr->mpi_comm_mygroup);
430     }
431 #else
432     int i;
433
434     if (nr > cr->mpb->dbuf_alloc) {
435         cr->mpb->dbuf_alloc = nr;
436         srenew(cr->mpb->dbuf,cr->mpb->dbuf_alloc);
437     }
438     if (cr->nc.bUse) {
439         /* Use two step summing */
440         MPI_Allreduce(r,cr->mpb->dbuf,nr,MPI_DOUBLE,MPI_SUM,cr->nc.comm_intra);
441         if (cr->nc.rank_intra == 0) {
442             /* Sum with the buffers reversed */
443             MPI_Allreduce(cr->mpb->dbuf,r,nr,MPI_DOUBLE,MPI_SUM, 
444                           cr->nc.comm_inter);
445         }
446         MPI_Bcast(r,nr,MPI_DOUBLE,0,cr->nc.comm_intra);
447     } else {
448         MPI_Allreduce(r,cr->mpb->dbuf,nr,MPI_DOUBLE,MPI_SUM,
449                       cr->mpi_comm_mygroup);
450         for(i=0; i<nr; i++)
451             r[i] = cr->mpb->dbuf[i];
452     }
453 #endif
454 #endif
455 }
456
457 void gmx_sumf(int nr,float r[],const t_commrec *cr)
458 {
459 #ifndef GMX_MPI
460     gmx_call("gmx_sumf");
461 #else
462 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREADS)
463     if (cr->nc.bUse) {
464         /* Use two step summing.  */
465         if (cr->nc.rank_intra == 0)
466         {
467             MPI_Reduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,0,
468                        cr->nc.comm_intra);
469             /* Sum the roots of the internal (intra) buffers */
470             MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,
471                           cr->nc.comm_inter);
472         }
473         else
474         {
475             /* This is here because of the silly MPI specification
476                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
477             MPI_Reduce(r,NULL,nr,MPI_FLOAT,MPI_SUM,0,cr->nc.comm_intra);
478         }
479         MPI_Bcast(r,nr,MPI_FLOAT,0,cr->nc.comm_intra);
480     } 
481     else 
482     {
483         MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,cr->mpi_comm_mygroup);
484     }
485 #else
486     int i;
487
488     if (nr > cr->mpb->fbuf_alloc) {
489         cr->mpb->fbuf_alloc = nr;
490         srenew(cr->mpb->fbuf,cr->mpb->fbuf_alloc);
491     }
492     if (cr->nc.bUse) {
493         /* Use two step summing */
494         MPI_Allreduce(r,cr->mpb->fbuf,nr,MPI_FLOAT,MPI_SUM,cr->nc.comm_intra);
495         if (cr->nc.rank_intra == 0) {
496             /* Sum with the buffers reversed */
497             MPI_Allreduce(cr->mpb->fbuf,r,nr,MPI_FLOAT,MPI_SUM, 
498                           cr->nc.comm_inter);
499         }
500         MPI_Bcast(r,nr,MPI_FLOAT,0,cr->nc.comm_intra);
501     } else {
502         MPI_Allreduce(r,cr->mpb->fbuf,nr,MPI_FLOAT,MPI_SUM,
503                       cr->mpi_comm_mygroup);
504         for(i=0; i<nr; i++)
505             r[i] = cr->mpb->fbuf[i];
506     }
507 #endif
508 #endif
509 }
510
511 void gmx_sumi(int nr,int r[],const t_commrec *cr)
512 {
513 #ifndef GMX_MPI
514     gmx_call("gmx_sumi");
515 #else
516 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREADS)
517     if (cr->nc.bUse) {
518         /* Use two step summing */
519         if (cr->nc.rank_intra == 0) 
520         {
521             MPI_Reduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,0,cr->nc.comm_intra);
522             /* Sum with the buffers reversed */
523             MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,cr->nc.comm_inter);
524         }
525         else
526         {
527             /* This is here because of the silly MPI specification
528                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
529             MPI_Reduce(r,NULL,nr,MPI_INT,MPI_SUM,0,cr->nc.comm_intra);
530         }
531         MPI_Bcast(r,nr,MPI_INT,0,cr->nc.comm_intra);
532     } 
533     else 
534     {
535         MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,cr->mpi_comm_mygroup);
536     }
537 #else
538     int i;
539
540     if (nr > cr->mpb->ibuf_alloc) {
541         cr->mpb->ibuf_alloc = nr;
542         srenew(cr->mpb->ibuf,cr->mpb->ibuf_alloc);
543     }
544     if (cr->nc.bUse) {
545         /* Use two step summing */
546         MPI_Allreduce(r,cr->mpb->ibuf,nr,MPI_INT,MPI_SUM,cr->nc.comm_intra);
547         if (cr->nc.rank_intra == 0) {
548             /* Sum with the buffers reversed */
549             MPI_Allreduce(cr->mpb->ibuf,r,nr,MPI_INT,MPI_SUM,cr->nc.comm_inter);
550         }
551         MPI_Bcast(r,nr,MPI_INT,0,cr->nc.comm_intra);
552     } else {
553         MPI_Allreduce(r,cr->mpb->ibuf,nr,MPI_INT,MPI_SUM,cr->mpi_comm_mygroup);
554         for(i=0; i<nr; i++)
555             r[i] = cr->mpb->ibuf[i];
556     }
557 #endif
558 #endif
559 }
560
561 #ifdef GMX_MPI
562 void gmx_sumd_comm(int nr,double r[],MPI_Comm mpi_comm)
563 {
564 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREADS)
565     MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM,mpi_comm);
566 #else
567     /* this function is only used in code that is not performance critical,
568        (during setup, when comm_rec is not the appropriate communication  
569        structure), so this isn't as bad as it looks. */
570     double *buf;
571     int i;
572
573     snew(buf, nr);
574     MPI_Allreduce(r,buf,nr,MPI_DOUBLE,MPI_SUM,mpi_comm);
575     for(i=0; i<nr; i++)
576         r[i] = buf[i];
577     sfree(buf);
578 #endif
579 }
580 #endif
581
582 #ifdef GMX_MPI
583 void gmx_sumf_comm(int nr,float r[],MPI_Comm mpi_comm)
584 {
585 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREADS)
586     MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,mpi_comm);
587 #else
588     /* this function is only used in code that is not performance critical,
589        (during setup, when comm_rec is not the appropriate communication  
590        structure), so this isn't as bad as it looks. */
591     float *buf;
592     int i;
593
594     snew(buf, nr);
595     MPI_Allreduce(r,buf,nr,MPI_FLOAT,MPI_SUM,mpi_comm);
596     for(i=0; i<nr; i++)
597         r[i] = buf[i];
598     sfree(buf);
599 #endif
600 }
601 #endif
602
603 void gmx_sumd_sim(int nr,double r[],const gmx_multisim_t *ms)
604 {
605 #ifndef GMX_MPI
606   gmx_call("gmx_sumd");
607 #else
608   gmx_sumd_comm(nr,r,ms->mpi_comm_masters);
609 #endif
610 }
611
612 void gmx_sumf_sim(int nr,float r[],const gmx_multisim_t *ms)
613 {
614 #ifndef GMX_MPI
615   gmx_call("gmx_sumf");
616 #else
617   gmx_sumf_comm(nr,r,ms->mpi_comm_masters);
618 #endif
619 }
620
621 void gmx_sumi_sim(int nr,int r[], const gmx_multisim_t *ms)
622 {
623 #ifndef GMX_MPI
624     gmx_call("gmx_sumd");
625 #else
626 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREADS)
627     MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,ms->mpi_comm_masters);
628 #else
629     /* this is thread-unsafe, but it will do for now: */
630     int i;
631
632     if (nr > ms->mpb->ibuf_alloc) {
633         ms->mpb->ibuf_alloc = nr;
634         srenew(ms->mpb->ibuf,ms->mpb->ibuf_alloc);
635     }
636     MPI_Allreduce(r,ms->mpb->ibuf,nr,MPI_INT,MPI_SUM,ms->mpi_comm_masters);
637     for(i=0; i<nr; i++)
638         r[i] = ms->mpb->ibuf[i];
639 #endif
640 #endif
641 }
642
643 void gmx_finalize(void)
644 {
645 #ifndef GMX_MPI
646   gmx_call("gmx_finalize");
647 #else
648   int ret;
649
650   /* just as a check; we don't want to finalize twice */
651   int finalized;
652   MPI_Finalized(&finalized);
653   if (finalized)
654       return;
655
656   /* We sync the processes here to try to avoid problems
657    * with buggy MPI implementations that could cause
658    * unfinished processes to terminate.
659    */
660   MPI_Barrier(MPI_COMM_WORLD);
661
662   /*
663   if (DOMAINDECOMP(cr)) {
664     if (cr->npmenodes > 0 || cr->dd->bCartesian) 
665       MPI_Comm_free(&cr->mpi_comm_mygroup);
666     if (cr->dd->bCartesian)
667       MPI_Comm_free(&cr->mpi_comm_mysim);
668   }
669   */
670
671   /* Apparently certain mpich implementations cause problems
672    * with MPI_Finalize. In that case comment out MPI_Finalize.
673    */
674   if (debug)
675     fprintf(debug,"Will call MPI_Finalize now\n");
676
677   ret = MPI_Finalize();
678   if (debug)
679     fprintf(debug,"Return code from MPI_Finalize = %d\n",ret);
680 #endif
681 }
682