src/gmxlib/network.c

   1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
   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_THREAD_MPI
  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 gmx_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
 258 int gmx_hostname_num()
 259 {
 260 #ifndef GMX_MPI
 261   return 0;
 262 #else
 263 #ifdef GMX_THREAD_MPI
 264   /* thread-MPI currently puts the thread number in the process name,
 265    * we might want to change this, as this is inconsistent with what
 266    * most MPI implementations would do when running on a single node.
 267    */
 268   return 0;
 269 #else
 270   int  resultlen,hostnum,i,j;
 271   char mpi_hostname[MPI_MAX_PROCESSOR_NAME],hostnum_str[MPI_MAX_PROCESSOR_NAME];
 272
 273   MPI_Get_processor_name(mpi_hostname,&resultlen);
 274   /* This procedure can only differentiate nodes with host names
 275    * that end on unique numbers.
 276    */
 277   i = 0;
 278   j = 0;
 279   /* Only parse the host name up to the first dot */
 280   while(i < resultlen && mpi_hostname[i] != '.') {
 281     if (isdigit(mpi_hostname[i])) {
 282       hostnum_str[j++] = mpi_hostname[i];
 283     }
 284     i++;
 285   }
 286   hostnum_str[j] = '\0';
 287   if (j == 0) {
 288     hostnum = 0;
 289   } else {
 290     /* Use only the last 9 decimals, so we don't overflow an int */
 291     hostnum = strtol(hostnum_str + max(0,j-9), NULL, 10);
 292   }
 293
 294   if (debug) {
 295     fprintf(debug,"In gmx_setup_nodecomm: hostname '%s', hostnum %d\n",
 296         mpi_hostname,hostnum);
 297   }
 298   return hostnum;
 299 #endif
 300 #endif
 301 }
 302
 303 void gmx_setup_nodecomm(FILE *fplog,t_commrec *cr)
 304 {
 305     gmx_nodecomm_t *nc;
 306     int  n,rank,hostnum,ng,ni;
 307
 308     /* Many MPI implementations do not optimize MPI_Allreduce
 309      * (and probably also other global communication calls)
 310      * for multi-core nodes connected by a network.
 311      * We can optimize such communication by using one MPI call
 312      * within each node and one between the nodes.
 313      * For MVAPICH2 and Intel MPI this reduces the time for
 314      * the global_stat communication by 25%
 315      * for 2x2-core 3 GHz Woodcrest connected by mixed DDR/SDR Infiniband.
 316      * B. Hess, November 2007
 317      */
 318
 319     nc = &cr->nc;
 320
 321     nc->bUse = FALSE;
 322 #ifndef GMX_THREAD_MPI
 323 #ifdef GMX_MPI
 324     MPI_Comm_size(cr->mpi_comm_mygroup,&n);
 325     MPI_Comm_rank(cr->mpi_comm_mygroup,&rank);
 326
 327     hostnum = gmx_hostname_num();
 328
 329     if (debug)
 330     {
 331         fprintf(debug,"In gmx_setup_nodecomm: splitting communicator of size %d\n",n);
 332     }
 333
 334
 335     /* The intra-node communicator, split on node number */
 336     MPI_Comm_split(cr->mpi_comm_mygroup,hostnum,rank,&nc->comm_intra);
 337     MPI_Comm_rank(nc->comm_intra,&nc->rank_intra);
 338     if (debug)
 339     {
 340         fprintf(debug,"In gmx_setup_nodecomm: node rank %d rank_intra %d\n",
 341                 rank,nc->rank_intra);
 342     }
 343     /* The inter-node communicator, split on rank_intra.
 344      * We actually only need the one for rank=0,
 345      * but it is easier to create them all.
 346      */
 347     MPI_Comm_split(cr->mpi_comm_mygroup,nc->rank_intra,rank,&nc->comm_inter);
 348     /* Check if this really created two step communication */
 349     MPI_Comm_size(nc->comm_inter,&ng);
 350     MPI_Comm_size(nc->comm_intra,&ni);
 351     if (debug)
 352     {
 353         fprintf(debug,"In gmx_setup_nodecomm: groups %d, my group size %d\n",
 354                 ng,ni);
 355     }
 356
 357     if (getenv("GMX_NO_NODECOMM") == NULL &&
 358         ((ng > 1 && ng < n) || (ni > 1 && ni < n)))
 359     {
 360         nc->bUse = TRUE;
 361         if (fplog)
 362         {
 363             fprintf(fplog,"Using two step summing over %d groups of on average %.1f processes\n\n",
 364                     ng,(real)n/(real)ng);
 365         }
 366         if (nc->rank_intra > 0)
 367         {
 368             MPI_Comm_free(&nc->comm_inter);
 369         }
 370     }
 371     else
 372     {
 373         /* One group or all processes in a separate group, use normal summing */
 374         MPI_Comm_free(&nc->comm_inter);
 375         MPI_Comm_free(&nc->comm_intra);
 376         if (debug)
 377         {
 378             fprintf(debug,"In gmx_setup_nodecomm: not unsing separate inter- and intra-node communicators.\n");
 379         }
 380     }
 381 #endif
 382 #else
 383     /* tMPI runs only on a single node so just use the nodeid */
 384     nc->rank_intra = cr->nodeid;
 385 #endif
 386 }
 387
 388 void gmx_init_intra_counters(t_commrec *cr)
 389 {
 390     /* counters for PP+PME and PP-only processes on my node */
 391     int nnodes, nnodes_pp, id_mynode=-1, id_mynode_group=-1, nproc_mynode, nproc_mynode_pp;
 392 #if defined GMX_MPI && !defined GMX_THREAD_MPI
 393     int i, mynum, *num, *num_s, *num_pp, *num_pp_s;
 394 #endif
 395
 396     nnodes    = cr->nnodes;
 397     nnodes_pp = nnodes - cr->npmenodes;
 398
 399 #if defined GMX_MPI && !defined GMX_THREAD_MPI
 400     /* We have MPI and can expect to have different compute nodes */
 401     mynum = gmx_hostname_num();
 402
 403     /* We can't rely on MPI_IN_PLACE, so we need send and receive buffers */
 404     snew(num,   nnodes);
 405     snew(num_s, nnodes);
 406     snew(num_pp,   nnodes_pp);
 407     snew(num_pp_s, nnodes_pp);
 408
 409     num_s[cr->sim_nodeid] = mynum;
 410     if (cr->duty & DUTY_PP)
 411     {
 412         num_pp_s[cr->nodeid] = mynum;
 413     }
 414
 415     MPI_Allreduce(num_s, num, nnodes, MPI_INT, MPI_SUM, cr->mpi_comm_mysim);
 416     MPI_Allreduce(num_pp_s, num_pp, nnodes_pp, MPI_INT, MPI_SUM, cr->mpi_comm_mygroup);
 417
 418     id_mynode       = 0;
 419     id_mynode_group = 0;
 420     nproc_mynode    = 0;
 421     nproc_mynode_pp = 0;
 422     for(i=0; i<nnodes; i++)
 423     {
 424         if (num[i] == mynum)
 425         {
 426             nproc_mynode++;
 427             if (i < cr->sim_nodeid)
 428             {
 429                 id_mynode++;
 430             }
 431             if (i < cr->nodeid)
 432             {
 433                 id_mynode_group++;
 434             }
 435         }
 436     }
 437     for(i=0; i<nnodes_pp; i++)
 438     {
 439         if (num_pp[i] == mynum)
 440         {
 441             nproc_mynode_pp++;
 442         }
 443     }
 444     sfree(num);
 445     sfree(num_s);
 446     sfree(num_pp);
 447     sfree(num_pp_s);
 448 #else
 449     /* Serial or thread-MPI code, we are running within a node */
 450     id_mynode       = cr->sim_nodeid;
 451     id_mynode_group = cr->nodeid;
 452     nproc_mynode    = cr->nnodes;
 453     nproc_mynode_pp = cr->nnodes - cr->npmenodes;
 454 #endif
 455
 456     if (debug)
 457     {
 458         char sbuf[STRLEN];
 459         if (cr->duty & DUTY_PP && cr->duty & DUTY_PME)
 460         {
 461             sprintf(sbuf, "PP+PME");
 462         }
 463         else
 464         {
 465             sprintf(sbuf, "%s", cr->duty & DUTY_PP ? "PP" : "PME");
 466         }
 467         fprintf(debug, "On %3s node %d: nodeid_intra=%d, nodeid_group_intra=%d, "
 468                 "nnodes_intra=%d, nnodes_pp_intra=%d\n", sbuf, cr->sim_nodeid,
 469                 id_mynode, id_mynode_group, nproc_mynode, nproc_mynode_pp);
 470     }
 471
 472     cr->nodeid_intra        = id_mynode;
 473     cr->nodeid_group_intra  = id_mynode_group;
 474     cr->nnodes_intra        = nproc_mynode;
 475     cr->nnodes_pp_intra     = nproc_mynode_pp;
 476 }
 477
 478
 479 void gmx_barrier(const t_commrec *cr)
 480 {
 481 #ifndef GMX_MPI
 482   gmx_call("gmx_barrier");
 483 #else
 484   MPI_Barrier(cr->mpi_comm_mygroup);
 485 #endif
 486 }
 487
 488 void gmx_abort(int noderank,int nnodes,int errorno)
 489 {
 490 #ifndef GMX_MPI
 491   gmx_call("gmx_abort");
 492 #else
 493 #ifdef GMX_THREAD_MPI
 494   fprintf(stderr,"Halting program %s\n",ShortProgram());
 495   thanx(stderr);
 496   exit(1);
 497 #else
 498   if (nnodes > 1)
 499   {
 500       fprintf(stderr,"Halting parallel program %s on CPU %d out of %d\n",
 501               ShortProgram(),noderank,nnodes);
 502   }
 503   else
 504   {
 505       fprintf(stderr,"Halting program %s\n",ShortProgram());
 506   }
 507
 508   thanx(stderr);
 509   MPI_Abort(MPI_COMM_WORLD,errorno);
 510   exit(1);
 511 #endif
 512 #endif
 513 }
 514
 515 void gmx_bcast(int nbytes,void *b,const t_commrec *cr)
 516 {
 517 #ifndef GMX_MPI
 518   gmx_call("gmx_bast");
 519 #else
 520   MPI_Bcast(b,nbytes,MPI_BYTE,MASTERRANK(cr),cr->mpi_comm_mygroup);
 521 #endif
 522 }
 523
 524 void gmx_bcast_sim(int nbytes,void *b,const t_commrec *cr)
 525 {
 526 #ifndef GMX_MPI
 527   gmx_call("gmx_bast");
 528 #else
 529   MPI_Bcast(b,nbytes,MPI_BYTE,MASTERRANK(cr),cr->mpi_comm_mysim);
 530 #endif
 531 }
 532
 533 void gmx_sumd(int nr,double r[],const t_commrec *cr)
 534 {
 535 #ifndef GMX_MPI
 536     gmx_call("gmx_sumd");
 537 #else
 538 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 539     if (cr->nc.bUse) {
 540         if (cr->nc.rank_intra == 0)
 541         {
 542             /* Use two step summing. */
 543             MPI_Reduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM,0,
 544                        cr->nc.comm_intra);
 545             /* Sum the roots of the internal (intra) buffers. */
 546             MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM,
 547                           cr->nc.comm_inter);
 548         }
 549         else
 550         {
 551             /* This is here because of the silly MPI specification
 552                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 553             MPI_Reduce(r,NULL,nr,MPI_DOUBLE,MPI_SUM,0,cr->nc.comm_intra);
 554         }
 555         MPI_Bcast(r,nr,MPI_DOUBLE,0,cr->nc.comm_intra);
 556     }
 557     else
 558     {
 559         MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM,
 560                       cr->mpi_comm_mygroup);
 561     }
 562 #else
 563     int i;
 564
 565     if (nr > cr->mpb->dbuf_alloc) {
 566         cr->mpb->dbuf_alloc = nr;
 567         srenew(cr->mpb->dbuf,cr->mpb->dbuf_alloc);
 568     }
 569     if (cr->nc.bUse) {
 570         /* Use two step summing */
 571         MPI_Allreduce(r,cr->mpb->dbuf,nr,MPI_DOUBLE,MPI_SUM,cr->nc.comm_intra);
 572         if (cr->nc.rank_intra == 0) {
 573             /* Sum with the buffers reversed */
 574             MPI_Allreduce(cr->mpb->dbuf,r,nr,MPI_DOUBLE,MPI_SUM,
 575                           cr->nc.comm_inter);
 576         }
 577         MPI_Bcast(r,nr,MPI_DOUBLE,0,cr->nc.comm_intra);
 578     } else {
 579         MPI_Allreduce(r,cr->mpb->dbuf,nr,MPI_DOUBLE,MPI_SUM,
 580                       cr->mpi_comm_mygroup);
 581         for(i=0; i<nr; i++)
 582             r[i] = cr->mpb->dbuf[i];
 583     }
 584 #endif
 585 #endif
 586 }
 587
 588 void gmx_sumf(int nr,float r[],const t_commrec *cr)
 589 {
 590 #ifndef GMX_MPI
 591     gmx_call("gmx_sumf");
 592 #else
 593 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 594     if (cr->nc.bUse) {
 595         /* Use two step summing.  */
 596         if (cr->nc.rank_intra == 0)
 597         {
 598             MPI_Reduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,0,
 599                        cr->nc.comm_intra);
 600             /* Sum the roots of the internal (intra) buffers */
 601             MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,
 602                           cr->nc.comm_inter);
 603         }
 604         else
 605         {
 606             /* This is here because of the silly MPI specification
 607                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 608             MPI_Reduce(r,NULL,nr,MPI_FLOAT,MPI_SUM,0,cr->nc.comm_intra);
 609         }
 610         MPI_Bcast(r,nr,MPI_FLOAT,0,cr->nc.comm_intra);
 611     }
 612     else
 613     {
 614         MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,cr->mpi_comm_mygroup);
 615     }
 616 #else
 617     int i;
 618
 619     if (nr > cr->mpb->fbuf_alloc) {
 620         cr->mpb->fbuf_alloc = nr;
 621         srenew(cr->mpb->fbuf,cr->mpb->fbuf_alloc);
 622     }
 623     if (cr->nc.bUse) {
 624         /* Use two step summing */
 625         MPI_Allreduce(r,cr->mpb->fbuf,nr,MPI_FLOAT,MPI_SUM,cr->nc.comm_intra);
 626         if (cr->nc.rank_intra == 0) {
 627             /* Sum with the buffers reversed */
 628             MPI_Allreduce(cr->mpb->fbuf,r,nr,MPI_FLOAT,MPI_SUM,
 629                           cr->nc.comm_inter);
 630         }
 631         MPI_Bcast(r,nr,MPI_FLOAT,0,cr->nc.comm_intra);
 632     } else {
 633         MPI_Allreduce(r,cr->mpb->fbuf,nr,MPI_FLOAT,MPI_SUM,
 634                       cr->mpi_comm_mygroup);
 635         for(i=0; i<nr; i++)
 636             r[i] = cr->mpb->fbuf[i];
 637     }
 638 #endif
 639 #endif
 640 }
 641
 642 void gmx_sumi(int nr,int r[],const t_commrec *cr)
 643 {
 644 #ifndef GMX_MPI
 645     gmx_call("gmx_sumi");
 646 #else
 647 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 648     if (cr->nc.bUse) {
 649         /* Use two step summing */
 650         if (cr->nc.rank_intra == 0)
 651         {
 652             MPI_Reduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,0,cr->nc.comm_intra);
 653             /* Sum with the buffers reversed */
 654             MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,cr->nc.comm_inter);
 655         }
 656         else
 657         {
 658             /* This is here because of the silly MPI specification
 659                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 660             MPI_Reduce(r,NULL,nr,MPI_INT,MPI_SUM,0,cr->nc.comm_intra);
 661         }
 662         MPI_Bcast(r,nr,MPI_INT,0,cr->nc.comm_intra);
 663     }
 664     else
 665     {
 666         MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,cr->mpi_comm_mygroup);
 667     }
 668 #else
 669     int i;
 670
 671     if (nr > cr->mpb->ibuf_alloc) {
 672         cr->mpb->ibuf_alloc = nr;
 673         srenew(cr->mpb->ibuf,cr->mpb->ibuf_alloc);
 674     }
 675     if (cr->nc.bUse) {
 676         /* Use two step summing */
 677         MPI_Allreduce(r,cr->mpb->ibuf,nr,MPI_INT,MPI_SUM,cr->nc.comm_intra);
 678         if (cr->nc.rank_intra == 0) {
 679             /* Sum with the buffers reversed */
 680             MPI_Allreduce(cr->mpb->ibuf,r,nr,MPI_INT,MPI_SUM,cr->nc.comm_inter);
 681         }
 682         MPI_Bcast(r,nr,MPI_INT,0,cr->nc.comm_intra);
 683     } else {
 684         MPI_Allreduce(r,cr->mpb->ibuf,nr,MPI_INT,MPI_SUM,cr->mpi_comm_mygroup);
 685         for(i=0; i<nr; i++)
 686             r[i] = cr->mpb->ibuf[i];
 687     }
 688 #endif
 689 #endif
 690 }
 691
 692 void gmx_sumli(int nr,gmx_large_int_t r[],const t_commrec *cr)
 693 {
 694 #ifndef GMX_MPI
 695     gmx_call("gmx_sumli");
 696 #else
 697 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 698     if (cr->nc.bUse) {
 699         /* Use two step summing */
 700         if (cr->nc.rank_intra == 0)
 701         {
 702             MPI_Reduce(MPI_IN_PLACE,r,nr,GMX_MPI_LARGE_INT,MPI_SUM,0,
 703                        cr->nc.comm_intra);
 704             /* Sum with the buffers reversed */
 705             MPI_Allreduce(MPI_IN_PLACE,r,nr,GMX_MPI_LARGE_INT,MPI_SUM,
 706                           cr->nc.comm_inter);
 707         }
 708         else
 709         {
 710             /* This is here because of the silly MPI specification
 711                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 712             MPI_Reduce(r,NULL,nr,GMX_MPI_LARGE_INT,MPI_SUM,0,cr->nc.comm_intra);
 713         }
 714         MPI_Bcast(r,nr,GMX_MPI_LARGE_INT,0,cr->nc.comm_intra);
 715     }
 716     else
 717     {
 718         MPI_Allreduce(MPI_IN_PLACE,r,nr,GMX_MPI_LARGE_INT,MPI_SUM,cr->mpi_comm_mygroup);
 719     }
 720 #else
 721     int i;
 722
 723     if (nr > cr->mpb->libuf_alloc) {
 724         cr->mpb->libuf_alloc = nr;
 725         srenew(cr->mpb->libuf,cr->mpb->libuf_alloc);
 726     }
 727     if (cr->nc.bUse) {
 728         /* Use two step summing */
 729         MPI_Allreduce(r,cr->mpb->libuf,nr,GMX_MPI_LARGE_INT,MPI_SUM,
 730                       cr->nc.comm_intra);
 731         if (cr->nc.rank_intra == 0) {
 732             /* Sum with the buffers reversed */
 733             MPI_Allreduce(cr->mpb->libuf,r,nr,GMX_MPI_LARGE_INT,MPI_SUM,
 734                           cr->nc.comm_inter);
 735         }
 736         MPI_Bcast(r,nr,GMX_MPI_LARGE_INT,0,cr->nc.comm_intra);
 737     } else {
 738         MPI_Allreduce(r,cr->mpb->libuf,nr,GMX_MPI_LARGE_INT,MPI_SUM,
 739                       cr->mpi_comm_mygroup);
 740         for(i=0; i<nr; i++)
 741             r[i] = cr->mpb->libuf[i];
 742     }
 743 #endif
 744 #endif
 745 }
 746
 747
 748
 749 #ifdef GMX_MPI
 750 void gmx_sumd_comm(int nr,double r[],MPI_Comm mpi_comm)
 751 {
 752 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 753     MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_DOUBLE,MPI_SUM,mpi_comm);
 754 #else
 755     /* this function is only used in code that is not performance critical,
 756        (during setup, when comm_rec is not the appropriate communication
 757        structure), so this isn't as bad as it looks. */
 758     double *buf;
 759     int i;
 760
 761     snew(buf, nr);
 762     MPI_Allreduce(r,buf,nr,MPI_DOUBLE,MPI_SUM,mpi_comm);
 763     for(i=0; i<nr; i++)
 764         r[i] = buf[i];
 765     sfree(buf);
 766 #endif
 767 }
 768 #endif
 769
 770 #ifdef GMX_MPI
 771 void gmx_sumf_comm(int nr,float r[],MPI_Comm mpi_comm)
 772 {
 773 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 774     MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_FLOAT,MPI_SUM,mpi_comm);
 775 #else
 776     /* this function is only used in code that is not performance critical,
 777        (during setup, when comm_rec is not the appropriate communication
 778        structure), so this isn't as bad as it looks. */
 779     float *buf;
 780     int i;
 781
 782     snew(buf, nr);
 783     MPI_Allreduce(r,buf,nr,MPI_FLOAT,MPI_SUM,mpi_comm);
 784     for(i=0; i<nr; i++)
 785         r[i] = buf[i];
 786     sfree(buf);
 787 #endif
 788 }
 789 #endif
 790
 791 void gmx_sumd_sim(int nr,double r[],const gmx_multisim_t *ms)
 792 {
 793 #ifndef GMX_MPI
 794   gmx_call("gmx_sumd_sim");
 795 #else
 796   gmx_sumd_comm(nr,r,ms->mpi_comm_masters);
 797 #endif
 798 }
 799
 800 void gmx_sumf_sim(int nr,float r[],const gmx_multisim_t *ms)
 801 {
 802 #ifndef GMX_MPI
 803   gmx_call("gmx_sumf_sim");
 804 #else
 805   gmx_sumf_comm(nr,r,ms->mpi_comm_masters);
 806 #endif
 807 }
 808
 809 void gmx_sumi_sim(int nr,int r[], const gmx_multisim_t *ms)
 810 {
 811 #ifndef GMX_MPI
 812     gmx_call("gmx_sumi_sim");
 813 #else
 814 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 815     MPI_Allreduce(MPI_IN_PLACE,r,nr,MPI_INT,MPI_SUM,ms->mpi_comm_masters);
 816 #else
 817     /* this is thread-unsafe, but it will do for now: */
 818     int i;
 819
 820     if (nr > ms->mpb->ibuf_alloc) {
 821         ms->mpb->ibuf_alloc = nr;
 822         srenew(ms->mpb->ibuf,ms->mpb->ibuf_alloc);
 823     }
 824     MPI_Allreduce(r,ms->mpb->ibuf,nr,MPI_INT,MPI_SUM,ms->mpi_comm_masters);
 825     for(i=0; i<nr; i++)
 826         r[i] = ms->mpb->ibuf[i];
 827 #endif
 828 #endif
 829 }
 830
 831 void gmx_sumli_sim(int nr,gmx_large_int_t r[], const gmx_multisim_t *ms)
 832 {
 833 #ifndef GMX_MPI
 834     gmx_call("gmx_sumli_sim");
 835 #else
 836 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 837     MPI_Allreduce(MPI_IN_PLACE,r,nr,GMX_MPI_LARGE_INT,MPI_SUM,
 838                   ms->mpi_comm_masters);
 839 #else
 840     /* this is thread-unsafe, but it will do for now: */
 841     int i;
 842
 843     if (nr > ms->mpb->libuf_alloc) {
 844         ms->mpb->libuf_alloc = nr;
 845         srenew(ms->mpb->libuf,ms->mpb->libuf_alloc);
 846     }
 847     MPI_Allreduce(r,ms->mpb->libuf,nr,GMX_MPI_LARGE_INT,MPI_SUM,
 848                   ms->mpi_comm_masters);
 849     for(i=0; i<nr; i++)
 850         r[i] = ms->mpb->libuf[i];
 851 #endif
 852 #endif
 853 }
 854
 855
 856 void gmx_finalize_par(void)
 857 {
 858 #ifndef GMX_MPI
 859     /* Compiled without MPI, no MPI finalizing needed */
 860     return;
 861 #else
 862     int initialized,finalized;
 863     int ret;
 864
 865     MPI_Initialized(&initialized);
 866     if (!initialized)
 867     {
 868         return;
 869     }
 870     /* just as a check; we don't want to finalize twice */
 871     MPI_Finalized(&finalized);
 872     if (finalized)
 873     {
 874         return;
 875     }
 876
 877   /* We sync the processes here to try to avoid problems
 878    * with buggy MPI implementations that could cause
 879    * unfinished processes to terminate.
 880    */
 881   MPI_Barrier(MPI_COMM_WORLD);
 882
 883   /*
 884   if (DOMAINDECOMP(cr)) {
 885     if (cr->npmenodes > 0 || cr->dd->bCartesian)
 886       MPI_Comm_free(&cr->mpi_comm_mygroup);
 887     if (cr->dd->bCartesian)
 888       MPI_Comm_free(&cr->mpi_comm_mysim);
 889   }
 890   */
 891
 892   /* Apparently certain mpich implementations cause problems
 893    * with MPI_Finalize. In that case comment out MPI_Finalize.
 894    */
 895   if (debug)
 896     fprintf(debug,"Will call MPI_Finalize now\n");
 897
 898   ret = MPI_Finalize();
 899   if (debug)
 900     fprintf(debug,"Return code from MPI_Finalize = %d\n",ret);
 901 #endif
 902 }
 903