From 5915b8359132652990ff9a269e21f28a2e42b6e8 Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Thu, 5 Sep 2019 21:29:47 +0200 Subject: [PATCH] Clean up PME rank variables Clarified the naming of several PME rank count variables and added a boolean to DDRankSetup that tells whether we use separate PME ranks. Fixes several misuses of npmenodes, all of which did not cause issues in pratice. Todo: Remove npmenodes from t_commrec. Change-Id: I03b4547c975e3ffa354b53b35abcd91b4f3f6b26 --- src/gromacs/domdec/domdec.cpp | 75 ++++++++++--------- src/gromacs/domdec/domdec.h | 2 +- src/gromacs/domdec/domdec_internal.h | 3 +- src/gromacs/domdec/domdec_setup.cpp | 21 +++--- src/gromacs/domdec/partition.cpp | 4 +- src/gromacs/mdlib/gmx_omp_nthreads.cpp | 2 +- .../taskassignment/resourcedivision.cpp | 2 +- 7 files changed, 56 insertions(+), 53 deletions(-) diff --git a/src/gromacs/domdec/domdec.cpp b/src/gromacs/domdec/domdec.cpp index 66ff10c65e..0fab6566ae 100644 --- a/src/gromacs/domdec/domdec.cpp +++ b/src/gromacs/domdec/domdec.cpp @@ -707,7 +707,7 @@ static int ddindex2pmeindex(const DDRankSetup &ddRankSetup, const int ddCellIndex) { const int npp = ddRankSetup.numPPRanks; - const int npme = ddRankSetup.npmenodes; + const int npme = ddRankSetup.numRanksDoingPme; /* Here we assign a PME node to communicate with this DD node * by assuming that the major index of both is x. @@ -721,7 +721,7 @@ static int *dd_interleaved_pme_ranks(const DDRankSetup &ddRankSetup) int *pme_rank; int n, i, p0, p1; - snew(pme_rank, ddRankSetup.npmenodes); + snew(pme_rank, ddRankSetup.numRanksDoingPme); n = 0; for (i = 0; i < ddRankSetup.numPPRanks; i++) { @@ -885,9 +885,9 @@ std::vector get_pme_ddranks(const t_commrec *cr, const int pmenodeid) { const DDRankSetup &ddRankSetup = cr->dd->comm->ddRankSetup; - GMX_RELEASE_ASSERT(ddRankSetup.npmenodes > 0, "This function should only be called when PME-only ranks are in use"); + GMX_RELEASE_ASSERT(ddRankSetup.usePmeOnlyRanks, "This function should only be called when PME-only ranks are in use"); std::vector ddranks; - ddranks.reserve((ddRankSetup.numPPRanks + ddRankSetup.npmenodes - 1)/ddRankSetup.npmenodes); + ddranks.reserve((ddRankSetup.numPPRanks + ddRankSetup.numRanksDoingPme - 1)/ddRankSetup.numRanksDoingPme); for (int x = 0; x < ddRankSetup.numPPCells[XX]; x++) { @@ -923,11 +923,11 @@ std::vector get_pme_ddranks(const t_commrec *cr, static gmx_bool receive_vir_ener(const gmx_domdec_t *dd, const t_commrec *cr) { - gmx_bool bReceive = TRUE; + gmx_bool bReceive = TRUE; - if (cr->npmenodes < dd->nnodes) + const DDRankSetup &ddRankSetup = dd->comm->ddRankSetup; + if (ddRankSetup.usePmeOnlyRanks) { - const DDRankSetup &ddRankSetup = dd->comm->ddRankSetup; if (ddRankSetup.bCartesianPP_PME) { #if GMX_MPI @@ -1010,7 +1010,7 @@ static void init_ddpme(gmx_domdec_t *dd, gmx_ddpme_t *ddpme, int dimind) return; } - const int nso = ddRankSetup.npmenodes/ddpme->nslab; + const int nso = ddRankSetup.numRanksDoingPme/ddpme->nslab; /* Determine for each PME slab the PP location range for dimension dim */ snew(ddpme->pp_min, ddpme->nslab); snew(ddpme->pp_max, ddpme->nslab); @@ -1435,7 +1435,7 @@ static void make_pp_communicator(const gmx::MDLogger &mdlog, } else if (ddRankSetup.bCartesianPP) { - if (cr->npmenodes == 0) + if (!ddRankSetup.usePmeOnlyRanks) { /* The PP communicator is also * the communicator for this simulation @@ -1538,7 +1538,7 @@ static void split_communicator(const gmx::MDLogger &mdlog, bool bDiv[DIM]; for (int i = 1; i < DIM; i++) { - bDiv[i] = ((cr->npmenodes*numDDCells[i]) % numDDCellsTot == 0); + bDiv[i] = ((ddRankSetup->numRanksDoingPme*numDDCells[i]) % numDDCellsTot == 0); } if (bDiv[YY] || bDiv[ZZ]) { @@ -1561,13 +1561,15 @@ static void split_communicator(const gmx::MDLogger &mdlog, ddRankSetup->cartpmedim = YY; } ddRankSetup->ntot[ddRankSetup->cartpmedim] - += (cr->npmenodes*numDDCells[ddRankSetup->cartpmedim])/numDDCellsTot; + += (ddRankSetup->numRanksDoingPme*numDDCells[ddRankSetup->cartpmedim])/numDDCellsTot; } else { GMX_LOG(mdlog.info).appendTextFormatted( "Number of PME-only ranks (%d) is not a multiple of nx*ny (%d*%d) or nx*nz (%d*%d)", - cr->npmenodes, numDDCells[XX], numDDCells[YY], numDDCells[XX], numDDCells[ZZ]); + ddRankSetup->numRanksDoingPme, + numDDCells[XX], numDDCells[YY], + numDDCells[XX], numDDCells[ZZ]); GMX_LOG(mdlog.info).appendText("Will not use a Cartesian communicator for PP <-> PME\n"); } } @@ -1611,7 +1613,7 @@ static void split_communicator(const gmx::MDLogger &mdlog, { cr->duty = DUTY_PP; } - if (cr->npmenodes == 0 || + if (!ddRankSetup->usePmeOnlyRanks || ddCellIndex[ddRankSetup->cartpmedim] >= numDDCells[ddRankSetup->cartpmedim]) { cr->duty = DUTY_PME; @@ -1672,7 +1674,6 @@ static void split_communicator(const gmx::MDLogger &mdlog, */ static void makeGroupCommunicators(const gmx::MDLogger &mdlog, const DDSettings &ddSettings, - const DDSetup &ddSetup, const DdRankOrder ddRankOrder, DDRankSetup *ddRankSetup, t_commrec *cr, @@ -1681,12 +1682,12 @@ static void makeGroupCommunicators(const gmx::MDLogger &mdlog, /* Initially we set ntot to the number of PP cells, * This will be increased with PME cells when using Cartesian communicators. */ - copy_ivec(ddSetup.numDomains, ddRankSetup->ntot); + copy_ivec(ddRankSetup->numPPCells, ddRankSetup->ntot); ddRankSetup->bCartesianPP = (ddRankOrder == DdRankOrder::cartesian); ddRankSetup->bCartesianPP_PME = FALSE; - if (ddSetup.numPmeRanks > 0) + if (ddRankSetup->usePmeOnlyRanks) { /* Split the communicator into a PP and PME part */ split_communicator(mdlog, cr, ddRankOrder, ddSettings.useCartesianReorder, @@ -2329,7 +2330,7 @@ getDDSetup(const gmx::MDLogger &mdlog, if (options.numPmeRanks >= 0) { - ddSetup.numPmeRanks = options.numPmeRanks; + ddSetup.numPmeOnlyRanks = options.numPmeRanks; } else { @@ -2337,7 +2338,7 @@ getDDSetup(const gmx::MDLogger &mdlog, * don't use PME ranks. We check later if the DD grid is * compatible with the total number of ranks. */ - ddSetup.numPmeRanks = 0; + ddSetup.numPmeOnlyRanks = 0; } } else @@ -2366,7 +2367,7 @@ getDDSetup(const gmx::MDLogger &mdlog, "There is no domain decomposition for %d ranks that is compatible with the given box and a minimum cell size of %g nm\n" "%s\n" "Look in the log file for details on the domain decomposition", - cr->nnodes - ddSetup.numPmeRanks, + cr->nnodes - ddSetup.numPmeOnlyRanks, ddSetup.cellsizeLimit, buf); } @@ -2388,16 +2389,16 @@ getDDSetup(const gmx::MDLogger &mdlog, } const int numPPRanks = ddSetup.numDomains[XX]*ddSetup.numDomains[YY]*ddSetup.numDomains[ZZ]; - if (cr->nnodes - numPPRanks != ddSetup.numPmeRanks) + if (cr->nnodes - numPPRanks != ddSetup.numPmeOnlyRanks) { gmx_fatal_collective(FARGS, cr->mpi_comm_mysim, MASTER(cr), "The size of the domain decomposition grid (%d) does not match the number of ranks (%d). The total number of ranks is %d", - numPPRanks, cr->nnodes - ddSetup.numPmeRanks, cr->nnodes); + numPPRanks, cr->nnodes - ddSetup.numPmeOnlyRanks, cr->nnodes); } - if (ddSetup.numPmeRanks > numPPRanks) + if (ddSetup.numPmeOnlyRanks > numPPRanks) { gmx_fatal_collective(FARGS, cr->mpi_comm_mysim, MASTER(cr), - "The number of separate PME ranks (%d) is larger than the number of PP ranks (%d), this is not supported.", ddSetup.numPmeRanks, numPPRanks); + "The number of separate PME ranks (%d) is larger than the number of PP ranks (%d), this is not supported.", ddSetup.numPmeOnlyRanks, numPPRanks); } ddSetup.numDDDimensions = set_dd_dim(ddSetup.numDomains, ddSettings, @@ -2416,20 +2417,21 @@ getDDRankSetup(const gmx::MDLogger &mdlog, GMX_LOG(mdlog.info).appendTextFormatted( "Domain decomposition grid %d x %d x %d, separate PME ranks %d", ddSetup.numDomains[XX], ddSetup.numDomains[YY], ddSetup.numDomains[ZZ], - ddSetup.numPmeRanks); + ddSetup.numPmeOnlyRanks); DDRankSetup ddRankSetup; - ddRankSetup.numPPRanks = cr->nnodes - cr->npmenodes; + ddRankSetup.numPPRanks = cr->nnodes - ddSetup.numPmeOnlyRanks; copy_ivec(ddSetup.numDomains, ddRankSetup.numPPCells); - if (cr->npmenodes > 0) + ddRankSetup.usePmeOnlyRanks = (ddSetup.numPmeOnlyRanks > 0); + if (ddRankSetup.usePmeOnlyRanks) { - ddRankSetup.npmenodes = cr->npmenodes; + ddRankSetup.numRanksDoingPme = ddSetup.numPmeOnlyRanks; } else { - ddRankSetup.npmenodes = ddSetup.numDomains[XX]*ddSetup.numDomains[YY]*ddSetup.numDomains[ZZ]; + ddRankSetup.numRanksDoingPme = ddSetup.numDomains[XX]*ddSetup.numDomains[YY]*ddSetup.numDomains[ZZ]; } if (EEL_PME(ir.coulombtype) || EVDW_PME(ir.vdwtype)) @@ -2445,13 +2447,13 @@ getDDRankSetup(const gmx::MDLogger &mdlog, if (ddSetup.numDDDimensions >= 2 && ddSetup.ddDimensions[0] == XX && ddSetup.ddDimensions[1] == YY && - ddRankSetup.npmenodes > ddSetup.numDomains[XX] && - ddRankSetup.npmenodes % ddSetup.numDomains[XX] == 0 && + ddRankSetup.numRanksDoingPme > ddSetup.numDomains[XX] && + ddRankSetup.numRanksDoingPme % ddSetup.numDomains[XX] == 0 && getenv("GMX_PMEONEDD") == nullptr) { ddRankSetup.npmedecompdim = 2; ddRankSetup.npmenodes_x = ddSetup.numDomains[XX]; - ddRankSetup.npmenodes_y = ddRankSetup.npmenodes/ddRankSetup.npmenodes_x; + ddRankSetup.npmenodes_y = ddRankSetup.numRanksDoingPme/ddRankSetup.npmenodes_x; } else { @@ -2462,11 +2464,11 @@ getDDRankSetup(const gmx::MDLogger &mdlog, if (ddSetup.ddDimensions[0] == YY) { ddRankSetup.npmenodes_x = 1; - ddRankSetup.npmenodes_y = ddRankSetup.npmenodes; + ddRankSetup.npmenodes_y = ddRankSetup.numRanksDoingPme; } else { - ddRankSetup.npmenodes_x = ddRankSetup.npmenodes; + ddRankSetup.npmenodes_x = ddRankSetup.numRanksDoingPme; ddRankSetup.npmenodes_y = 1; } } @@ -2532,7 +2534,6 @@ static void set_dd_limits(const gmx::MDLogger &mdlog, comm->cgs_gl = gmx_mtop_global_cgs(mtop); /* Set the DD setup given by ddSetup */ - cr->npmenodes = ddSetup.numPmeRanks; copy_ivec(ddSetup.numDomains, dd->nc); dd->ndim = ddSetup.numDDDimensions; copy_ivec(ddSetup.ddDimensions, dd->dim); @@ -2937,7 +2938,7 @@ static void set_ddgrid_parameters(const gmx::MDLogger &mdlog, } else { - ddRankSetup.npmenodes = 0; + ddRankSetup.numRanksDoingPme = 0; if (dd->pme_nodeid >= 0) { gmx_fatal_collective(FARGS, dd->mpi_comm_all, DDMASTER(dd), @@ -3049,12 +3050,12 @@ gmx_domdec_t *init_domain_decomposition(const gmx::MDLogger &mdlog, DDSetup ddSetup = getDDSetup(mdlog, cr, options, ddSettings, systemInfo, mtop, ir, box, xGlobal, &ddbox); - cr->npmenodes = ddSetup.numPmeRanks; + cr->npmenodes = ddSetup.numPmeOnlyRanks; DDRankSetup ddRankSetup = getDDRankSetup(mdlog, cr, ddSetup, *ir); ivec ddCellIndex = { 0, 0, 0 }; - makeGroupCommunicators(mdlog, ddSettings, ddSetup, options.rankOrder, + makeGroupCommunicators(mdlog, ddSettings, options.rankOrder, &ddRankSetup, cr, ddCellIndex); gmx_domdec_t *dd = new gmx_domdec_t(*ir); diff --git a/src/gromacs/domdec/domdec.h b/src/gromacs/domdec/domdec.h index 2ef0ade747..52036a02bf 100644 --- a/src/gromacs/domdec/domdec.h +++ b/src/gromacs/domdec/domdec.h @@ -328,7 +328,7 @@ real comm_box_frac(const ivec dd_nc, real cutoff, const gmx_ddbox_t *ddbox); struct DDSetup { - int numPmeRanks = 0; + int numPmeOnlyRanks = 0; ivec numDomains = { 0, 0, 0 }; real cellsizeLimit = 0; int numDDDimensions = 0; diff --git a/src/gromacs/domdec/domdec_internal.h b/src/gromacs/domdec/domdec_internal.h index a56fc42fc9..0d67456f15 100644 --- a/src/gromacs/domdec/domdec_internal.h +++ b/src/gromacs/domdec/domdec_internal.h @@ -519,10 +519,11 @@ struct DDRankSetup ivec numPPCells = { 0, 0, 0 }; /* PME and Cartesian communicator stuff */ + bool usePmeOnlyRanks = false; /**< The number of decomposition dimensions for PME, 0: no PME */ int npmedecompdim = 0; /**< The number of ranks doing PME (PP/PME or only PME) */ - int npmenodes = 0; + int numRanksDoingPme = 0; /**< The number of PME ranks/domains along x */ int npmenodes_x = 0; /**< The number of PME ranks/domains along y */ diff --git a/src/gromacs/domdec/domdec_setup.cpp b/src/gromacs/domdec/domdec_setup.cpp index 2c5698b152..b72e8a11dd 100644 --- a/src/gromacs/domdec/domdec_setup.cpp +++ b/src/gromacs/domdec/domdec_setup.cpp @@ -746,7 +746,7 @@ dd_choose_grid(const gmx::MDLogger &mdlog, } else { - ddSetup.numPmeRanks = 0; + ddSetup.numPmeOnlyRanks = 0; } if (nnodes_div > 12) @@ -767,25 +767,26 @@ dd_choose_grid(const gmx::MDLogger &mdlog, /* Use PME nodes when the number of nodes is more than 16 */ if (cr->nnodes <= 18) { - ddSetup.numPmeRanks = 0; + ddSetup.numPmeOnlyRanks = 0; GMX_LOG(mdlog.info).appendTextFormatted( "Using %d separate PME ranks, as there are too few total\n" " ranks for efficient splitting", - cr->npmenodes); + ddSetup.numPmeOnlyRanks); } else { - ddSetup.numPmeRanks = guess_npme(mdlog, mtop, ir, box, cr->nnodes); + ddSetup.numPmeOnlyRanks = guess_npme(mdlog, mtop, ir, box, cr->nnodes); GMX_LOG(mdlog.info).appendTextFormatted( - "Using %d separate PME ranks, as guessed by mdrun", ddSetup.numPmeRanks); + "Using %d separate PME ranks, as guessed by mdrun", + ddSetup.numPmeOnlyRanks); } } else { /* We checked above that nPmeRanks is a valid number */ - ddSetup.numPmeRanks = numPmeRanksRequested; + ddSetup.numPmeOnlyRanks = numPmeRanksRequested; GMX_LOG(mdlog.info).appendTextFormatted( - "Using %d separate PME ranks", ddSetup.numPmeRanks); + "Using %d separate PME ranks", ddSetup.numPmeOnlyRanks); // TODO: there was a ", per user request" note here, but it's not correct anymore, // as with GPUs decision about nPmeRanks can be made in runner() as well. // Consider a single spot for setting nPmeRanks. @@ -793,7 +794,7 @@ dd_choose_grid(const gmx::MDLogger &mdlog, } ddSetup.cellsizeLimit = - optimize_ncells(mdlog, cr->nnodes, ddSetup.numPmeRanks, + optimize_ncells(mdlog, cr->nnodes, ddSetup.numPmeOnlyRanks, bDynLoadBal, dlb_scale, mtop, box, ddbox, ir, systemInfo, @@ -804,11 +805,11 @@ dd_choose_grid(const gmx::MDLogger &mdlog, gmx_bcast(sizeof(ddSetup.numDomains), ddSetup.numDomains, cr); if (EEL_PME(ir->coulombtype)) { - gmx_bcast(sizeof(ddSetup.numPmeRanks), &ddSetup.numPmeRanks, cr); + gmx_bcast(sizeof(ddSetup.numPmeOnlyRanks), &ddSetup.numPmeOnlyRanks, cr); } else { - ddSetup.numPmeRanks = 0; + ddSetup.numPmeOnlyRanks = 0; } return ddSetup; diff --git a/src/gromacs/domdec/partition.cpp b/src/gromacs/domdec/partition.cpp index a508fece1a..08462f1b26 100644 --- a/src/gromacs/domdec/partition.cpp +++ b/src/gromacs/domdec/partition.cpp @@ -1057,7 +1057,7 @@ static void print_dd_load_av(FILE *fplog, gmx_domdec_t *dd) char buf[STRLEN]; int numPpRanks = dd->nnodes; - int numPmeRanks = (dd->pme_nodeid >= 0) ? comm->ddRankSetup.npmenodes : 0; + int numPmeRanks = (comm->ddRankSetup.usePmeOnlyRanks ? comm->ddRankSetup.numRanksDoingPme : 0); int numRanks = numPpRanks + numPmeRanks; float lossFraction = 0; @@ -2921,7 +2921,7 @@ void dd_partition_system(FILE *fplog, * cost on the PME ranks, which will then surely result * in lower total performance. */ - if (cr->npmenodes > 0 && + if (comm->ddRankSetup.usePmeOnlyRanks && dd_pme_f_ratio(dd) > 1 - DD_PERF_LOSS_DLB_ON) { turnOnDlb = FALSE; diff --git a/src/gromacs/mdlib/gmx_omp_nthreads.cpp b/src/gromacs/mdlib/gmx_omp_nthreads.cpp index 5879179907..b2a0e46066 100644 --- a/src/gromacs/mdlib/gmx_omp_nthreads.cpp +++ b/src/gromacs/mdlib/gmx_omp_nthreads.cpp @@ -373,7 +373,7 @@ reportOpenmpSettings(const gmx::MDLogger &mdlog, } #if GMX_MPI - if (cr->nnodes + cr->npmenodes > 1) + if (cr->nnodes > 1) { /* Get the min and max thread counts over the MPI ranks */ int buf_in[4], buf_out[4]; diff --git a/src/gromacs/taskassignment/resourcedivision.cpp b/src/gromacs/taskassignment/resourcedivision.cpp index 90a40414c7..452a771616 100644 --- a/src/gromacs/taskassignment/resourcedivision.cpp +++ b/src/gromacs/taskassignment/resourcedivision.cpp @@ -563,7 +563,7 @@ void check_resource_division_efficiency(const gmx_hw_info_t *hwinfo, bool anyRankIsUsingGpus = willUsePhysicalGpu; /* Thread-MPI seems to have a bug with reduce on 1 node, so use a cond. */ - if (cr->nnodes + cr->npmenodes > 1) + if (cr->nnodes > 1) { int count[3], count_max[3]; -- 2.22.0