src/gmxlib/network.c

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
   5  * Copyright (c) 2001-2004, The GROMACS development team,
   6  * check out http://www.gromacs.org for more information.
   7  * Copyright (c) 2012,2013, by the GROMACS development team, led by
   8  * David van der Spoel, Berk Hess, Erik Lindahl, and including many
   9  * others, as listed in the AUTHORS file in the top-level source
  10  * directory and at http://www.gromacs.org.
  11  *
  12  * GROMACS is free software; you can redistribute it and/or
  13  * modify it under the terms of the GNU Lesser General Public License
  14  * as published by the Free Software Foundation; either version 2.1
  15  * of the License, or (at your option) any later version.
  16  *
  17  * GROMACS is distributed in the hope that it will be useful,
  18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  20  * Lesser General Public License for more details.
  21  *
  22  * You should have received a copy of the GNU Lesser General Public
  23  * License along with GROMACS; if not, see
  24  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  25  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  26  *
  27  * If you want to redistribute modifications to GROMACS, please
  28  * consider that scientific software is very special. Version
  29  * control is crucial - bugs must be traceable. We will be happy to
  30  * consider code for inclusion in the official distribution, but
  31  * derived work must not be called official GROMACS. Details are found
  32  * in the README & COPYING files - if they are missing, get the
  33  * official version at http://www.gromacs.org.
  34  *
  35  * To help us fund GROMACS development, we humbly ask that you cite
  36  * the research papers on the package. Check out http://www.gromacs.org.
  37  */
  38 #ifdef HAVE_CONFIG_H
  39 #include <config.h>
  40 #endif
  41
  42 #include <string.h>
  43 #include "gmx_fatal.h"
  44 #include "main.h"
  45 #include "smalloc.h"
  46 #include "network.h"
  47 #include "copyrite.h"
  48 #include "statutil.h"
  49 #include <ctype.h>
  50 #include "macros.h"
  51
  52 #ifdef GMX_LIB_MPI
  53 #include <mpi.h>
  54 #endif
  55
  56 #ifdef GMX_THREAD_MPI
  57 #include "tmpi.h"
  58 #endif
  59
  60 #include "mpelogging.h"
  61
  62 /* The source code in this file should be thread-safe.
  63       Please keep it that way. */
  64
  65 gmx_bool gmx_mpi_initialized(void)
  66 {
  67     int n;
  68 #ifndef GMX_MPI
  69     return 0;
  70 #else
  71     MPI_Initialized(&n);
  72
  73     return n;
  74 #endif
  75 }
  76
  77 int gmx_setup(int *argc, char **argv, int *nnodes)
  78 {
  79 #ifndef GMX_MPI
  80     gmx_call("gmx_setup");
  81     return 0;
  82 #else
  83     char   buf[256];
  84     int    resultlen;             /* actual length of node name      */
  85     int    i, flag;
  86     int    mpi_num_nodes;
  87     int    mpi_my_rank;
  88     char   mpi_hostname[MPI_MAX_PROCESSOR_NAME];
  89
  90     /* Call the MPI routines */
  91 #ifdef GMX_LIB_MPI
  92 #ifdef GMX_FAHCORE
  93     (void) fah_MPI_Init(argc, &argv);
  94 #else
  95     (void) MPI_Init(argc, &argv);
  96 #endif
  97 #endif
  98     (void) MPI_Comm_size( MPI_COMM_WORLD, &mpi_num_nodes );
  99     (void) MPI_Comm_rank( MPI_COMM_WORLD, &mpi_my_rank );
 100     (void) MPI_Get_processor_name( mpi_hostname, &resultlen );
 101
 102
 103 #ifdef USE_MPE
 104     /* MPE logging routines. Get event IDs from MPE: */
 105     /* General events */
 106     ev_timestep1               = MPE_Log_get_event_number( );
 107     ev_timestep2               = MPE_Log_get_event_number( );
 108     ev_force_start             = MPE_Log_get_event_number( );
 109     ev_force_finish            = MPE_Log_get_event_number( );
 110     ev_do_fnbf_start           = MPE_Log_get_event_number( );
 111     ev_do_fnbf_finish          = MPE_Log_get_event_number( );
 112     ev_ns_start                = MPE_Log_get_event_number( );
 113     ev_ns_finish               = MPE_Log_get_event_number( );
 114     ev_calc_bonds_start        = MPE_Log_get_event_number( );
 115     ev_calc_bonds_finish       = MPE_Log_get_event_number( );
 116     ev_global_stat_start       = MPE_Log_get_event_number( );
 117     ev_global_stat_finish      = MPE_Log_get_event_number( );
 118     ev_virial_start            = MPE_Log_get_event_number( );
 119     ev_virial_finish           = MPE_Log_get_event_number( );
 120
 121     /* Shift related events */
 122     ev_shift_start             = MPE_Log_get_event_number( );
 123     ev_shift_finish            = MPE_Log_get_event_number( );
 124     ev_unshift_start           = MPE_Log_get_event_number( );
 125     ev_unshift_finish          = MPE_Log_get_event_number( );
 126     ev_mk_mshift_start         = MPE_Log_get_event_number( );
 127     ev_mk_mshift_finish        = MPE_Log_get_event_number( );
 128
 129     /* PME related events */
 130     ev_pme_start                = MPE_Log_get_event_number( );
 131     ev_pme_finish               = MPE_Log_get_event_number( );
 132     ev_spread_on_grid_start     = MPE_Log_get_event_number( );
 133     ev_spread_on_grid_finish    = MPE_Log_get_event_number( );
 134     ev_sum_qgrid_start          = MPE_Log_get_event_number( );
 135     ev_sum_qgrid_finish         = MPE_Log_get_event_number( );
 136     ev_gmxfft3d_start           = MPE_Log_get_event_number( );
 137     ev_gmxfft3d_finish          = MPE_Log_get_event_number( );
 138     ev_solve_pme_start          = MPE_Log_get_event_number( );
 139     ev_solve_pme_finish         = MPE_Log_get_event_number( );
 140     ev_gather_f_bsplines_start  = MPE_Log_get_event_number( );
 141     ev_gather_f_bsplines_finish = MPE_Log_get_event_number( );
 142     ev_reduce_start             = MPE_Log_get_event_number( );
 143     ev_reduce_finish            = MPE_Log_get_event_number( );
 144     ev_rscatter_start           = MPE_Log_get_event_number( );
 145     ev_rscatter_finish          = MPE_Log_get_event_number( );
 146     ev_alltoall_start           = MPE_Log_get_event_number( );
 147     ev_alltoall_finish          = MPE_Log_get_event_number( );
 148     ev_pmeredist_start          = MPE_Log_get_event_number( );
 149     ev_pmeredist_finish         = MPE_Log_get_event_number( );
 150     ev_init_pme_start           = MPE_Log_get_event_number( );
 151     ev_init_pme_finish          = MPE_Log_get_event_number( );
 152     ev_send_coordinates_start   = MPE_Log_get_event_number( );
 153     ev_send_coordinates_finish  = MPE_Log_get_event_number( );
 154     ev_update_fr_start          = MPE_Log_get_event_number( );
 155     ev_update_fr_finish         = MPE_Log_get_event_number( );
 156     ev_clear_rvecs_start        = MPE_Log_get_event_number( );
 157     ev_clear_rvecs_finish       = MPE_Log_get_event_number( );
 158     ev_update_start             = MPE_Log_get_event_number( );
 159     ev_update_finish            = MPE_Log_get_event_number( );
 160     ev_output_start             = MPE_Log_get_event_number( );
 161     ev_output_finish            = MPE_Log_get_event_number( );
 162     ev_sum_lrforces_start       = MPE_Log_get_event_number( );
 163     ev_sum_lrforces_finish      = MPE_Log_get_event_number( );
 164     ev_sort_start               = MPE_Log_get_event_number( );
 165     ev_sort_finish              = MPE_Log_get_event_number( );
 166     ev_sum_qgrid_start          = MPE_Log_get_event_number( );
 167     ev_sum_qgrid_finish         = MPE_Log_get_event_number( );
 168
 169     /* Essential dynamics related events */
 170     ev_edsam_start             = MPE_Log_get_event_number( );
 171     ev_edsam_finish            = MPE_Log_get_event_number( );
 172     ev_get_coords_start        = MPE_Log_get_event_number( );
 173     ev_get_coords_finish       = MPE_Log_get_event_number( );
 174     ev_ed_apply_cons_start     = MPE_Log_get_event_number( );
 175     ev_ed_apply_cons_finish    = MPE_Log_get_event_number( );
 176     ev_fit_to_reference_start  = MPE_Log_get_event_number( );
 177     ev_fit_to_reference_finish = MPE_Log_get_event_number( );
 178
 179     /* describe events: */
 180     if (mpi_my_rank == 0)
 181     {
 182         /* General events */
 183         MPE_Describe_state(ev_timestep1,               ev_timestep2,                "timestep START",  "magenta" );
 184         MPE_Describe_state(ev_force_start,             ev_force_finish,             "force",           "cornflower blue" );
 185         MPE_Describe_state(ev_do_fnbf_start,           ev_do_fnbf_finish,           "do_fnbf",         "navy" );
 186         MPE_Describe_state(ev_ns_start,                ev_ns_finish,                "neighbor search", "tomato" );
 187         MPE_Describe_state(ev_calc_bonds_start,        ev_calc_bonds_finish,        "bonded forces",   "slate blue" );
 188         MPE_Describe_state(ev_global_stat_start,       ev_global_stat_finish,       "global stat",     "firebrick3");
 189         MPE_Describe_state(ev_update_fr_start,         ev_update_fr_finish,         "update forcerec", "goldenrod");
 190         MPE_Describe_state(ev_clear_rvecs_start,       ev_clear_rvecs_finish,       "clear rvecs",     "bisque");
 191         MPE_Describe_state(ev_update_start,            ev_update_finish,            "update",          "cornsilk");
 192         MPE_Describe_state(ev_output_start,            ev_output_finish,            "output",          "black");
 193         MPE_Describe_state(ev_virial_start,            ev_virial_finish,            "calc_virial",     "thistle4");
 194
 195         /* PME related events */
 196         MPE_Describe_state(ev_pme_start,               ev_pme_finish,               "doing PME",       "grey" );
 197         MPE_Describe_state(ev_spread_on_grid_start,    ev_spread_on_grid_finish,    "spread",          "dark orange" );
 198         MPE_Describe_state(ev_sum_qgrid_start,         ev_sum_qgrid_finish,         "sum qgrid",       "slate blue");
 199         MPE_Describe_state(ev_gmxfft3d_start,          ev_gmxfft3d_finish,          "fft3d",           "snow2" );
 200         MPE_Describe_state(ev_solve_pme_start,         ev_solve_pme_finish,         "solve PME",       "indian red" );
 201         MPE_Describe_state(ev_gather_f_bsplines_start, ev_gather_f_bsplines_finish, "bsplines",        "light sea green" );
 202         MPE_Describe_state(ev_reduce_start,            ev_reduce_finish,            "reduce",          "cyan1" );
 203         MPE_Describe_state(ev_rscatter_start,          ev_rscatter_finish,          "rscatter",        "cyan3" );
 204         MPE_Describe_state(ev_alltoall_start,          ev_alltoall_finish,          "alltoall",        "LightCyan4" );
 205         MPE_Describe_state(ev_pmeredist_start,         ev_pmeredist_finish,         "pmeredist",       "thistle" );
 206         MPE_Describe_state(ev_init_pme_start,          ev_init_pme_finish,          "init PME",        "snow4");
 207         MPE_Describe_state(ev_send_coordinates_start,  ev_send_coordinates_finish,  "send_coordinates", "blue");
 208         MPE_Describe_state(ev_sum_lrforces_start,      ev_sum_lrforces_finish,      "sum_LRforces",    "lime green");
 209         MPE_Describe_state(ev_sort_start,              ev_sort_finish,              "sort pme atoms",  "brown");
 210         MPE_Describe_state(ev_sum_qgrid_start,         ev_sum_qgrid_finish,         "sum charge grid", "medium orchid");
 211
 212         /* Shift related events */
 213         MPE_Describe_state(ev_shift_start,             ev_shift_finish,             "shift",           "orange");
 214         MPE_Describe_state(ev_unshift_start,           ev_unshift_finish,           "unshift",         "dark orange");
 215         MPE_Describe_state(ev_mk_mshift_start,         ev_mk_mshift_finish,         "mk_mshift",       "maroon");
 216
 217         /* Essential dynamics related events */
 218         MPE_Describe_state(ev_edsam_start,             ev_edsam_finish,             "EDSAM",           "deep sky blue");
 219         MPE_Describe_state(ev_get_coords_start,        ev_get_coords_finish,        "ED get coords",   "steel blue");
 220         MPE_Describe_state(ev_ed_apply_cons_start,     ev_ed_apply_cons_finish,     "ED apply constr", "forest green");
 221         MPE_Describe_state(ev_fit_to_reference_start,  ev_fit_to_reference_finish,  "ED fit to ref",   "lavender");
 222
 223     }
 224     MPE_Init_log();
 225 #endif
 226
 227 #ifdef GMX_LIB_MPI
 228     if (debug)
 229     {
 230         fprintf(debug, "NNODES=%d, MYRANK=%d, HOSTNAME=%s\n",
 231                 mpi_num_nodes, mpi_my_rank, mpi_hostname);
 232     }
 233 #endif
 234
 235     *nnodes = mpi_num_nodes;
 236
 237     return mpi_my_rank;
 238 #endif
 239 }
 240
 241 int  gmx_node_num(void)
 242 {
 243 #ifndef GMX_MPI
 244     return 1;
 245 #else
 246     int i;
 247     (void) MPI_Comm_size(MPI_COMM_WORLD, &i);
 248     return i;
 249 #endif
 250 }
 251
 252 int gmx_node_rank(void)
 253 {
 254 #ifndef GMX_MPI
 255     return 0;
 256 #else
 257     int i;
 258     (void) MPI_Comm_rank(MPI_COMM_WORLD, &i);
 259     return i;
 260 #endif
 261 }
 262
 263
 264 int gmx_hostname_num()
 265 {
 266 #ifndef GMX_MPI
 267     return 0;
 268 #else
 269 #ifdef GMX_THREAD_MPI
 270     /* thread-MPI currently puts the thread number in the process name,
 271      * we might want to change this, as this is inconsistent with what
 272      * most MPI implementations would do when running on a single node.
 273      */
 274     return 0;
 275 #else
 276     int  resultlen, hostnum, i, j;
 277     char mpi_hostname[MPI_MAX_PROCESSOR_NAME], hostnum_str[MPI_MAX_PROCESSOR_NAME];
 278
 279     MPI_Get_processor_name(mpi_hostname, &resultlen);
 280     /* This procedure can only differentiate nodes with host names
 281      * that end on unique numbers.
 282      */
 283     i = 0;
 284     j = 0;
 285     /* Only parse the host name up to the first dot */
 286     while (i < resultlen && mpi_hostname[i] != '.')
 287     {
 288         if (isdigit(mpi_hostname[i]))
 289         {
 290             hostnum_str[j++] = mpi_hostname[i];
 291         }
 292         i++;
 293     }
 294     hostnum_str[j] = '\0';
 295     if (j == 0)
 296     {
 297         hostnum = 0;
 298     }
 299     else
 300     {
 301         /* Use only the last 9 decimals, so we don't overflow an int */
 302         hostnum = strtol(hostnum_str + max(0, j-9), NULL, 10);
 303     }
 304
 305     if (debug)
 306     {
 307         fprintf(debug, "In gmx_setup_nodecomm: hostname '%s', hostnum %d\n",
 308                 mpi_hostname, hostnum);
 309     }
 310     return hostnum;
 311 #endif
 312 #endif
 313 }
 314
 315 void gmx_setup_nodecomm(FILE *fplog, t_commrec *cr)
 316 {
 317     gmx_nodecomm_t *nc;
 318     int             n, rank, hostnum, ng, ni;
 319
 320     /* Many MPI implementations do not optimize MPI_Allreduce
 321      * (and probably also other global communication calls)
 322      * for multi-core nodes connected by a network.
 323      * We can optimize such communication by using one MPI call
 324      * within each node and one between the nodes.
 325      * For MVAPICH2 and Intel MPI this reduces the time for
 326      * the global_stat communication by 25%
 327      * for 2x2-core 3 GHz Woodcrest connected by mixed DDR/SDR Infiniband.
 328      * B. Hess, November 2007
 329      */
 330
 331     nc = &cr->nc;
 332
 333     nc->bUse = FALSE;
 334 #ifndef GMX_THREAD_MPI
 335 #ifdef GMX_MPI
 336     MPI_Comm_size(cr->mpi_comm_mygroup, &n);
 337     MPI_Comm_rank(cr->mpi_comm_mygroup, &rank);
 338
 339     hostnum = gmx_hostname_num();
 340
 341     if (debug)
 342     {
 343         fprintf(debug, "In gmx_setup_nodecomm: splitting communicator of size %d\n", n);
 344     }
 345
 346
 347     /* The intra-node communicator, split on node number */
 348     MPI_Comm_split(cr->mpi_comm_mygroup, hostnum, rank, &nc->comm_intra);
 349     MPI_Comm_rank(nc->comm_intra, &nc->rank_intra);
 350     if (debug)
 351     {
 352         fprintf(debug, "In gmx_setup_nodecomm: node rank %d rank_intra %d\n",
 353                 rank, nc->rank_intra);
 354     }
 355     /* The inter-node communicator, split on rank_intra.
 356      * We actually only need the one for rank=0,
 357      * but it is easier to create them all.
 358      */
 359     MPI_Comm_split(cr->mpi_comm_mygroup, nc->rank_intra, rank, &nc->comm_inter);
 360     /* Check if this really created two step communication */
 361     MPI_Comm_size(nc->comm_inter, &ng);
 362     MPI_Comm_size(nc->comm_intra, &ni);
 363     if (debug)
 364     {
 365         fprintf(debug, "In gmx_setup_nodecomm: groups %d, my group size %d\n",
 366                 ng, ni);
 367     }
 368
 369     if (getenv("GMX_NO_NODECOMM") == NULL &&
 370         ((ng > 1 && ng < n) || (ni > 1 && ni < n)))
 371     {
 372         nc->bUse = TRUE;
 373         if (fplog)
 374         {
 375             fprintf(fplog, "Using two step summing over %d groups of on average %.1f processes\n\n",
 376                     ng, (real)n/(real)ng);
 377         }
 378         if (nc->rank_intra > 0)
 379         {
 380             MPI_Comm_free(&nc->comm_inter);
 381         }
 382     }
 383     else
 384     {
 385         /* One group or all processes in a separate group, use normal summing */
 386         MPI_Comm_free(&nc->comm_inter);
 387         MPI_Comm_free(&nc->comm_intra);
 388         if (debug)
 389         {
 390             fprintf(debug, "In gmx_setup_nodecomm: not unsing separate inter- and intra-node communicators.\n");
 391         }
 392     }
 393 #endif
 394 #else
 395     /* tMPI runs only on a single node so just use the nodeid */
 396     nc->rank_intra = cr->nodeid;
 397 #endif
 398 }
 399
 400 void gmx_init_intranode_counters(t_commrec *cr)
 401 {
 402     /* counters for PP+PME and PP-only processes on my physical node */
 403     int nrank_intranode, rank_intranode;
 404     int nrank_pp_intranode, rank_pp_intranode;
 405     /* thread-MPI is not initialized when not running in parallel */
 406 #if defined GMX_MPI && !defined GMX_THREAD_MPI
 407     int nrank_world, rank_world;
 408     int i, mynum, *num, *num_s, *num_pp, *num_pp_s;
 409
 410     MPI_Comm_size(MPI_COMM_WORLD, &nrank_world);
 411     MPI_Comm_rank(MPI_COMM_WORLD, &rank_world);
 412
 413     /* Get the node number from the hostname to identify the nodes */
 414     mynum = gmx_hostname_num();
 415
 416     /* We can't rely on MPI_IN_PLACE, so we need send and receive buffers */
 417     snew(num,   nrank_world);
 418     snew(num_s, nrank_world);
 419     snew(num_pp,   nrank_world);
 420     snew(num_pp_s, nrank_world);
 421
 422     num_s[rank_world]    = mynum;
 423     num_pp_s[rank_world] = (cr->duty & DUTY_PP) ? mynum : -1;
 424
 425     MPI_Allreduce(num_s,    num,    nrank_world, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
 426     MPI_Allreduce(num_pp_s, num_pp, nrank_world, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
 427
 428     nrank_intranode    = 0;
 429     rank_intranode     = 0;
 430     nrank_pp_intranode = 0;
 431     rank_pp_intranode  = 0;
 432     for (i = 0; i < nrank_world; i++)
 433     {
 434         if (num[i] == mynum)
 435         {
 436             nrank_intranode++;
 437             if (i < rank_world)
 438             {
 439                 rank_intranode++;
 440             }
 441         }
 442         if ((cr->duty & DUTY_PP) && num_pp[i] == mynum)
 443         {
 444             nrank_pp_intranode++;
 445             if (i < rank_world)
 446             {
 447                 rank_pp_intranode++;
 448             }
 449         }
 450     }
 451     sfree(num);
 452     sfree(num_s);
 453     sfree(num_pp);
 454     sfree(num_pp_s);
 455 #else
 456     /* Serial or thread-MPI code: we run within a single physical node */
 457     nrank_intranode    = cr->nnodes;
 458     rank_intranode     = cr->sim_nodeid;
 459     nrank_pp_intranode = cr->nnodes - cr->npmenodes;
 460     rank_pp_intranode  = cr->nodeid;
 461 #endif
 462
 463     if (debug)
 464     {
 465         char sbuf[STRLEN];
 466         if (cr->duty & DUTY_PP && cr->duty & DUTY_PME)
 467         {
 468             sprintf(sbuf, "PP+PME");
 469         }
 470         else
 471         {
 472             sprintf(sbuf, "%s", cr->duty & DUTY_PP ? "PP" : "PME");
 473         }
 474         fprintf(debug, "On %3s node %d: nrank_intranode=%d, rank_intranode=%d, "
 475                 "nrank_pp_intranode=%d, rank_pp_intranode=%d\n",
 476                 sbuf, cr->sim_nodeid,
 477                 nrank_intranode, rank_intranode,
 478                 nrank_pp_intranode, rank_pp_intranode);
 479     }
 480
 481     cr->nrank_intranode    = nrank_intranode;
 482     cr->rank_intranode     = rank_intranode;
 483     cr->nrank_pp_intranode = nrank_pp_intranode;
 484     cr->rank_pp_intranode  = rank_pp_intranode;
 485 }
 486
 487
 488 void gmx_barrier(const t_commrec *cr)
 489 {
 490 #ifndef GMX_MPI
 491     gmx_call("gmx_barrier");
 492 #else
 493     MPI_Barrier(cr->mpi_comm_mygroup);
 494 #endif
 495 }
 496
 497 void gmx_abort(int noderank, int nnodes, int errorno)
 498 {
 499 #ifndef GMX_MPI
 500     gmx_call("gmx_abort");
 501 #else
 502 #ifdef GMX_THREAD_MPI
 503     fprintf(stderr, "Halting program %s\n", ShortProgram());
 504     thanx(stderr);
 505     exit(1);
 506 #else
 507     if (nnodes > 1)
 508     {
 509         fprintf(stderr, "Halting parallel program %s on CPU %d out of %d\n",
 510                 ShortProgram(), noderank, nnodes);
 511     }
 512     else
 513     {
 514         fprintf(stderr, "Halting program %s\n", ShortProgram());
 515     }
 516
 517     thanx(stderr);
 518     MPI_Abort(MPI_COMM_WORLD, errorno);
 519     exit(1);
 520 #endif
 521 #endif
 522 }
 523
 524 void gmx_bcast(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_mygroup);
 530 #endif
 531 }
 532
 533 void gmx_bcast_sim(int nbytes, void *b, const t_commrec *cr)
 534 {
 535 #ifndef GMX_MPI
 536     gmx_call("gmx_bast");
 537 #else
 538     MPI_Bcast(b, nbytes, MPI_BYTE, MASTERRANK(cr), cr->mpi_comm_mysim);
 539 #endif
 540 }
 541
 542 void gmx_sumd(int nr, double r[], const t_commrec *cr)
 543 {
 544 #ifndef GMX_MPI
 545     gmx_call("gmx_sumd");
 546 #else
 547 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 548     if (cr->nc.bUse)
 549     {
 550         if (cr->nc.rank_intra == 0)
 551         {
 552             /* Use two step summing. */
 553             MPI_Reduce(MPI_IN_PLACE, r, nr, MPI_DOUBLE, MPI_SUM, 0,
 554                        cr->nc.comm_intra);
 555             /* Sum the roots of the internal (intra) buffers. */
 556             MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_DOUBLE, MPI_SUM,
 557                           cr->nc.comm_inter);
 558         }
 559         else
 560         {
 561             /* This is here because of the silly MPI specification
 562                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 563             MPI_Reduce(r, NULL, nr, MPI_DOUBLE, MPI_SUM, 0, cr->nc.comm_intra);
 564         }
 565         MPI_Bcast(r, nr, MPI_DOUBLE, 0, cr->nc.comm_intra);
 566     }
 567     else
 568     {
 569         MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_DOUBLE, MPI_SUM,
 570                       cr->mpi_comm_mygroup);
 571     }
 572 #else
 573     int i;
 574
 575     if (nr > cr->mpb->dbuf_alloc)
 576     {
 577         cr->mpb->dbuf_alloc = nr;
 578         srenew(cr->mpb->dbuf, cr->mpb->dbuf_alloc);
 579     }
 580     if (cr->nc.bUse)
 581     {
 582         /* Use two step summing */
 583         MPI_Allreduce(r, cr->mpb->dbuf, nr, MPI_DOUBLE, MPI_SUM, cr->nc.comm_intra);
 584         if (cr->nc.rank_intra == 0)
 585         {
 586             /* Sum with the buffers reversed */
 587             MPI_Allreduce(cr->mpb->dbuf, r, nr, MPI_DOUBLE, MPI_SUM,
 588                           cr->nc.comm_inter);
 589         }
 590         MPI_Bcast(r, nr, MPI_DOUBLE, 0, cr->nc.comm_intra);
 591     }
 592     else
 593     {
 594         MPI_Allreduce(r, cr->mpb->dbuf, nr, MPI_DOUBLE, MPI_SUM,
 595                       cr->mpi_comm_mygroup);
 596         for (i = 0; i < nr; i++)
 597         {
 598             r[i] = cr->mpb->dbuf[i];
 599         }
 600     }
 601 #endif
 602 #endif
 603 }
 604
 605 void gmx_sumf(int nr, float r[], const t_commrec *cr)
 606 {
 607 #ifndef GMX_MPI
 608     gmx_call("gmx_sumf");
 609 #else
 610 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 611     if (cr->nc.bUse)
 612     {
 613         /* Use two step summing.  */
 614         if (cr->nc.rank_intra == 0)
 615         {
 616             MPI_Reduce(MPI_IN_PLACE, r, nr, MPI_FLOAT, MPI_SUM, 0,
 617                        cr->nc.comm_intra);
 618             /* Sum the roots of the internal (intra) buffers */
 619             MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_FLOAT, MPI_SUM,
 620                           cr->nc.comm_inter);
 621         }
 622         else
 623         {
 624             /* This is here because of the silly MPI specification
 625                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 626             MPI_Reduce(r, NULL, nr, MPI_FLOAT, MPI_SUM, 0, cr->nc.comm_intra);
 627         }
 628         MPI_Bcast(r, nr, MPI_FLOAT, 0, cr->nc.comm_intra);
 629     }
 630     else
 631     {
 632         MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_FLOAT, MPI_SUM, cr->mpi_comm_mygroup);
 633     }
 634 #else
 635     int i;
 636
 637     if (nr > cr->mpb->fbuf_alloc)
 638     {
 639         cr->mpb->fbuf_alloc = nr;
 640         srenew(cr->mpb->fbuf, cr->mpb->fbuf_alloc);
 641     }
 642     if (cr->nc.bUse)
 643     {
 644         /* Use two step summing */
 645         MPI_Allreduce(r, cr->mpb->fbuf, nr, MPI_FLOAT, MPI_SUM, cr->nc.comm_intra);
 646         if (cr->nc.rank_intra == 0)
 647         {
 648             /* Sum with the buffers reversed */
 649             MPI_Allreduce(cr->mpb->fbuf, r, nr, MPI_FLOAT, MPI_SUM,
 650                           cr->nc.comm_inter);
 651         }
 652         MPI_Bcast(r, nr, MPI_FLOAT, 0, cr->nc.comm_intra);
 653     }
 654     else
 655     {
 656         MPI_Allreduce(r, cr->mpb->fbuf, nr, MPI_FLOAT, MPI_SUM,
 657                       cr->mpi_comm_mygroup);
 658         for (i = 0; i < nr; i++)
 659         {
 660             r[i] = cr->mpb->fbuf[i];
 661         }
 662     }
 663 #endif
 664 #endif
 665 }
 666
 667 void gmx_sumi(int nr, int r[], const t_commrec *cr)
 668 {
 669 #ifndef GMX_MPI
 670     gmx_call("gmx_sumi");
 671 #else
 672 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 673     if (cr->nc.bUse)
 674     {
 675         /* Use two step summing */
 676         if (cr->nc.rank_intra == 0)
 677         {
 678             MPI_Reduce(MPI_IN_PLACE, r, nr, MPI_INT, MPI_SUM, 0, cr->nc.comm_intra);
 679             /* Sum with the buffers reversed */
 680             MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_INT, MPI_SUM, cr->nc.comm_inter);
 681         }
 682         else
 683         {
 684             /* This is here because of the silly MPI specification
 685                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 686             MPI_Reduce(r, NULL, nr, MPI_INT, MPI_SUM, 0, cr->nc.comm_intra);
 687         }
 688         MPI_Bcast(r, nr, MPI_INT, 0, cr->nc.comm_intra);
 689     }
 690     else
 691     {
 692         MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_INT, MPI_SUM, cr->mpi_comm_mygroup);
 693     }
 694 #else
 695     int i;
 696
 697     if (nr > cr->mpb->ibuf_alloc)
 698     {
 699         cr->mpb->ibuf_alloc = nr;
 700         srenew(cr->mpb->ibuf, cr->mpb->ibuf_alloc);
 701     }
 702     if (cr->nc.bUse)
 703     {
 704         /* Use two step summing */
 705         MPI_Allreduce(r, cr->mpb->ibuf, nr, MPI_INT, MPI_SUM, cr->nc.comm_intra);
 706         if (cr->nc.rank_intra == 0)
 707         {
 708             /* Sum with the buffers reversed */
 709             MPI_Allreduce(cr->mpb->ibuf, r, nr, MPI_INT, MPI_SUM, cr->nc.comm_inter);
 710         }
 711         MPI_Bcast(r, nr, MPI_INT, 0, cr->nc.comm_intra);
 712     }
 713     else
 714     {
 715         MPI_Allreduce(r, cr->mpb->ibuf, nr, MPI_INT, MPI_SUM, cr->mpi_comm_mygroup);
 716         for (i = 0; i < nr; i++)
 717         {
 718             r[i] = cr->mpb->ibuf[i];
 719         }
 720     }
 721 #endif
 722 #endif
 723 }
 724
 725 void gmx_sumli(int nr, gmx_large_int_t r[], const t_commrec *cr)
 726 {
 727 #ifndef GMX_MPI
 728     gmx_call("gmx_sumli");
 729 #else
 730 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 731     if (cr->nc.bUse)
 732     {
 733         /* Use two step summing */
 734         if (cr->nc.rank_intra == 0)
 735         {
 736             MPI_Reduce(MPI_IN_PLACE, r, nr, GMX_MPI_LARGE_INT, MPI_SUM, 0,
 737                        cr->nc.comm_intra);
 738             /* Sum with the buffers reversed */
 739             MPI_Allreduce(MPI_IN_PLACE, r, nr, GMX_MPI_LARGE_INT, MPI_SUM,
 740                           cr->nc.comm_inter);
 741         }
 742         else
 743         {
 744             /* This is here because of the silly MPI specification
 745                 that MPI_IN_PLACE should be put in sendbuf instead of recvbuf */
 746             MPI_Reduce(r, NULL, nr, GMX_MPI_LARGE_INT, MPI_SUM, 0, cr->nc.comm_intra);
 747         }
 748         MPI_Bcast(r, nr, GMX_MPI_LARGE_INT, 0, cr->nc.comm_intra);
 749     }
 750     else
 751     {
 752         MPI_Allreduce(MPI_IN_PLACE, r, nr, GMX_MPI_LARGE_INT, MPI_SUM, cr->mpi_comm_mygroup);
 753     }
 754 #else
 755     int i;
 756
 757     if (nr > cr->mpb->libuf_alloc)
 758     {
 759         cr->mpb->libuf_alloc = nr;
 760         srenew(cr->mpb->libuf, cr->mpb->libuf_alloc);
 761     }
 762     if (cr->nc.bUse)
 763     {
 764         /* Use two step summing */
 765         MPI_Allreduce(r, cr->mpb->libuf, nr, GMX_MPI_LARGE_INT, MPI_SUM,
 766                       cr->nc.comm_intra);
 767         if (cr->nc.rank_intra == 0)
 768         {
 769             /* Sum with the buffers reversed */
 770             MPI_Allreduce(cr->mpb->libuf, r, nr, GMX_MPI_LARGE_INT, MPI_SUM,
 771                           cr->nc.comm_inter);
 772         }
 773         MPI_Bcast(r, nr, GMX_MPI_LARGE_INT, 0, cr->nc.comm_intra);
 774     }
 775     else
 776     {
 777         MPI_Allreduce(r, cr->mpb->libuf, nr, GMX_MPI_LARGE_INT, MPI_SUM,
 778                       cr->mpi_comm_mygroup);
 779         for (i = 0; i < nr; i++)
 780         {
 781             r[i] = cr->mpb->libuf[i];
 782         }
 783     }
 784 #endif
 785 #endif
 786 }
 787
 788
 789
 790 #ifdef GMX_MPI
 791 void gmx_sumd_comm(int nr, double r[], MPI_Comm mpi_comm)
 792 {
 793 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 794     MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_DOUBLE, MPI_SUM, mpi_comm);
 795 #else
 796     /* this function is only used in code that is not performance critical,
 797        (during setup, when comm_rec is not the appropriate communication
 798        structure), so this isn't as bad as it looks. */
 799     double *buf;
 800     int     i;
 801
 802     snew(buf, nr);
 803     MPI_Allreduce(r, buf, nr, MPI_DOUBLE, MPI_SUM, mpi_comm);
 804     for (i = 0; i < nr; i++)
 805     {
 806         r[i] = buf[i];
 807     }
 808     sfree(buf);
 809 #endif
 810 }
 811 #endif
 812
 813 #ifdef GMX_MPI
 814 void gmx_sumf_comm(int nr, float r[], MPI_Comm mpi_comm)
 815 {
 816 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 817     MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_FLOAT, MPI_SUM, mpi_comm);
 818 #else
 819     /* this function is only used in code that is not performance critical,
 820        (during setup, when comm_rec is not the appropriate communication
 821        structure), so this isn't as bad as it looks. */
 822     float *buf;
 823     int    i;
 824
 825     snew(buf, nr);
 826     MPI_Allreduce(r, buf, nr, MPI_FLOAT, MPI_SUM, mpi_comm);
 827     for (i = 0; i < nr; i++)
 828     {
 829         r[i] = buf[i];
 830     }
 831     sfree(buf);
 832 #endif
 833 }
 834 #endif
 835
 836 void gmx_sumd_sim(int nr, double r[], const gmx_multisim_t *ms)
 837 {
 838 #ifndef GMX_MPI
 839     gmx_call("gmx_sumd_sim");
 840 #else
 841     gmx_sumd_comm(nr, r, ms->mpi_comm_masters);
 842 #endif
 843 }
 844
 845 void gmx_sumf_sim(int nr, float r[], const gmx_multisim_t *ms)
 846 {
 847 #ifndef GMX_MPI
 848     gmx_call("gmx_sumf_sim");
 849 #else
 850     gmx_sumf_comm(nr, r, ms->mpi_comm_masters);
 851 #endif
 852 }
 853
 854 void gmx_sumi_sim(int nr, int r[], const gmx_multisim_t *ms)
 855 {
 856 #ifndef GMX_MPI
 857     gmx_call("gmx_sumi_sim");
 858 #else
 859 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 860     MPI_Allreduce(MPI_IN_PLACE, r, nr, MPI_INT, MPI_SUM, ms->mpi_comm_masters);
 861 #else
 862     /* this is thread-unsafe, but it will do for now: */
 863     int i;
 864
 865     if (nr > ms->mpb->ibuf_alloc)
 866     {
 867         ms->mpb->ibuf_alloc = nr;
 868         srenew(ms->mpb->ibuf, ms->mpb->ibuf_alloc);
 869     }
 870     MPI_Allreduce(r, ms->mpb->ibuf, nr, MPI_INT, MPI_SUM, ms->mpi_comm_masters);
 871     for (i = 0; i < nr; i++)
 872     {
 873         r[i] = ms->mpb->ibuf[i];
 874     }
 875 #endif
 876 #endif
 877 }
 878
 879 void gmx_sumli_sim(int nr, gmx_large_int_t r[], const gmx_multisim_t *ms)
 880 {
 881 #ifndef GMX_MPI
 882     gmx_call("gmx_sumli_sim");
 883 #else
 884 #if defined(MPI_IN_PLACE_EXISTS) || defined(GMX_THREAD_MPI)
 885     MPI_Allreduce(MPI_IN_PLACE, r, nr, GMX_MPI_LARGE_INT, MPI_SUM,
 886                   ms->mpi_comm_masters);
 887 #else
 888     /* this is thread-unsafe, but it will do for now: */
 889     int i;
 890
 891     if (nr > ms->mpb->libuf_alloc)
 892     {
 893         ms->mpb->libuf_alloc = nr;
 894         srenew(ms->mpb->libuf, ms->mpb->libuf_alloc);
 895     }
 896     MPI_Allreduce(r, ms->mpb->libuf, nr, GMX_MPI_LARGE_INT, MPI_SUM,
 897                   ms->mpi_comm_masters);
 898     for (i = 0; i < nr; i++)
 899     {
 900         r[i] = ms->mpb->libuf[i];
 901     }
 902 #endif
 903 #endif
 904 }
 905
 906
 907 void gmx_finalize_par(void)
 908 {
 909 #ifndef GMX_MPI
 910     /* Compiled without MPI, no MPI finalizing needed */
 911     return;
 912 #else
 913     int initialized, finalized;
 914     int ret;
 915
 916     MPI_Initialized(&initialized);
 917     if (!initialized)
 918     {
 919         return;
 920     }
 921     /* just as a check; we don't want to finalize twice */
 922     MPI_Finalized(&finalized);
 923     if (finalized)
 924     {
 925         return;
 926     }
 927
 928     /* We sync the processes here to try to avoid problems
 929      * with buggy MPI implementations that could cause
 930      * unfinished processes to terminate.
 931      */
 932     MPI_Barrier(MPI_COMM_WORLD);
 933
 934     /*
 935        if (DOMAINDECOMP(cr)) {
 936        if (cr->npmenodes > 0 || cr->dd->bCartesian)
 937         MPI_Comm_free(&cr->mpi_comm_mygroup);
 938        if (cr->dd->bCartesian)
 939         MPI_Comm_free(&cr->mpi_comm_mysim);
 940        }
 941      */
 942
 943     /* Apparently certain mpich implementations cause problems
 944      * with MPI_Finalize. In that case comment out MPI_Finalize.
 945      */
 946     if (debug)
 947     {
 948         fprintf(debug, "Will call MPI_Finalize now\n");
 949     }
 950
 951     ret = MPI_Finalize();
 952     if (debug)
 953     {
 954         fprintf(debug, "Return code from MPI_Finalize = %d\n", ret);
 955     }
 956 #endif
 957 }