From 5ec4eb24d033403f4a051365f1d117f5f6cb6a64 Mon Sep 17 00:00:00 2001 From: Aleksei Iupinov Date: Thu, 26 Oct 2017 17:04:34 +0200 Subject: [PATCH] Simplify PME data handling in runner Differing ownership of the PME data for PME-only and other ranks is now hidden behind a reference. gmx_pme_init() now returns a pointer to the allocated structure. Change-Id: Ia9c5117a0db43a6564298dd621cf9254f0423acf --- src/gromacs/ewald/pme.cpp | 40 +++++++++++++++---------------- src/gromacs/ewald/pme.h | 26 ++++++++++---------- src/gromacs/ewald/tests/pmetestcommon.cpp | 5 ++-- src/programs/mdrun/runner.cpp | 30 +++++++++++------------ 4 files changed, 50 insertions(+), 51 deletions(-) diff --git a/src/gromacs/ewald/pme.cpp b/src/gromacs/ewald/pme.cpp index 7698e90eef..3fb2517dbd 100644 --- a/src/gromacs/ewald/pme.cpp +++ b/src/gromacs/ewald/pme.cpp @@ -505,22 +505,21 @@ static int div_round_up(int enumerator, int denominator) return (enumerator + denominator - 1)/denominator; } -void gmx_pme_init(struct gmx_pme_t **pmedata, - t_commrec *cr, - int nnodes_major, - int nnodes_minor, - const t_inputrec *ir, - int homenr, - gmx_bool bFreeEnergy_q, - gmx_bool bFreeEnergy_lj, - gmx_bool bReproducible, - real ewaldcoeff_q, - real ewaldcoeff_lj, - int nthread, - PmeRunMode runMode, - pme_gpu_t *pmeGPU, - gmx_device_info_t *gpuInfo, - const gmx::MDLogger &mdlog) +gmx_pme_t *gmx_pme_init(const t_commrec *cr, + int nnodes_major, + int nnodes_minor, + const t_inputrec *ir, + int homenr, + gmx_bool bFreeEnergy_q, + gmx_bool bFreeEnergy_lj, + gmx_bool bReproducible, + real ewaldcoeff_q, + real ewaldcoeff_lj, + int nthread, + PmeRunMode runMode, + pme_gpu_t *pmeGPU, + gmx_device_info_t *gpuInfo, + const gmx::MDLogger &mdlog) { int use_threads, sum_use_threads, i; ivec ndata; @@ -859,7 +858,7 @@ void gmx_pme_init(struct gmx_pme_t **pmedata, pme_init_all_work(&pme->solve_work, pme->nthread, pme->nkx); // no exception was thrown during the init, so we hand over the PME structure handle - *pmedata = pme.release(); + return pme.release(); } void gmx_pme_reinit(struct gmx_pme_t **pmedata, @@ -902,9 +901,10 @@ void gmx_pme_reinit(struct gmx_pme_t **pmedata, // This is reinit which is currently only changing grid size/coefficients, // so we don't expect the actual logging. // TODO: when PME is an object, it should take reference to mdlog on construction and save it. - gmx_pme_init(pmedata, cr, pme_src->nnodes_major, pme_src->nnodes_minor, - &irc, homenr, pme_src->bFEP_q, pme_src->bFEP_lj, FALSE, ewaldcoeff_q, ewaldcoeff_lj, - pme_src->nthread, pme_src->runMode, pme_src->gpu, nullptr, dummyLogger); + GMX_ASSERT(pmedata, "Invalid PME pointer"); + *pmedata = gmx_pme_init(cr, pme_src->nnodes_major, pme_src->nnodes_minor, + &irc, homenr, pme_src->bFEP_q, pme_src->bFEP_lj, FALSE, ewaldcoeff_q, ewaldcoeff_lj, + pme_src->nthread, pme_src->runMode, pme_src->gpu, nullptr, dummyLogger); //TODO this is mostly passing around current values } GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR; diff --git a/src/gromacs/ewald/pme.h b/src/gromacs/ewald/pme.h index 7f03067052..af9fd7c77b 100644 --- a/src/gromacs/ewald/pme.h +++ b/src/gromacs/ewald/pme.h @@ -61,6 +61,7 @@ struct t_nrnb; struct pme_gpu_t; struct gmx_wallclock_gpu_pme_t; struct gmx_device_info_t; +struct gmx_pme_t; namespace gmx { @@ -92,21 +93,22 @@ bool gmx_pme_check_restrictions(int pme_order, bool useThreads, bool errorsAreFatal); -/*! \brief Initialize \p pmedata +/*! \brief Construct PME data * * \throws gmx::InconsistentInputError if input grid sizes/PME order are inconsistent. + * \returns Pointer to newly allocated and initialized PME data. */ -void gmx_pme_init(struct gmx_pme_t **pmedata, struct t_commrec *cr, - int nnodes_major, int nnodes_minor, - const t_inputrec *ir, int homenr, - gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj, - gmx_bool bReproducible, - real ewaldcoeff_q, real ewaldcoeff_lj, - int nthread, - PmeRunMode runMode, - pme_gpu_t *pmeGPU, - gmx_device_info_t *gpuInfo, - const gmx::MDLogger &mdlog); +gmx_pme_t *gmx_pme_init(const t_commrec *cr, + int nnodes_major, int nnodes_minor, + const t_inputrec *ir, int homenr, + gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj, + gmx_bool bReproducible, + real ewaldcoeff_q, real ewaldcoeff_lj, + int nthread, + PmeRunMode runMode, + pme_gpu_t *pmeGPU, + gmx_device_info_t *gpuInfo, + const gmx::MDLogger &mdlog); /*! \brief Destroys the PME data structure.*/ void gmx_pme_destroy(gmx_pme_t *pme); diff --git a/src/gromacs/ewald/tests/pmetestcommon.cpp b/src/gromacs/ewald/tests/pmetestcommon.cpp index dc3af72015..0cd03e2da2 100644 --- a/src/gromacs/ewald/tests/pmetestcommon.cpp +++ b/src/gromacs/ewald/tests/pmetestcommon.cpp @@ -102,12 +102,11 @@ static PmeSafePointer pmeInitInternal(const t_inputrec *inputRec, real ewaldCoeff_lj = 1.0f ) { - gmx_pme_t *pmeDataRaw = nullptr; const MDLogger dummyLogger; const auto runMode = (mode == CodePath::CPU) ? PmeRunMode::CPU : PmeRunMode::GPU; t_commrec dummyCommrec = {0}; - gmx_pme_init(&pmeDataRaw, &dummyCommrec, 1, 1, inputRec, atomCount, false, false, true, - ewaldCoeff_q, ewaldCoeff_lj, 1, runMode, nullptr, gpuInfo, dummyLogger); + gmx_pme_t *pmeDataRaw = gmx_pme_init(&dummyCommrec, 1, 1, inputRec, atomCount, false, false, true, + ewaldCoeff_q, ewaldCoeff_lj, 1, runMode, nullptr, gpuInfo, dummyLogger); PmeSafePointer pme(pmeDataRaw); // taking ownership // TODO get rid of this with proper matrix type diff --git a/src/programs/mdrun/runner.cpp b/src/programs/mdrun/runner.cpp index edd1f1b3f0..7d617b97fd 100644 --- a/src/programs/mdrun/runner.cpp +++ b/src/programs/mdrun/runner.cpp @@ -441,7 +441,6 @@ int Mdrunner::mdrunner() t_fcdata *fcd = nullptr; real ewaldcoeff_q = 0; real ewaldcoeff_lj = 0; - struct gmx_pme_t **pmedata = nullptr; gmx_vsite_t *vsite = nullptr; gmx_constr_t constr; int nChargePerturbed = -1, nTypePerturbed = 0; @@ -1013,11 +1012,6 @@ int Mdrunner::mdrunner() { ewaldcoeff_q = fr->ic->ewaldcoeff_q; ewaldcoeff_lj = fr->ic->ewaldcoeff_lj; - pmedata = &fr->pmedata; - } - else - { - pmedata = nullptr; } } else @@ -1028,9 +1022,13 @@ int Mdrunner::mdrunner() ewaldcoeff_q = calc_ewaldcoeff_q(inputrec->rcoulomb, inputrec->ewald_rtol); ewaldcoeff_lj = calc_ewaldcoeff_lj(inputrec->rvdw, inputrec->ewald_rtol_lj); - snew(pmedata, 1); } + gmx_pme_t *sepPmeData = nullptr; + // This reference hides the fact that PME data is owned by runner on PME-only ranks and by forcerec on other ranks + GMX_ASSERT(thisRankHasDuty(cr, DUTY_PP) == (fr != nullptr), "Double-checking that only PME-only ranks have no forcerec"); + gmx_pme_t * &pmedata = fr ? fr->pmedata : sepPmeData; + if (hw_opt.thread_affinity != threadaffOFF) { /* Before setting affinity, check whether the affinity has changed @@ -1080,12 +1078,12 @@ int Mdrunner::mdrunner() try { gmx_device_info_t *pmeGpuInfo = nullptr; - gmx_pme_init(pmedata, cr, npme_major, npme_minor, inputrec, - mtop ? mtop->natoms : 0, nChargePerturbed, nTypePerturbed, - mdrunOptions.reproducible, - ewaldcoeff_q, ewaldcoeff_lj, - nthreads_pme, - pmeRunMode, nullptr, pmeGpuInfo, mdlog); + pmedata = gmx_pme_init(cr, npme_major, npme_minor, inputrec, + mtop ? mtop->natoms : 0, nChargePerturbed, nTypePerturbed, + mdrunOptions.reproducible, + ewaldcoeff_q, ewaldcoeff_lj, + nthreads_pme, + pmeRunMode, nullptr, pmeGpuInfo, mdlog); } GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR; } @@ -1172,7 +1170,7 @@ int Mdrunner::mdrunner() GMX_RELEASE_ASSERT(pmedata, "pmedata was NULL while cr->duty was not DUTY_PP"); /* do PME only */ walltime_accounting = walltime_accounting_init(gmx_omp_nthreads_get(emntPME)); - gmx_pmeonly(*pmedata, cr, nrnb, wcycle, walltime_accounting, inputrec, pmeRunMode); + gmx_pmeonly(pmedata, cr, nrnb, wcycle, walltime_accounting, inputrec, pmeRunMode); } wallcycle_stop(wcycle, ewcRUN); @@ -1183,13 +1181,13 @@ int Mdrunner::mdrunner() finish_run(fplog, mdlog, cr, inputrec, nrnb, wcycle, walltime_accounting, fr ? fr->nbv : nullptr, - fr ? fr->pmedata : nullptr, + pmedata, EI_DYNAMICS(inputrec->eI) && !MULTISIM(cr)); // Free PME data if (pmedata) { - gmx_pme_destroy(*pmedata); // TODO: pmedata is always a single element list, refactor + gmx_pme_destroy(pmedata); pmedata = nullptr; } -- 2.11.4.GIT