From a5d4bf712ef319814efd90ccd50c492f0b9c07fc Mon Sep 17 00:00:00 2001 From: Artem Zhmurov Date: Sat, 20 Feb 2021 22:48:03 +0000 Subject: [PATCH] Set temperature scaling groups upon construction of integrator The temperature scaling parameters do not change, so can be set upon construction of the integrator object. Closes #3794 --- src/gromacs/mdlib/leapfrog_gpu.cu | 34 +++++++++---------- src/gromacs/mdlib/leapfrog_gpu.h | 19 +++++------ src/gromacs/mdlib/leapfrog_gpu_sycl.cpp | 30 ++++++++-------- .../mdlib/tests/leapfrogtestrunners_gpu.cpp | 10 +++--- src/gromacs/mdlib/update_constrain_gpu.h | 24 ++++++------- .../mdlib/update_constrain_gpu_impl.cpp | 6 ++-- .../mdlib/update_constrain_gpu_impl.cu | 18 +++++----- src/gromacs/mdlib/update_constrain_gpu_impl.h | 26 +++++++------- src/gromacs/mdrun/md.cpp | 4 +-- 9 files changed, 81 insertions(+), 90 deletions(-) diff --git a/src/gromacs/mdlib/leapfrog_gpu.cu b/src/gromacs/mdlib/leapfrog_gpu.cu index d47619109e..2f5b589870 100644 --- a/src/gromacs/mdlib/leapfrog_gpu.cu +++ b/src/gromacs/mdlib/leapfrog_gpu.cu @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2019,2020, by the GROMACS development team, led by + * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -307,9 +307,12 @@ void LeapFrogGpu::integrate(const DeviceBuffer d_x, return; } -LeapFrogGpu::LeapFrogGpu(const DeviceContext& deviceContext, const DeviceStream& deviceStream) : +LeapFrogGpu::LeapFrogGpu(const DeviceContext& deviceContext, + const DeviceStream& deviceStream, + const int numTempScaleValues) : deviceContext_(deviceContext), - deviceStream_(deviceStream) + deviceStream_(deviceStream), + numTempScaleValues_(numTempScaleValues) { numAtoms_ = 0; @@ -319,6 +322,14 @@ LeapFrogGpu::LeapFrogGpu(const DeviceContext& deviceContext, const DeviceStream& kernelLaunchConfig_.blockSize[1] = 1; kernelLaunchConfig_.blockSize[2] = 1; kernelLaunchConfig_.sharedMemorySize = 0; + + // If the temperature coupling is enabled, we need to make space for scaling factors + if (numTempScaleValues_ > 0) + { + h_lambdas_.resize(numTempScaleValues_); + reallocateDeviceBuffer( + &d_lambdas_, numTempScaleValues_, &numLambdas_, &numLambdasAlloc_, deviceContext_); + } } LeapFrogGpu::~LeapFrogGpu() @@ -326,37 +337,24 @@ LeapFrogGpu::~LeapFrogGpu() freeDeviceBuffer(&d_inverseMasses_); } -void LeapFrogGpu::set(const int numAtoms, - const real* inverseMasses, - const int numTempScaleValues, - const unsigned short* tempScaleGroups) +void LeapFrogGpu::set(const int numAtoms, const real* inverseMasses, const unsigned short* tempScaleGroups) { numAtoms_ = numAtoms; kernelLaunchConfig_.gridSize[0] = (numAtoms_ + c_threadsPerBlock - 1) / c_threadsPerBlock; - numTempScaleValues_ = numTempScaleValues; - reallocateDeviceBuffer( &d_inverseMasses_, numAtoms_, &numInverseMasses_, &numInverseMassesAlloc_, deviceContext_); copyToDeviceBuffer( &d_inverseMasses_, (float*)inverseMasses, 0, numAtoms_, deviceStream_, GpuApiCallBehavior::Sync, nullptr); // Temperature scale group map only used if there are more then one group - if (numTempScaleValues > 1) + if (numTempScaleValues_ > 1) { reallocateDeviceBuffer( &d_tempScaleGroups_, numAtoms_, &numTempScaleGroups_, &numTempScaleGroupsAlloc_, deviceContext_); copyToDeviceBuffer( &d_tempScaleGroups_, tempScaleGroups, 0, numAtoms_, deviceStream_, GpuApiCallBehavior::Sync, nullptr); } - - // If the temperature coupling is enabled, we need to make space for scaling factors - if (numTempScaleValues_ > 0) - { - h_lambdas_.resize(numTempScaleValues); - reallocateDeviceBuffer( - &d_lambdas_, numTempScaleValues_, &numLambdas_, &numLambdasAlloc_, deviceContext_); - } } } // namespace gmx diff --git a/src/gromacs/mdlib/leapfrog_gpu.h b/src/gromacs/mdlib/leapfrog_gpu.h index e802d4c8e3..108259f1c6 100644 --- a/src/gromacs/mdlib/leapfrog_gpu.h +++ b/src/gromacs/mdlib/leapfrog_gpu.h @@ -101,10 +101,11 @@ class LeapFrogGpu public: /*! \brief Constructor. * - * \param[in] deviceContext Device context (dummy in CUDA). - * \param[in] deviceStream Device stream to use. + * \param[in] deviceContext Device context (dummy in CUDA). + * \param[in] deviceStream Device stream to use. + * \param[in] numTempScaleValues Number of temperature scale groups. */ - LeapFrogGpu(const DeviceContext& deviceContext, const DeviceStream& deviceStream); + LeapFrogGpu(const DeviceContext& deviceContext, const DeviceStream& deviceStream, int numTempScaleValues); ~LeapFrogGpu(); /*! \brief Integrate @@ -140,15 +141,11 @@ public: * and temperature coupling groups. Copies inverse masses and temperature coupling groups * to the GPU. * - * \param[in] numAtoms Number of atoms in the system. - * \param[in] inverseMasses Inverse masses of atoms. - * \param[in] numTempScaleValues Number of temperature scale groups. - * \param[in] tempScaleGroups Maps the atom index to temperature scale value. + * \param[in] numAtoms Number of atoms in the system. + * \param[in] inverseMasses Inverse masses of atoms. + * \param[in] tempScaleGroups Maps the atom index to temperature scale value. */ - void set(const int numAtoms, - const real* inverseMasses, - int numTempScaleValues, - const unsigned short* tempScaleGroups); + void set(const int numAtoms, const real* inverseMasses, const unsigned short* tempScaleGroups); /*! \brief Class with hardware-specific interfaces and implementations.*/ class Impl; diff --git a/src/gromacs/mdlib/leapfrog_gpu_sycl.cpp b/src/gromacs/mdlib/leapfrog_gpu_sycl.cpp index d82c9c12eb..9afb0320fb 100644 --- a/src/gromacs/mdlib/leapfrog_gpu_sycl.cpp +++ b/src/gromacs/mdlib/leapfrog_gpu_sycl.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2020, by the GROMACS development team, led by + * Copyright (c) 2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -273,11 +273,20 @@ void LeapFrogGpu::integrate(DeviceBuffer d_x, prVelocityScalingMatrixDiagonal_); } -LeapFrogGpu::LeapFrogGpu(const DeviceContext& deviceContext, const DeviceStream& deviceStream) : +LeapFrogGpu::LeapFrogGpu(const DeviceContext& deviceContext, + const DeviceStream& deviceStream, + const int numTempScaleValues) : deviceContext_(deviceContext), deviceStream_(deviceStream), - numAtoms_(0) + numAtoms_(0), + numTempScaleValues_(numTempScaleValues) { + // If the temperature coupling is enabled, we need to make space for scaling factors + if (numTempScaleValues_ > 0) + { + reallocateDeviceBuffer( + &d_lambdas_, numTempScaleValues_, &numLambdas_, &numLambdasAlloc_, deviceContext_); + } } LeapFrogGpu::~LeapFrogGpu() @@ -285,13 +294,9 @@ LeapFrogGpu::~LeapFrogGpu() freeDeviceBuffer(&d_inverseMasses_); } -void LeapFrogGpu::set(const int numAtoms, - const real* inverseMasses, - const int numTempScaleValues, - const unsigned short* tempScaleGroups) +void LeapFrogGpu::set(const int numAtoms, const real* inverseMasses, const unsigned short* tempScaleGroups) { - numAtoms_ = numAtoms; - numTempScaleValues_ = numTempScaleValues; + numAtoms_ = numAtoms; reallocateDeviceBuffer( &d_inverseMasses_, numAtoms_, &numInverseMasses_, &numInverseMassesAlloc_, deviceContext_); @@ -306,13 +311,6 @@ void LeapFrogGpu::set(const int numAtoms, copyToDeviceBuffer( &d_tempScaleGroups_, tempScaleGroups, 0, numAtoms_, deviceStream_, GpuApiCallBehavior::Sync, nullptr); } - - // If the temperature coupling is enabled, we need to make space for scaling factors - if (numTempScaleValues_ > 0) - { - reallocateDeviceBuffer( - &d_lambdas_, numTempScaleValues_, &numLambdas_, &numLambdasAlloc_, deviceContext_); - } } } // namespace gmx diff --git a/src/gromacs/mdlib/tests/leapfrogtestrunners_gpu.cpp b/src/gromacs/mdlib/tests/leapfrogtestrunners_gpu.cpp index f88ae219ce..b85ed94b52 100644 --- a/src/gromacs/mdlib/tests/leapfrogtestrunners_gpu.cpp +++ b/src/gromacs/mdlib/tests/leapfrogtestrunners_gpu.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2019,2020, by the GROMACS development team, led by + * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -95,12 +95,10 @@ void LeapFrogDeviceTestRunner::integrate(LeapFrogTestData* testData, int numStep copyToDeviceBuffer(&d_v, h_v, 0, numAtoms, deviceStream, GpuApiCallBehavior::Sync, nullptr); copyToDeviceBuffer(&d_f, h_f, 0, numAtoms, deviceStream, GpuApiCallBehavior::Sync, nullptr); - auto integrator = std::make_unique(deviceContext, deviceStream); + auto integrator = + std::make_unique(deviceContext, deviceStream, testData->numTCoupleGroups_); - integrator->set(testData->numAtoms_, - testData->inverseMasses_.data(), - testData->numTCoupleGroups_, - testData->mdAtoms_.cTC); + integrator->set(testData->numAtoms_, testData->inverseMasses_.data(), testData->mdAtoms_.cTC); bool doTempCouple = testData->numTCoupleGroups_ > 0; for (int step = 0; step < numSteps; step++) diff --git a/src/gromacs/mdlib/update_constrain_gpu.h b/src/gromacs/mdlib/update_constrain_gpu.h index 3e77d8d3ce..ce23e1dfa0 100644 --- a/src/gromacs/mdlib/update_constrain_gpu.h +++ b/src/gromacs/mdlib/update_constrain_gpu.h @@ -76,18 +76,20 @@ public: * any) consumers of the updated coordinates. The \p xUpdatedOnDevice also can not be a nullptr * because the markEvent(...) method is called unconditionally. * - * \param[in] ir Input record data: LINCS takes number of iterations and order of - * projection from it. - * \param[in] mtop Topology of the system: SETTLE gets the masses for O and H atoms - * and target O-H and H-H distances from this object. - * \param[in] deviceContext GPU device context. - * \param[in] deviceStream GPU stream to use. - * \param[in] xUpdatedOnDevice The event synchronizer to use to mark that update is done - * on the GPU. - * \param[in] wcycle The wallclock counter + * \param[in] ir Input record data: LINCS takes number of iterations and order of + * projection from it. + * \param[in] mtop Topology of the system: SETTLE gets the masses for O and H atoms + * and target O-H and H-H distances from this object. + * \param[in] numTempScaleValues Number of temperature scaling groups. Zero for no temperature scaling. + * \param[in] deviceContext GPU device context. + * \param[in] deviceStream GPU stream to use. + * \param[in] xUpdatedOnDevice The event synchronizer to use to mark that update is done + * on the GPU. + * \param[in] wcycle The wallclock counter */ UpdateConstrainGpu(const t_inputrec& ir, const gmx_mtop_t& mtop, + int numTempScaleValues, const DeviceContext& deviceContext, const DeviceStream& deviceStream, GpuEventSynchronizer* xUpdatedOnDevice, @@ -147,14 +149,12 @@ public: * \param[in] d_f Device buffer with forces. * \param[in] idef System topology * \param[in] md Atoms data. - * \param[in] numTempScaleValues Number of temperature scaling groups. Zero for no temperature scaling. */ void set(DeviceBuffer d_x, DeviceBuffer d_v, DeviceBuffer d_f, const InteractionDefinitions& idef, - const t_mdatoms& md, - int numTempScaleValues); + const t_mdatoms& md); /*! \brief * Update PBC data. diff --git a/src/gromacs/mdlib/update_constrain_gpu_impl.cpp b/src/gromacs/mdlib/update_constrain_gpu_impl.cpp index dc2b0421c3..b0f3872860 100644 --- a/src/gromacs/mdlib/update_constrain_gpu_impl.cpp +++ b/src/gromacs/mdlib/update_constrain_gpu_impl.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2019,2020, by the GROMACS development team, led by + * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -58,6 +58,7 @@ class UpdateConstrainGpu::Impl UpdateConstrainGpu::UpdateConstrainGpu(const t_inputrec& /* ir */, const gmx_mtop_t& /* mtop */, + const int /* numTempScaleValues */, const DeviceContext& /* deviceContext */, const DeviceStream& /* deviceStream */, GpuEventSynchronizer* /* xUpdatedOnDevice */, @@ -101,8 +102,7 @@ void UpdateConstrainGpu::set(DeviceBuffer /* d_x */, DeviceBuffer /* d_v */, const DeviceBuffer /* d_f */, const InteractionDefinitions& /* idef */, - const t_mdatoms& /* md */, - const int /* numTempScaleValues */) + const t_mdatoms& /* md */) { GMX_ASSERT(!impl_, "A CPU stub for UpdateConstrain was called instead of the correct implementation."); diff --git a/src/gromacs/mdlib/update_constrain_gpu_impl.cu b/src/gromacs/mdlib/update_constrain_gpu_impl.cu index ba27a8a58c..f03ab778c3 100644 --- a/src/gromacs/mdlib/update_constrain_gpu_impl.cu +++ b/src/gromacs/mdlib/update_constrain_gpu_impl.cu @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2019,2020, by the GROMACS development team, led by + * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -215,6 +215,7 @@ void UpdateConstrainGpu::Impl::scaleVelocities(const matrix scalingMatrix) UpdateConstrainGpu::Impl::Impl(const t_inputrec& ir, const gmx_mtop_t& mtop, + const int numTempScaleValues, const DeviceContext& deviceContext, const DeviceStream& deviceStream, GpuEventSynchronizer* xUpdatedOnDevice, @@ -227,7 +228,7 @@ UpdateConstrainGpu::Impl::Impl(const t_inputrec& ir, GMX_ASSERT(xUpdatedOnDevice != nullptr, "The event synchronizer can not be nullptr."); - integrator_ = std::make_unique(deviceContext_, deviceStream_); + integrator_ = std::make_unique(deviceContext_, deviceStream_, numTempScaleValues); lincsGpu_ = std::make_unique(ir.nLincsIter, ir.nProjOrder, deviceContext_, deviceStream_); settleGpu_ = std::make_unique(mtop, deviceContext_, deviceStream_); @@ -243,8 +244,7 @@ void UpdateConstrainGpu::Impl::set(DeviceBuffer d_x, DeviceBuffer d_v, const DeviceBuffer d_f, const InteractionDefinitions& idef, - const t_mdatoms& md, - const int numTempScaleValues) + const t_mdatoms& md) { // TODO wallcycle wallcycle_start_nocount(wcycle_, ewcLAUNCH_GPU); @@ -266,7 +266,7 @@ void UpdateConstrainGpu::Impl::set(DeviceBuffer d_x, &d_inverseMasses_, numAtoms_, &numInverseMasses_, &numInverseMassesAlloc_, deviceContext_); // Integrator should also update something, but it does not even have a method yet - integrator_->set(numAtoms_, md.invmass, numTempScaleValues, md.cTC); + integrator_->set(numAtoms_, md.invmass, md.cTC); lincsGpu_->set(idef, numAtoms_, md.invmass); settleGpu_->set(idef); @@ -290,11 +290,12 @@ GpuEventSynchronizer* UpdateConstrainGpu::Impl::getCoordinatesReadySync() UpdateConstrainGpu::UpdateConstrainGpu(const t_inputrec& ir, const gmx_mtop_t& mtop, + const int numTempScaleValues, const DeviceContext& deviceContext, const DeviceStream& deviceStream, GpuEventSynchronizer* xUpdatedOnDevice, gmx_wallcycle* wcycle) : - impl_(new Impl(ir, mtop, deviceContext, deviceStream, xUpdatedOnDevice, wcycle)) + impl_(new Impl(ir, mtop, numTempScaleValues, deviceContext, deviceStream, xUpdatedOnDevice, wcycle)) { } @@ -337,10 +338,9 @@ void UpdateConstrainGpu::set(DeviceBuffer d_x, DeviceBuffer d_v, const DeviceBuffer d_f, const InteractionDefinitions& idef, - const t_mdatoms& md, - const int numTempScaleValues) + const t_mdatoms& md) { - impl_->set(d_x, d_v, d_f, idef, md, numTempScaleValues); + impl_->set(d_x, d_v, d_f, idef, md); } void UpdateConstrainGpu::setPbc(const PbcType pbcType, const matrix box) diff --git a/src/gromacs/mdlib/update_constrain_gpu_impl.h b/src/gromacs/mdlib/update_constrain_gpu_impl.h index 7453a98105..8e101b8cd0 100644 --- a/src/gromacs/mdlib/update_constrain_gpu_impl.h +++ b/src/gromacs/mdlib/update_constrain_gpu_impl.h @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2019,2020, by the GROMACS development team, led by + * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -71,18 +71,20 @@ public: * any) consumers of the updated coordinates. The \p xUpdatedOnDevice also can not be a nullptr * because the markEvent(...) method is called unconditionally. * - * \param[in] ir Input record data: LINCS takes number of iterations and order of - * projection from it. - * \param[in] mtop Topology of the system: SETTLE gets the masses for O and H atoms - * and target O-H and H-H distances from this object. - * \param[in] deviceContext GPU device context. - * \param[in] deviceStream GPU stream to use. - * \param[in] xUpdatedOnDevice The event synchronizer to use to mark that - * update is done on the GPU. - * \param[in] wcycle The wallclock counter + * \param[in] ir Input record data: LINCS takes number of iterations and order of + * projection from it. + * \param[in] mtop Topology of the system: SETTLE gets the masses for O and H atoms + * and target O-H and H-H distances from this object. + * \param[in] numTempScaleValues Number of temperature scaling groups. Set zero for no temperature coupling. + * \param[in] deviceContext GPU device context. + * \param[in] deviceStream GPU stream to use. + * \param[in] xUpdatedOnDevice The event synchronizer to use to mark that + * update is done on the GPU. + * \param[in] wcycle The wallclock counter */ Impl(const t_inputrec& ir, const gmx_mtop_t& mtop, + int numTempScaleValues, const DeviceContext& deviceContext, const DeviceStream& deviceStream, GpuEventSynchronizer* xUpdatedOnDevice, @@ -147,14 +149,12 @@ public: * \param[in] d_f Device buffer with forces. * \param[in] idef System topology * \param[in] md Atoms data. - * \param[in] numTempScaleValues Number of temperature scaling groups. Set zero for no temperature coupling. */ void set(DeviceBuffer d_x, DeviceBuffer d_v, const DeviceBuffer d_f, const InteractionDefinitions& idef, - const t_mdatoms& md, - const int numTempScaleValues); + const t_mdatoms& md); /*! \brief * Update PBC data. diff --git a/src/gromacs/mdrun/md.cpp b/src/gromacs/mdrun/md.cpp index cf64bf19c6..f161798d84 100644 --- a/src/gromacs/mdrun/md.cpp +++ b/src/gromacs/mdrun/md.cpp @@ -473,6 +473,7 @@ void gmx::LegacySimulator::do_md() integrator = std::make_unique( *ir, *top_global, + ekind->ngtc, fr->deviceStreamManager->context(), fr->deviceStreamManager->stream(gmx::DeviceStreamType::UpdateAndConstraints), stateGpu->xUpdatedOnDevice(), @@ -1481,8 +1482,7 @@ void gmx::LegacySimulator::do_md() stateGpu->getVelocities(), stateGpu->getForces(), top.idef, - *mdatoms, - ekind->ngtc); + *mdatoms); // Copy data to the GPU after buffers might have being reinitialized stateGpu->copyVelocitiesToGpu(state->v, AtomLocality::Local); -- 2.22.0