From 51cb3a47c7a2c62e5b5f1df5baeb5fda0d4162f9 Mon Sep 17 00:00:00 2001 From: Andrey Alekseenko Date: Wed, 14 Oct 2020 08:45:57 +0000 Subject: [PATCH] Revert "CUDA nightly fail: Dirty hack" This reverts commit c8b6d3b269730ad613b1810aa255e05395d7e8fe. --- src/gromacs/mdrun/md.cpp | 15 +++++++-------- src/gromacs/mdrun/mimic.cpp | 5 +++-- src/gromacs/mdrun/minimize.cpp | 5 +++-- src/gromacs/mdrun/rerun.cpp | 5 +++-- src/gromacs/mdrun/shellfc.cpp | 17 ++++++++++++++++- src/gromacs/mdrun/shellfc.h | 15 ++++++++++++--- src/gromacs/modularsimulator/forceelement.cpp | 4 +++- 7 files changed, 47 insertions(+), 19 deletions(-) diff --git a/src/gromacs/mdrun/md.cpp b/src/gromacs/mdrun/md.cpp index b95eb0b84e..f010e6a10a 100644 --- a/src/gromacs/mdrun/md.cpp +++ b/src/gromacs/mdrun/md.cpp @@ -300,9 +300,15 @@ void gmx::LegacySimulator::do_md() gstat = global_stat_init(ir); + const auto& simulationWork = runScheduleWork->simulationWork; + const bool useGpuForPme = simulationWork.useGpuPme; + const bool useGpuForNonbonded = simulationWork.useGpuNonbonded; + const bool useGpuForBufferOps = simulationWork.useGpuBufferOps; + const bool useGpuForUpdate = simulationWork.useGpuUpdate; + /* Check for polarizable models and flexible constraints */ shellfc = init_shell_flexcon(fplog, top_global, constr ? constr->numFlexibleConstraints() : 0, - ir->nstcalcenergy, DOMAINDECOMP(cr)); + ir->nstcalcenergy, DOMAINDECOMP(cr), useGpuForPme); { double io = compute_io(ir, top_global->natoms, *groups, energyOutput.numEnergyTerms(), 1); @@ -316,17 +322,10 @@ void gmx::LegacySimulator::do_md() std::unique_ptr stateInstance; t_state* state; - gmx_localtop_t top(top_global->ffparams); auto mdatoms = mdAtoms->mdatoms(); - const auto& simulationWork = runScheduleWork->simulationWork; - const bool useGpuForPme = simulationWork.useGpuPme; - const bool useGpuForNonbonded = simulationWork.useGpuNonbonded; - const bool useGpuForBufferOps = simulationWork.useGpuBufferOps; - const bool useGpuForUpdate = simulationWork.useGpuUpdate; - ForceBuffers f(fr->useMts, ((useGpuForNonbonded && useGpuForBufferOps) || useGpuForUpdate) ? PinningPolicy::PinnedIfSupported : PinningPolicy::CannotBePinned); diff --git a/src/gromacs/mdrun/mimic.cpp b/src/gromacs/mdrun/mimic.cpp index d232621817..e5ccee65e2 100644 --- a/src/gromacs/mdrun/mimic.cpp +++ b/src/gromacs/mdrun/mimic.cpp @@ -114,6 +114,7 @@ #include "gromacs/mdtypes/mdatom.h" #include "gromacs/mdtypes/mdrunoptions.h" #include "gromacs/mdtypes/observableshistory.h" +#include "gromacs/mdtypes/simulation_workload.h" #include "gromacs/mdtypes/state.h" #include "gromacs/mimic/communicator.h" #include "gromacs/mimic/utilities.h" @@ -131,7 +132,6 @@ #include "gromacs/utility/fatalerror.h" #include "gromacs/utility/logger.h" #include "gromacs/utility/real.h" -#include "gromacs/utility/smalloc.h" #include "legacysimulator.h" #include "replicaexchange.h" @@ -238,7 +238,8 @@ void gmx::LegacySimulator::do_mimic() /* Check for polarizable models and flexible constraints */ shellfc = init_shell_flexcon(fplog, top_global, constr ? constr->numFlexibleConstraints() : 0, - ir->nstcalcenergy, DOMAINDECOMP(cr)); + ir->nstcalcenergy, DOMAINDECOMP(cr), + runScheduleWork->simulationWork.useGpuPme); { double io = compute_io(ir, top_global->natoms, *groups, energyOutput.numEnergyTerms(), 1); diff --git a/src/gromacs/mdrun/minimize.cpp b/src/gromacs/mdrun/minimize.cpp index f03445b605..44276f96b1 100644 --- a/src/gromacs/mdrun/minimize.cpp +++ b/src/gromacs/mdrun/minimize.cpp @@ -399,8 +399,9 @@ static void init_em(FILE* fplog, { GMX_ASSERT(shellfc != nullptr, "With NM we always support shells"); - *shellfc = init_shell_flexcon(stdout, top_global, constr ? constr->numFlexibleConstraints() : 0, - ir->nstcalcenergy, DOMAINDECOMP(cr)); + *shellfc = + init_shell_flexcon(stdout, top_global, constr ? constr->numFlexibleConstraints() : 0, + ir->nstcalcenergy, DOMAINDECOMP(cr), thisRankHasDuty(cr, DUTY_PME)); } else { diff --git a/src/gromacs/mdrun/rerun.cpp b/src/gromacs/mdrun/rerun.cpp index 1f0f7d40ce..8b83a92f5a 100644 --- a/src/gromacs/mdrun/rerun.cpp +++ b/src/gromacs/mdrun/rerun.cpp @@ -112,6 +112,7 @@ #include "gromacs/mdtypes/mdatom.h" #include "gromacs/mdtypes/mdrunoptions.h" #include "gromacs/mdtypes/observableshistory.h" +#include "gromacs/mdtypes/simulation_workload.h" #include "gromacs/mdtypes/state.h" #include "gromacs/mimic/utilities.h" #include "gromacs/pbcutil/pbc.h" @@ -129,7 +130,6 @@ #include "gromacs/utility/fatalerror.h" #include "gromacs/utility/logger.h" #include "gromacs/utility/real.h" -#include "gromacs/utility/smalloc.h" #include "legacysimulator.h" #include "replicaexchange.h" @@ -294,7 +294,8 @@ void gmx::LegacySimulator::do_rerun() /* Check for polarizable models and flexible constraints */ shellfc = init_shell_flexcon(fplog, top_global, constr ? constr->numFlexibleConstraints() : 0, - ir->nstcalcenergy, DOMAINDECOMP(cr)); + ir->nstcalcenergy, DOMAINDECOMP(cr), + runScheduleWork->simulationWork.useGpuPme); { double io = compute_io(ir, top_global->natoms, *groups, energyOutput.numEnergyTerms(), 1); diff --git a/src/gromacs/mdrun/shellfc.cpp b/src/gromacs/mdrun/shellfc.cpp index a2c8b86f08..ccc3c3e3ca 100644 --- a/src/gromacs/mdrun/shellfc.cpp +++ b/src/gromacs/mdrun/shellfc.cpp @@ -253,7 +253,12 @@ static void predict_shells(FILE* fplog, } } -gmx_shellfc_t* init_shell_flexcon(FILE* fplog, const gmx_mtop_t* mtop, int nflexcon, int nstcalcenergy, bool usingDomainDecomposition) +gmx_shellfc_t* init_shell_flexcon(FILE* fplog, + const gmx_mtop_t* mtop, + int nflexcon, + int nstcalcenergy, + bool usingDomainDecomposition, + bool usingPmeOnGpu) { gmx_shellfc_t* shfc; @@ -535,6 +540,16 @@ gmx_shellfc_t* init_shell_flexcon(FILE* fplog, const gmx_mtop_t* mtop, int nflex } } + /* shfc->x is used as a coordinate buffer for the sim_util's `do_force` function, and + * when using PME it must be pinned. */ + if (usingPmeOnGpu) + { + for (i = 0; i < 2; i++) + { + changePinningPolicy(&shfc->x[i], gmx::PinningPolicy::PinnedIfSupported); + } + } + return shfc; } diff --git a/src/gromacs/mdrun/shellfc.h b/src/gromacs/mdrun/shellfc.h index d14646b646..f5145cfb8d 100644 --- a/src/gromacs/mdrun/shellfc.h +++ b/src/gromacs/mdrun/shellfc.h @@ -72,14 +72,23 @@ class MdrunScheduleWorkload; class VirtualSitesHandler; } // namespace gmx -/* Initialization function, also predicts the initial shell postions. - * Returns a pointer to an initialized shellfc object. +/*! \brief Initialization function, also predicts the initial shell positions. + * + * \param fplog Pointer to the log stream. Can be set to \c nullptr to disable verbose log. + * \param mtop Pointer to a global system topology object. + * \param nflexcon Number of flexible constraints. + * \param nstcalcenergy How often are energies calculated. Must be provided for sanity check. + * \param usingDomainDecomposition Whether domain decomposition is used. Must be provided for sanity check. + * \param usingPmeOnGpu Set to true if GPU will be used for PME calculations. Necessary for proper buffer initialization. + * + * \returns a pointer to an initialized \c shellfc object. */ gmx_shellfc_t* init_shell_flexcon(FILE* fplog, const gmx_mtop_t* mtop, int nflexcon, int nstcalcenergy, - bool usingDomainDecomposition); + bool usingDomainDecomposition, + bool usingPmeOnGpu); /* Optimize shell positions */ void relax_shell_flexcon(FILE* log, diff --git a/src/gromacs/modularsimulator/forceelement.cpp b/src/gromacs/modularsimulator/forceelement.cpp index 4d43376706..3aae7e3756 100644 --- a/src/gromacs/modularsimulator/forceelement.cpp +++ b/src/gromacs/modularsimulator/forceelement.cpp @@ -54,6 +54,7 @@ #include "gromacs/mdtypes/inputrec.h" #include "gromacs/mdtypes/mdatom.h" #include "gromacs/mdtypes/mdrunoptions.h" +#include "gromacs/mdtypes/simulation_workload.h" #include "gromacs/pbcutil/pbc.h" #include "energydata.h" @@ -93,7 +94,8 @@ ForceElement::ForceElement(StatePropagatorData* statePropagatorData, globalTopology, constr ? constr->numFlexibleConstraints() : 0, inputrec->nstcalcenergy, - DOMAINDECOMP(cr))), + DOMAINDECOMP(cr), + runScheduleWork->simulationWork.useGpuPme)), doShellFC_(shellfc_ != nullptr), nextNSStep_(-1), nextEnergyCalculationStep_(-1), -- 2.22.0