From 6385f296e79e04ff3737d637ae4558296da519cc Mon Sep 17 00:00:00 2001 From: Artem Zhmurov Date: Fri, 17 May 2019 16:08:24 +0200 Subject: [PATCH] Remove PImpl scaffolding from CUDA version of LINCS The CUDA implementation of LINCS was initially introduced as a stand-alone feature. This required hiding CUDA-specific variables and subroutines into the private implementation subclass. Since the LINCS is not a part of Update and Constraints module, this is no longer required and can be removed. Refs #2816, #2888 Change-Id: I9698224d4702dfb8d99106999335c62e83a511df --- src/gromacs/mdlib/CMakeLists.txt | 2 +- .../mdlib/{lincs_cuda_impl.cu => lincs_cuda.cu} | 99 +------ .../mdlib/{lincs_cuda_impl.h => lincs_cuda.cuh} | 61 +--- src/gromacs/mdlib/lincs_cuda.h | 138 --------- src/gromacs/mdlib/lincs_cuda_impl.cpp | 96 ------ src/gromacs/mdlib/tests/CMakeLists.txt | 16 +- src/gromacs/mdlib/tests/constr.cpp | 326 ++++----------------- src/gromacs/mdlib/tests/constr_impl.cpp | 186 ++++++++++++ src/gromacs/mdlib/tests/constr_impl.cu | 113 +++++++ src/gromacs/mdlib/tests/constr_impl.h | 237 +++++++++++++++ src/gromacs/mdlib/update_constrain_cuda_impl.cu | 4 +- src/gromacs/mdlib/update_constrain_cuda_impl.h | 4 +- src/gromacs/pbcutil/pbc_aiuc.h | 7 +- 13 files changed, 639 insertions(+), 650 deletions(-) rename src/gromacs/mdlib/{lincs_cuda_impl.cu => lincs_cuda.cu} (92%) rename src/gromacs/mdlib/{lincs_cuda_impl.h => lincs_cuda.cuh} (76%) delete mode 100644 src/gromacs/mdlib/lincs_cuda.h delete mode 100644 src/gromacs/mdlib/lincs_cuda_impl.cpp create mode 100644 src/gromacs/mdlib/tests/constr_impl.cpp create mode 100644 src/gromacs/mdlib/tests/constr_impl.cu create mode 100644 src/gromacs/mdlib/tests/constr_impl.h diff --git a/src/gromacs/mdlib/CMakeLists.txt b/src/gromacs/mdlib/CMakeLists.txt index d543f04651..9b37d64041 100644 --- a/src/gromacs/mdlib/CMakeLists.txt +++ b/src/gromacs/mdlib/CMakeLists.txt @@ -41,7 +41,7 @@ endif() if(GMX_USE_CUDA) gmx_add_libgromacs_sources( leapfrog_cuda_impl.cu - lincs_cuda_impl.cu + lincs_cuda.cu settle_cuda_impl.cu update_constrain_cuda_impl.cu ) diff --git a/src/gromacs/mdlib/lincs_cuda_impl.cu b/src/gromacs/mdlib/lincs_cuda.cu similarity index 92% rename from src/gromacs/mdlib/lincs_cuda_impl.cu rename to src/gromacs/mdlib/lincs_cuda.cu index eb03a48214..df018e9e91 100644 --- a/src/gromacs/mdlib/lincs_cuda_impl.cu +++ b/src/gromacs/mdlib/lincs_cuda.cu @@ -52,7 +52,7 @@ */ #include "gmxpre.h" -#include "lincs_cuda_impl.h" +#include "lincs_cuda.cuh" #include #include @@ -68,7 +68,6 @@ #include "gromacs/gpu_utils/vectype_ops.cuh" #include "gromacs/math/vec.h" #include "gromacs/mdlib/constr.h" -#include "gromacs/mdlib/lincs_cuda.h" #include "gromacs/pbcutil/pbc.h" #include "gromacs/pbcutil/pbc_aiuc_cuda.cuh" #include "gromacs/topology/ifunc.h" @@ -431,13 +430,13 @@ inline auto getLincsKernelPtr(const bool updateVelocities, return kernelPtr; } -void LincsCuda::Impl::apply(const float3 *d_x, - float3 *d_xp, - const bool updateVelocities, - float3 *d_v, - const real invdt, - const bool computeVirial, - tensor virialScaled) +void LincsCuda::apply(const float3 *d_x, + float3 *d_xp, + const bool updateVelocities, + float3 *d_v, + const real invdt, + const bool computeVirial, + tensor virialScaled) { ensureNoPendingCudaError("In CUDA version of LINCS"); @@ -519,45 +518,8 @@ void LincsCuda::Impl::apply(const float3 *d_x, return; } -void LincsCuda::Impl::copyApplyCopy(const int numAtoms, - const rvec *h_x, - rvec *h_xp, - bool updateVelocities, - rvec *h_v, - real invdt, - bool computeVirial, - tensor virialScaled) -{ - - float3 *d_x, *d_xp, *d_v; - - allocateDeviceBuffer(&d_x, numAtoms, nullptr); - allocateDeviceBuffer(&d_xp, numAtoms, nullptr); - allocateDeviceBuffer(&d_v, numAtoms, nullptr); - - copyToDeviceBuffer(&d_x, (float3*)h_x, 0, numAtoms, stream_, GpuApiCallBehavior::Sync, nullptr); - copyToDeviceBuffer(&d_xp, (float3*)h_xp, 0, numAtoms, stream_, GpuApiCallBehavior::Sync, nullptr); - if (updateVelocities) - { - copyToDeviceBuffer(&d_v, (float3*)h_v, 0, numAtoms, stream_, GpuApiCallBehavior::Sync, nullptr); - } - apply(d_x, d_xp, - updateVelocities, d_v, invdt, - computeVirial, virialScaled); - - copyFromDeviceBuffer((float3*)h_xp, &d_xp, 0, numAtoms, stream_, GpuApiCallBehavior::Sync, nullptr); - if (updateVelocities) - { - copyFromDeviceBuffer((float3*)h_v, &d_v, 0, numAtoms, stream_, GpuApiCallBehavior::Sync, nullptr); - } - - freeDeviceBuffer(&d_x); - freeDeviceBuffer(&d_xp); - freeDeviceBuffer(&d_v); -} - -LincsCuda::Impl::Impl(int numIterations, - int expansionOrder) +LincsCuda::LincsCuda(int numIterations, + int expansionOrder) { kernelParams_.numIterations = numIterations; kernelParams_.expansionOrder = expansionOrder; @@ -579,7 +541,7 @@ LincsCuda::Impl::Impl(int numIterations, } -LincsCuda::Impl::~Impl() +LincsCuda::~LincsCuda() { freeDeviceBuffer(&kernelParams_.d_virialScaled); @@ -635,8 +597,8 @@ inline int countCoupled(int a, std::vector *spaceNeeded, return counted; } -void LincsCuda::Impl::set(const t_idef &idef, - const t_mdatoms &md) +void LincsCuda::set(const t_idef &idef, + const t_mdatoms &md) { int numAtoms = md.nr; // List of constrained atoms (CPU memory) @@ -891,42 +853,9 @@ void LincsCuda::Impl::set(const t_idef &idef, } -void LincsCuda::Impl::setPbc(const t_pbc *pbc) -{ - setPbcAiuc(pbc->ndim_ePBC, pbc->box, &kernelParams_.pbcAiuc); -} - -LincsCuda::LincsCuda(const int numIterations, - const int expansionOrder) - : impl_(new Impl(numIterations, expansionOrder)) -{ -} - -LincsCuda::~LincsCuda() = default; - -void LincsCuda::copyApplyCopy(const int numAtoms, - const rvec *h_x, - rvec *h_xp, - const bool updateVelocities, - rvec *h_v, - const real invdt, - const bool computeVirial, - tensor virialScaled) -{ - impl_->copyApplyCopy(numAtoms, h_x, h_xp, - updateVelocities, h_v, invdt, - computeVirial, virialScaled); -} - void LincsCuda::setPbc(const t_pbc *pbc) { - impl_->setPbc(pbc); -} - -void LincsCuda::set(const t_idef &idef, - const t_mdatoms &md) -{ - impl_->set(idef, md); + setPbcAiuc(pbc->ndim_ePBC, pbc->box, &kernelParams_.pbcAiuc); } } // namespace gmx diff --git a/src/gromacs/mdlib/lincs_cuda_impl.h b/src/gromacs/mdlib/lincs_cuda.cuh similarity index 76% rename from src/gromacs/mdlib/lincs_cuda_impl.h rename to src/gromacs/mdlib/lincs_cuda.cuh index ff21ae88ae..1fb9eada23 100644 --- a/src/gromacs/mdlib/lincs_cuda_impl.h +++ b/src/gromacs/mdlib/lincs_cuda.cuh @@ -32,26 +32,23 @@ * To help us fund GROMACS development, we humbly ask that you cite * the research papers on the package. Check out http://www.gromacs.org. */ -/*! \internal \file +/*! \libinternal \file * - * \brief Declares CUDA implementation class for LINCS - * - * This header file is needed to include from both the device-side - * kernels file, and the host-side management code. + * \brief Declaration of high-level functions of CUDA implementation of LINCS. * * \author Artem Zhmurov * * \ingroup module_mdlib + * \inlibraryapi */ -#ifndef GMX_MDLIB_LINCS_CUDA_IMPL_H -#define GMX_MDLIB_LINCS_CUDA_IMPL_H - +#ifndef GMX_MDLIB_LINCS_CUDA_CUH +#define GMX_MDLIB_LINCS_CUDA_CUH #include "gromacs/mdlib/constr.h" -#include "gromacs/mdlib/lincs_cuda.h" #include "gromacs/mdtypes/mdatom.h" -#include "gromacs/pbcutil/pbc_aiuc_cuda.cuh" +#include "gromacs/pbcutil/pbc_aiuc.h" #include "gromacs/topology/idef.h" +#include "gromacs/utility/classhelpers.h" namespace gmx { @@ -96,7 +93,7 @@ struct LincsCudaKernelParameters }; /*! \internal \brief Class with interfaces and data for CUDA version of LINCS. */ -class LincsCuda::Impl +class LincsCuda { public: @@ -106,10 +103,10 @@ class LincsCuda::Impl * \param[in] numIterations Number of iteration for the correction of the projection. * \param[in] expansionOrder Order of the matrix inversion algorithm. */ - Impl(int numIterations, - int expansionOrder); + LincsCuda(int numIterations, + int expansionOrder); /*! \brief Destructor.*/ - ~Impl(); + ~LincsCuda(); /*! \brief Apply LINCS. * @@ -137,38 +134,6 @@ class LincsCuda::Impl const bool computeVirial, tensor virialScaled); - /*! \brief Apply LINCS to the coordinates/velocities stored in CPU memory. - * - * This method should not be used in any code-path, where performance is of any value. - * Only suitable for test and will be removed in future patch sets. - * Allocates GPU memory, copies data from CPU, applies LINCS to coordinates and, - * if requested, to velocities, copies the results back, frees GPU memory. - * Method uses this class data structures which should be filled with set() and setPbc() - * methods. - * - * \todo Remove this method - * - * \param[in] numAtoms Number of atoms - * \param[in] h_x Coordinates before timestep (in CPU memory) - * \param[in,out] h_xp Coordinates after timestep (in CPU memory). The - * resulting constrained coordinates will be saved here. - * \param[in] updateVelocities If the velocities should be updated. - * \param[in,out] h_v Velocities to update (in CPU memory, can be nullptr - * if not updated) - * \param[in] invdt Reciprocal timestep (to scale Lagrange - * multipliers when velocities are updated) - * \param[in] computeVirial If virial should be updated. - * \param[in,out] virialScaled Scaled virial tensor to be updated. - */ - void copyApplyCopy(const int numAtoms, - const rvec *h_x, - rvec *h_xp, - const bool updateVelocities, - rvec *h_v, - const real invdt, - const bool computeVirial, - tensor virialScaled); - /*! \brief * Update data-structures (e.g. after NB search step). * @@ -234,6 +199,6 @@ class LincsCuda::Impl }; -} // namespace gmx +} // namespace gmx -#endif +#endif // GMX_MDLIB_LINCS_CUDA_CUH diff --git a/src/gromacs/mdlib/lincs_cuda.h b/src/gromacs/mdlib/lincs_cuda.h deleted file mode 100644 index 76014939bb..0000000000 --- a/src/gromacs/mdlib/lincs_cuda.h +++ /dev/null @@ -1,138 +0,0 @@ -/* - * This file is part of the GROMACS molecular simulation package. - * - * Copyright (c) 2019, 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. - * - * GROMACS is free software; you can redistribute it and/or - * modify it under the terms of the GNU Lesser General Public License - * as published by the Free Software Foundation; either version 2.1 - * of the License, or (at your option) any later version. - * - * GROMACS is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU - * Lesser General Public License for more details. - * - * You should have received a copy of the GNU Lesser General Public - * License along with GROMACS; if not, see - * http://www.gnu.org/licenses, or write to the Free Software Foundation, - * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. - * - * If you want to redistribute modifications to GROMACS, please - * consider that scientific software is very special. Version - * control is crucial - bugs must be traceable. We will be happy to - * consider code for inclusion in the official distribution, but - * derived work must not be called official GROMACS. Details are found - * in the README & COPYING files - if they are missing, get the - * official version at http://www.gromacs.org. - * - * To help us fund GROMACS development, we humbly ask that you cite - * the research papers on the package. Check out http://www.gromacs.org. - */ -/*! \libinternal \file - * - * \brief Declaration of high-level functions of CUDA implementation of LINCS. - * - * \todo This should only list interfaces needed for libgromacs clients. - * - * \author Artem Zhmurov - * - * \ingroup module_mdlib - * \inlibraryapi - */ -#ifndef GMX_MDLIB_LINCS_CUDA_H -#define GMX_MDLIB_LINCS_CUDA_H - - -#include "gromacs/mdlib/constr.h" -#include "gromacs/mdtypes/mdatom.h" -#include "gromacs/topology/idef.h" -#include "gromacs/utility/classhelpers.h" - -namespace gmx -{ - -class LincsCuda -{ - - public: - /*! \brief Constructor. - * - * \param[in] numIterations Number of iteration for the correction of the projection. - * \param[in] expansionOrder Order of the matrix inversion algorithm. - */ - LincsCuda(int numIterations, - int expansionOrder); - ~LincsCuda(); - - /*! \brief Apply LINCS to the coordinates/velocities stored in CPU memory. - * - * This method should not be used in any code-path, where performance is of any value. - * Only suitable for test and will be removed in future patch sets. - * Allocates GPU memory, copies data from CPU, applies LINCS to coordinates and, - * if requested, to velocities, copies the results back, frees GPU memory. - * Method uses this class data structures which should be filled with set() and setPbc() - * methods. - * - * \todo Remove this method - * - * \param[in] numAtoms Number of atoms - * \param[in] h_x Coordinates before timestep (in CPU memory) - * \param[in,out] h_xp Coordinates after timestep (in CPU memory). The - * resulting constrained coordinates will be saved here. - * \param[in] updateVelocities If the velocities should be updated. - * \param[in,out] h_v Velocities to update (in CPU memory, can be nullptr - * if not updated) - * \param[in] invdt Reciprocal timestep (to scale Lagrange - * multipliers when velocities are updated) - * \param[in] computeVirial If virial should be updated. - * \param[in,out] virialScaled Scaled virial tensor to be updated. - */ - void copyApplyCopy(int numAtoms, - const rvec *h_x, - rvec *h_xp, - bool updateVelocities, - rvec *h_v, - real invdt, - bool computeVirial, - tensor virialScaled); - - /*! \brief - * Update data-structures (e.g. after NB search step). - * - * Updates the constraints data. Should be called if the particles were sorted, - * redistributed between domains, etc. - * - * Information about constraints should be taken from: - * idef.il[F_CONSTR].iatoms --- type (T) of constraint and two atom indexes (i1, i2) - * idef.iparams[T].constr.dA --- target length for constraint of type T - * From t_mdatom, the code should take: - * md.invmass --- array of inverse square root of masses for each atom in the system. - * - * \param[in] idef Local topology data to get information on constraints from. - * \param[in] md Atoms data to get atom masses from. - */ - void set(const t_idef &idef, - const t_mdatoms &md); - - /*! \brief - * Update PBC data. - * - * \param[in] pbc The PBC data in t_pbc format. - */ - void setPbc(const t_pbc *pbc); - - /*! \brief Class with hardware-specific interfaces and implementations.*/ - class Impl; - - private: - gmx::PrivateImplPointer impl_; - -}; - -} // namespace gmx - -#endif diff --git a/src/gromacs/mdlib/lincs_cuda_impl.cpp b/src/gromacs/mdlib/lincs_cuda_impl.cpp deleted file mode 100644 index 39ab17d48a..0000000000 --- a/src/gromacs/mdlib/lincs_cuda_impl.cpp +++ /dev/null @@ -1,96 +0,0 @@ -/* - * This file is part of the GROMACS molecular simulation package. - * - * Copyright (c) 2019, 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. - * - * GROMACS is free software; you can redistribute it and/or - * modify it under the terms of the GNU Lesser General Public License - * as published by the Free Software Foundation; either version 2.1 - * of the License, or (at your option) any later version. - * - * GROMACS is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU - * Lesser General Public License for more details. - * - * You should have received a copy of the GNU Lesser General Public - * License along with GROMACS; if not, see - * http://www.gnu.org/licenses, or write to the Free Software Foundation, - * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. - * - * If you want to redistribute modifications to GROMACS, please - * consider that scientific software is very special. Version - * control is crucial - bugs must be traceable. We will be happy to - * consider code for inclusion in the official distribution, but - * derived work must not be called official GROMACS. Details are found - * in the README & COPYING files - if they are missing, get the - * official version at http://www.gromacs.org. - * - * To help us fund GROMACS development, we humbly ask that you cite - * the research papers on the package. Check out http://www.gromacs.org. - */ -/*! \internal \file - * - * \brief May be used to implement LINCS CUDA interfaces for non-GPU builds. - * - * Currently, reports and exits if any of the interfaces are called. - * May replace CPU LINCS module in future. Currently used to satisfy - * compiler on systems, where CUDA is not available. - * - * \author Artem Zhmurov - * - * \ingroup module_mdlib - */ -#include "gmxpre.h" - -#include "config.h" - -#include "gromacs/mdlib/lincs_cuda.h" - -#if GMX_GPU != GMX_GPU_CUDA - -namespace gmx -{ - -class LincsCuda::Impl -{ -}; - -LincsCuda::LincsCuda(gmx_unused int numIterations, - gmx_unused int expansionOrder) - : impl_(nullptr) -{ - GMX_ASSERT(false, "A CPU stub for LINCS was called insted of the correct implementation."); -} - -LincsCuda::~LincsCuda() = default; - -void LincsCuda::copyApplyCopy(gmx_unused int numAtoms, - gmx_unused const rvec *h_x, - gmx_unused rvec *h_xp, - gmx_unused bool updateVelocities, - gmx_unused rvec *h_v, - gmx_unused real invdt, - gmx_unused bool computeVirial, - gmx_unused tensor virialScaled) -{ - GMX_ASSERT(false, "A CPU stub for LINCS was called insted of the correct implementation."); -} - -void LincsCuda::set(gmx_unused const t_idef &idef, - gmx_unused const t_mdatoms &md) -{ - GMX_ASSERT(false, "A CPU stub for LINCS was called insted of the correct implementation."); -} - -void LincsCuda::setPbc(gmx_unused const t_pbc *pbc) -{ - GMX_ASSERT(false, "A CPU stub for LINCS was called insted of the correct implementation."); -} - -} // namespace gmx - -#endif /* GMX_GPU != GMX_GPU_CUDA */ diff --git a/src/gromacs/mdlib/tests/CMakeLists.txt b/src/gromacs/mdlib/tests/CMakeLists.txt index fc0c714c1a..e3ac4b1cfe 100644 --- a/src/gromacs/mdlib/tests/CMakeLists.txt +++ b/src/gromacs/mdlib/tests/CMakeLists.txt @@ -32,6 +32,16 @@ # To help us fund GROMACS development, we humbly ask that you cite # the research papers on the package. Check out http://www.gromacs.org. +file(GLOB MDLIB_TEST_SOURCES + constr_impl.cpp + ) + +if (GMX_USE_CUDA) + file(GLOB MDLIB_TEST_CUDA_SOURCES + constr_impl.cu + ) +endif() + gmx_add_unit_test(MdlibUnitTest mdlib-test calc_verletbuf.cpp constr.cpp @@ -42,4 +52,8 @@ gmx_add_unit_test(MdlibUnitTest mdlib-test shake.cpp simulationsignal.cpp updategroups.cpp - updategroupscog.cpp) + updategroupscog.cpp + ${MDLIB_TEST_SOURCES}) + +# TODO: Make CUDA source to compile inside the testing framework +gmx_add_libgromacs_sources(${MDLIB_TEST_CUDA_SOURCES}) \ No newline at end of file diff --git a/src/gromacs/mdlib/tests/constr.cpp b/src/gromacs/mdlib/tests/constr.cpp index 69a0260ed8..2f02a2ca4c 100644 --- a/src/gromacs/mdlib/tests/constr.cpp +++ b/src/gromacs/mdlib/tests/constr.cpp @@ -63,14 +63,10 @@ #include "gromacs/fileio/gmxfio.h" #include "gromacs/gmxlib/nrnb.h" -#include "gromacs/gmxlib/nonbonded/nonbonded.h" -#include "gromacs/gpu_utils/gpu_utils.h" #include "gromacs/math/paddedvector.h" #include "gromacs/math/vec.h" #include "gromacs/math/vectypes.h" -#include "gromacs/mdlib/gmx_omp_nthreads.h" #include "gromacs/mdlib/lincs.h" -#include "gromacs/mdlib/lincs_cuda.h" #include "gromacs/mdlib/shake.h" #include "gromacs/mdrunutility/multisim.h" #include "gromacs/mdtypes/commrec.h" @@ -88,128 +84,25 @@ #include "testutils/refdata.h" #include "testutils/testasserts.h" +#include "constr_impl.h" + namespace gmx { namespace test { -/*! \brief - * Constraints test data structure. - * - * Structure to collect all the necessary data, including system coordinates and topology, - * constraints information, etc. The structure can be reset and reused. - */ -struct ConstraintsTestData +ConstraintsTestData::ConstraintsTestData(const std::string &title, + int numAtoms, std::vector masses, + std::vector constraints, std::vector constraintsR0, + bool computeVirial, tensor virialScaledRef, + bool compute_dHdLambda, float dHdLambdaRef, + real initialTime, real timestep, + const std::vector &x, const std::vector &xPrime, const std::vector &v, + real shakeTolerance, gmx_bool shakeUseSOR, + int lincsNumIterations, int lincsExpansionOrder, real lincsWarnAngle) { - public: - //! Human-friendly name for a system - std::string title_; - //! Number of atoms - int numAtoms_; - //! Topology - gmx_mtop_t mtop_; - //! Masses - std::vector masses_; - //! Inverse masses - std::vector invmass_; - //! Communication record - t_commrec cr_; - //! Input record (info that usually in .mdp file) - t_inputrec ir_; - //! Local topology - t_idef idef_; - //! MD atoms - t_mdatoms md_; - //! Multisim data - gmx_multisim_t ms_; - //! Computational time array (normally used to benchmark performance) - t_nrnb nrnb_; - - //! Inverse timestep - real invdt_; - //! Number of flexible constraints - int nflexcon_ = 0; - //! Whether the virial should be computed - bool computeVirial_; - //! Scaled virial - tensor virialScaled_; - //! Scaled virial (reference values) - tensor virialScaledRef_; - //! If the free energy is computed - bool compute_dHdLambda_; - //! For free energy computation - real dHdLambda_; - //! For free energy computation (reference value) - real dHdLambdaRef_; - - //! Coordinates before the timestep - PaddedVector x_; - //! Coordinates after timestep, output for the constraints - PaddedVector xPrime_; - //! Backup for coordinates (for reset) - PaddedVector xPrime0_; - //! Intermediate set of coordinates (normally used for projection correction) - PaddedVector xPrime2_; - //! Velocities - PaddedVector v_; - //! Backup for velocities (for reset) - PaddedVector v0_; - - //! Constraints data (type1-i1-j1-type2-i2-j2-...) - std::vector constraints_; - //! Target lengths for all constraint types - std::vector constraintsR0_; - - /*! \brief - * Constructor for the object with all parameters and variables needed by constraints algorithms. - * - * This constructor assembles stubs for all the data structures, required to initialize - * and apply LINCS and SHAKE constraints. The coordinates and velocities before constraining - * are saved to allow for reset. The constraints data are stored for testing after constraints - * were applied. - * - * \param[in] title Human-friendly name of the system. - * \param[in] numAtoms Number of atoms in the system. - * \param[in] masses Atom masses. Size of this vector should be equal to numAtoms. - * \param[in] constraints List of constraints, organized in triples of integers. - * First integer is the index of type for a constraint, second - * and third are the indices of constrained atoms. The types - * of constraints should be sequential but not necessarily - * start from zero (which is the way they normally are in - * GROMACS). - * \param[in] constraintsR0 Target values for bond lengths for bonds of each type. The - * size of this vector should be equal to the total number of - * unique types in constraints vector. - * \param[in] computeVirial Whether the virial should be computed. - * \param[in] virialScaledRef Reference values for scaled virial tensor. - * \param[in] compute_dHdLambda Whether free energy should be computed. - * \param[in] dHdLambdaRef Reference value for dHdLambda. - * \param[in] initialTime Initial time. - * \param[in] timestep Timestep. - * \param[in] x Coordinates before integration step. - * \param[in] xPrime Coordinates after integration step, but before constraining. - * \param[in] v Velocities before constraining. - * \param[in] shakeTolerance Target tolerance for SHAKE. - * \param[in] shakeUseSOR Use successive over-relaxation method for SHAKE iterations. - * The general formula is: - * x_n+1 = (1-omega)*x_n + omega*f(x_n), - * where omega = 1 if SOR is off and may be < 1 if SOR is on. - * \param[in] lincsNumIterations Number of iterations used to compute the inverse matrix. - * \param[in] lincsExpansionOrder The order for algorithm that adjusts the direction of the - * bond after constraints are applied. - * \param[in] lincsWarnAngle The threshold value for the change in bond angle. When - * exceeded the program will issue a warning. - * - */ - ConstraintsTestData(const std::string &title, - int numAtoms, std::vector masses, - std::vector constraints, std::vector constraintsR0, - bool computeVirial, tensor virialScaledRef, - bool compute_dHdLambda, float dHdLambdaRef, - real initialTime, real timestep, - const std::vector &x, const std::vector &xPrime, const std::vector &v, - real shakeTolerance, gmx_bool shakeUseSOR, - int lincsNumIterations, int lincsExpansionOrder, real lincsWarnAngle) + // This is to trick Gerrit + { { title_ = title; // Human-friendly name of the system numAtoms_ = numAtoms; // Number of atoms @@ -361,43 +254,46 @@ struct ConstraintsTestData ir_.nProjOrder = lincsExpansionOrder; ir_.LincsWarnAngle = lincsWarnAngle; } + } +} - /*! \brief - * Reset the data structure so it can be reused. - * - * Set the coordinates and velocities back to their values before - * constraining. The scaled virial tensor and dHdLambda are zeroed. - * - */ - void reset() - { - xPrime_ = xPrime0_; - xPrime2_ = xPrime0_; - v_ = v0_; +/*! \brief + * Reset the data structure so it can be reused. + * + * Set the coordinates and velocities back to their values before + * constraining. The scaled virial tensor and dHdLambda are zeroed. + * + */ +void ConstraintsTestData::reset() +{ + xPrime_ = xPrime0_; + xPrime2_ = xPrime0_; + v_ = v0_; - if (computeVirial_) + if (computeVirial_) + { + for (int i = 0; i < DIM; i++) + { + for (int j = 0; j < DIM; j++) { - for (int i = 0; i < DIM; i++) - { - for (int j = 0; j < DIM; j++) - { - virialScaled_[i][j] = 0; - } - } + virialScaled_[i][j] = 0; } - dHdLambda_ = 0; } + } + dHdLambda_ = 0; +} - /*! \brief - * Cleaning up the memory. - */ - ~ConstraintsTestData() - { - sfree(idef_.il[F_CONSTR].iatoms); - sfree(idef_.iparams); - } +/*! \brief + * Cleaning up the memory. + */ +ConstraintsTestData::~ConstraintsTestData() +{ + sfree(idef_.il[F_CONSTR].iatoms); + sfree(idef_.iparams); +} -}; +namespace +{ /*! \brief The two-dimensional parameter space for test. * @@ -408,14 +304,11 @@ struct ConstraintsTestData */ typedef std::tuple ConstraintsTestParameters; -/*! \brief Names of all availible algorithms - */ -static std::vector algorithmsNames; +//! Names of all availible algorithms +std::vector algorithmsNames; -/*! \brief Method that fills and returns algorithmNames to the test macros. - * - */ -static std::vector getAlgorithmsNames() +//! Method that fills and returns algorithmNames to the test macros. +std::vector getAlgorithmsNames() { algorithmsNames.emplace_back("SHAKE"); algorithmsNames.emplace_back("LINCS"); @@ -484,127 +377,11 @@ class ConstraintsTest : public ::testing::TestWithParamidef_, testData->md_); - bool success = constrain_shake( - nullptr, - shaked, - testData->invmass_.data(), - testData->idef_, - testData->ir_, - as_rvec_array(testData->x_.data()), - as_rvec_array(testData->xPrime_.data()), - as_rvec_array(testData->xPrime2_.data()), - &testData->nrnb_, - testData->md_.lambda, - &testData->dHdLambda_, - testData->invdt_, - as_rvec_array(testData->v_.data()), - testData->computeVirial_, - testData->virialScaled_, - false, - gmx::ConstraintVariable::Positions); - EXPECT_TRUE(success) << "Test failed with a false return value in SHAKE."; - done_shake(shaked); - } - - /*! \brief - * Initialize and apply LINCS constraints. - * - * \param[in] testData Test data structure. - * \param[in] pbc Periodic boundary data. - */ - static void applyLincs(ConstraintsTestData *testData, t_pbc pbc) - { - - Lincs *lincsd; - int maxwarn = 100; - int warncount_lincs = 0; - gmx_omp_nthreads_set(emntLINCS, 1); - - // Make blocka structure for faster LINCS setup - std::vector at2con_mt; - at2con_mt.reserve(testData->mtop_.moltype.size()); - for (const gmx_moltype_t &moltype : testData->mtop_.moltype) - { - // This function is in constr.cpp - at2con_mt.push_back(make_at2con(moltype, - testData->mtop_.ffparams.iparams, - flexibleConstraintTreatment(EI_DYNAMICS(testData->ir_.eI)))); - } - // Initialize LINCS - lincsd = init_lincs(nullptr, testData->mtop_, - testData->nflexcon_, at2con_mt, - false, - testData->ir_.nLincsIter, testData->ir_.nProjOrder); - set_lincs(testData->idef_, testData->md_, EI_DYNAMICS(testData->ir_.eI), &testData->cr_, lincsd); - - // Evaluate constraints - bool sucess = constrain_lincs(false, testData->ir_, 0, lincsd, testData->md_, - &testData->cr_, - &testData->ms_, - as_rvec_array(testData->x_.data()), - as_rvec_array(testData->xPrime_.data()), - as_rvec_array(testData->xPrime2_.data()), - pbc.box, &pbc, - testData->md_.lambda, &testData->dHdLambda_, - testData->invdt_, - as_rvec_array(testData->v_.data()), - testData->computeVirial_, testData->virialScaled_, - gmx::ConstraintVariable::Positions, - &testData->nrnb_, - maxwarn, &warncount_lincs); - EXPECT_TRUE(sucess) << "Test failed with a false return value in LINCS."; - EXPECT_EQ(warncount_lincs, 0) << "There were warnings in LINCS."; - for (unsigned int i = 0; i < testData->mtop_.moltype.size(); i++) - { - sfree(at2con_mt.at(i).index); - sfree(at2con_mt.at(i).a); - } - done_lincs(lincsd); - } - - /*! \brief - * Initialize and apply LINCS constraints on CUDA-enabled GPU. - * - * \param[in] testData Test data structure. - * \param[in] pbc Periodic boundary data. - */ - static void applyLincsCuda(ConstraintsTestData *testData, t_pbc pbc) - { - std::string errorMessage; - GMX_RELEASE_ASSERT(GMX_GPU == GMX_GPU_CUDA, - "The test for CUDA version of LINCS was called non-CUDA build."); - GMX_RELEASE_ASSERT(canDetectGpus(&errorMessage), - "The test for CUDA version of LINCS was called on host that is not CUDA capable."); - - auto lincsCuda = std::make_unique(testData->ir_.nLincsIter, - testData->ir_.nProjOrder); - lincsCuda->set(testData->idef_, testData->md_); - lincsCuda->setPbc(&pbc); - lincsCuda->copyApplyCopy(testData->numAtoms_, - as_rvec_array(testData->x_.data()), - as_rvec_array(testData->xPrime_.data()), - true, - as_rvec_array(testData->v_.data()), - testData->invdt_, - testData->computeVirial_, - testData->virialScaled_); - } - - /*! \brief * The test on the final length of constrained bonds. * * Goes through all the constraints and checks if the final length of all the constraints is equal @@ -1172,5 +949,6 @@ INSTANTIATE_TEST_CASE_P(WithParameters, ConstraintsTest, ::testing::Combine(::testing::Values("PBCNone", "PBCXYZ"), ::testing::ValuesIn(getAlgorithmsNames()))); +} // namespace } // namespace test } // namespace gmx diff --git a/src/gromacs/mdlib/tests/constr_impl.cpp b/src/gromacs/mdlib/tests/constr_impl.cpp new file mode 100644 index 0000000000..452fbd71e1 --- /dev/null +++ b/src/gromacs/mdlib/tests/constr_impl.cpp @@ -0,0 +1,186 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2018,2019, 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. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief Function to run SHAKE and LINCS on CPU. + * + * Functions used in the test to apply constraints on the test data: + * CPU-based implementation and a stub for GPU-based implementation. + * + * \author Artem Zhmurov + * \ingroup module_mdlib + */ + +#include "gmxpre.h" + +#include "constr_impl.h" + +#include "config.h" + +#include + +#include + +#include +#include +#include + +#include + +#include "gromacs/math/paddedvector.h" +#include "gromacs/math/vec.h" +#include "gromacs/math/vectypes.h" +#include "gromacs/mdlib/constr.h" +#include "gromacs/mdlib/lincs.h" +#include "gromacs/mdlib/shake.h" +#include "gromacs/mdtypes/commrec.h" +#include "gromacs/mdtypes/inputrec.h" +#include "gromacs/mdtypes/mdatom.h" +#include "gromacs/pbcutil/pbc.h" +#include "gromacs/topology/block.h" +#include "gromacs/topology/idef.h" +#include "gromacs/topology/ifunc.h" +#include "gromacs/topology/topology.h" +#include "gromacs/utility/unique_cptr.h" + +#include "testutils/testasserts.h" + +namespace gmx +{ +namespace test +{ + +/*! \brief + * Initialize and apply SHAKE constraints. + * + * \param[in] testData Test data structure. + * \param[in] pbc Periodic boundary data. + */ +void applyShake(ConstraintsTestData *testData, t_pbc gmx_unused pbc) +{ + shakedata* shaked = shake_init(); + make_shake_sblock_serial(shaked, &testData->idef_, testData->md_); + bool success = constrain_shake( + nullptr, + shaked, + testData->invmass_.data(), + testData->idef_, + testData->ir_, + as_rvec_array(testData->x_.data()), + as_rvec_array(testData->xPrime_.data()), + as_rvec_array(testData->xPrime2_.data()), + &testData->nrnb_, + testData->md_.lambda, + &testData->dHdLambda_, + testData->invdt_, + as_rvec_array(testData->v_.data()), + testData->computeVirial_, + testData->virialScaled_, + false, + gmx::ConstraintVariable::Positions); + EXPECT_TRUE(success) << "Test failed with a false return value in SHAKE."; + done_shake(shaked); +} + +/*! \brief + * Initialize and apply LINCS constraints. + * + * \param[in] testData Test data structure. + * \param[in] pbc Periodic boundary data. + */ +void applyLincs(ConstraintsTestData *testData, t_pbc pbc) +{ + + Lincs *lincsd; + int maxwarn = 100; + int warncount_lincs = 0; + gmx_omp_nthreads_set(emntLINCS, 1); + + // Make blocka structure for faster LINCS setup + std::vector at2con_mt; + at2con_mt.reserve(testData->mtop_.moltype.size()); + for (const gmx_moltype_t &moltype : testData->mtop_.moltype) + { + // This function is in constr.cpp + at2con_mt.push_back(make_at2con(moltype, + testData->mtop_.ffparams.iparams, + flexibleConstraintTreatment(EI_DYNAMICS(testData->ir_.eI)))); + } + // Initialize LINCS + lincsd = init_lincs(nullptr, testData->mtop_, + testData->nflexcon_, at2con_mt, + false, + testData->ir_.nLincsIter, testData->ir_.nProjOrder); + set_lincs(testData->idef_, testData->md_, EI_DYNAMICS(testData->ir_.eI), &testData->cr_, lincsd); + + // Evaluate constraints + bool success = constrain_lincs(false, testData->ir_, 0, lincsd, testData->md_, + &testData->cr_, + &testData->ms_, + as_rvec_array(testData->x_.data()), + as_rvec_array(testData->xPrime_.data()), + as_rvec_array(testData->xPrime2_.data()), + pbc.box, &pbc, + testData->md_.lambda, &testData->dHdLambda_, + testData->invdt_, + as_rvec_array(testData->v_.data()), + testData->computeVirial_, testData->virialScaled_, + gmx::ConstraintVariable::Positions, + &testData->nrnb_, + maxwarn, &warncount_lincs); + EXPECT_TRUE(success) << "Test failed with a false return value in LINCS."; + EXPECT_EQ(warncount_lincs, 0) << "There were warnings in LINCS."; + for (unsigned int i = 0; i < testData->mtop_.moltype.size(); i++) + { + sfree(at2con_mt.at(i).index); + sfree(at2con_mt.at(i).a); + } + done_lincs(lincsd); +} + +#if GMX_GPU != GMX_GPU_CUDA +/*! \brief + * Stub for LINCS constraints on CUDA-enabled GPU to satisfy compiler. + * + * \param[in] testData Test data structure. + * \param[in] pbc Periodic boundary data. + */ +void applyLincsCuda(ConstraintsTestData gmx_unused *testData, t_pbc gmx_unused pbc) +{ + FAIL() << "Dummy LINCS CUDA function was called instead of the real one."; +} +#endif + +} // namespace test +} // namespace gmx diff --git a/src/gromacs/mdlib/tests/constr_impl.cu b/src/gromacs/mdlib/tests/constr_impl.cu new file mode 100644 index 0000000000..45027bdbb8 --- /dev/null +++ b/src/gromacs/mdlib/tests/constr_impl.cu @@ -0,0 +1,113 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2018,2019, 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. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief Subroutines to run LINCS on GPU + * + * Copies data to GPU, runs LINCS and copies the results back. + * + * \author Artem Zhmurov + * \ingroup module_mdlib + */ +#include "gmxpre.h" + +#include "constr_impl.h" + +#include + +#include + +#include +#include +#include + +#include "gromacs/gpu_utils/devicebuffer.cuh" +#include "gromacs/gpu_utils/gpu_utils.h" +#include "gromacs/math/vec.h" +#include "gromacs/math/vectypes.h" +#include "gromacs/mdlib/constr.h" +#include "gromacs/mdlib/lincs_cuda.cuh" +#include "gromacs/pbcutil/pbc.h" +#include "gromacs/utility/unique_cptr.h" + +namespace gmx +{ +namespace test +{ + +/*! \brief + * Initialize and apply LINCS constraints on CUDA-enabled GPU. + * + * \param[in] testData Test data structure. + * \param[in] pbc Periodic boundary data. + */ +void applyLincsCuda(ConstraintsTestData *testData, t_pbc pbc) +{ + auto lincsCuda = std::make_unique(testData->ir_.nLincsIter, + testData->ir_.nProjOrder); + + bool updateVelocities = true; + int numAtoms = testData->numAtoms_; + float3 *d_x, *d_xp, *d_v; + + lincsCuda->set(testData->idef_, testData->md_); + lincsCuda->setPbc(&pbc); + + allocateDeviceBuffer(&d_x, numAtoms, nullptr); + allocateDeviceBuffer(&d_xp, numAtoms, nullptr); + allocateDeviceBuffer(&d_v, numAtoms, nullptr); + + copyToDeviceBuffer(&d_x, (float3*)(testData->x_.data()), 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr); + copyToDeviceBuffer(&d_xp, (float3*)(testData->xPrime_.data()), 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr); + if (updateVelocities) + { + copyToDeviceBuffer(&d_v, (float3*)(testData->v_.data()), 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr); + } + lincsCuda->apply(d_x, d_xp, + updateVelocities, d_v, testData->invdt_, + testData->computeVirial_, testData->virialScaled_); + + copyFromDeviceBuffer((float3*)(testData->xPrime_.data()), &d_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr); + if (updateVelocities) + { + copyFromDeviceBuffer((float3*)(testData->v_.data()), &d_v, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr); + } + + freeDeviceBuffer(&d_x); + freeDeviceBuffer(&d_xp); + freeDeviceBuffer(&d_v); +} + +} // namespace test +} // namespace gmx diff --git a/src/gromacs/mdlib/tests/constr_impl.h b/src/gromacs/mdlib/tests/constr_impl.h new file mode 100644 index 0000000000..c46d00e2f4 --- /dev/null +++ b/src/gromacs/mdlib/tests/constr_impl.h @@ -0,0 +1,237 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2018,2019, 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. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief SHAKE and LINCS tests header. + * + * Contains description and constructor for the test data accumulating object, + * declares CPU- and GPU-based functions used to apply SHAKE or LINCS on the + * test data. + * + * \author Artem Zhmurov + * \ingroup module_mdlib + */ + +#ifndef GMX_MDLIB_TESTS_CONSTR_IMPL_H +#define GMX_MDLIB_TESTS_CONSTR_IMPL_H + +#include "gmxpre.h" + +#include + +#include + +#include +#include +#include + +#include "gromacs/fileio/gmxfio.h" +#include "gromacs/gmxlib/nrnb.h" +#include "gromacs/gmxlib/nonbonded/nonbonded.h" +#include "gromacs/gpu_utils/gpu_utils.h" +#include "gromacs/math/paddedvector.h" +#include "gromacs/math/vec.h" +#include "gromacs/math/vectypes.h" +#include "gromacs/mdlib/constr.h" +#include "gromacs/mdlib/gmx_omp_nthreads.h" +#include "gromacs/mdlib/lincs.h" +#include "gromacs/mdlib/shake.h" +#include "gromacs/mdrunutility/multisim.h" +#include "gromacs/mdtypes/commrec.h" +#include "gromacs/mdtypes/inputrec.h" +#include "gromacs/mdtypes/mdatom.h" +#include "gromacs/pbcutil/pbc.h" +#include "gromacs/topology/block.h" +#include "gromacs/topology/idef.h" +#include "gromacs/topology/ifunc.h" +#include "gromacs/topology/topology.h" +#include "gromacs/utility/smalloc.h" +#include "gromacs/utility/stringutil.h" +#include "gromacs/utility/unique_cptr.h" + +namespace gmx +{ +namespace test +{ + +/* \brief + * Constraints test data structure. + * + * Structure to collect all the necessary data, including system coordinates and topology, + * constraints information, etc. The structure can be reset and reused. + */ +class ConstraintsTestData +{ + public: + //! Human-friendly name for a system + std::string title_; + //! Number of atoms + int numAtoms_; + //! Topology + gmx_mtop_t mtop_; + //! Masses + std::vector masses_; + //! Inverse masses + std::vector invmass_; + //! Communication record + t_commrec cr_; + //! Input record (info that usually in .mdp file) + t_inputrec ir_; + //! Local topology + t_idef idef_; + //! MD atoms + t_mdatoms md_; + //! Multisim data + gmx_multisim_t ms_; + //! Computational time array (normally used to benchmark performance) + t_nrnb nrnb_; + + //! Inverse timestep + real invdt_; + //! Number of flexible constraints + int nflexcon_ = 0; + //! Whether the virial should be computed + bool computeVirial_; + //! Scaled virial + tensor virialScaled_; + //! Scaled virial (reference values) + tensor virialScaledRef_; + //! If the free energy is computed + bool compute_dHdLambda_; + //! For free energy computation + real dHdLambda_; + //! For free energy computation (reference value) + real dHdLambdaRef_; + + //! Coordinates before the timestep + PaddedVector x_; + //! Coordinates after timestep, output for the constraints + PaddedVector xPrime_; + //! Backup for coordinates (for reset) + PaddedVector xPrime0_; + //! Intermediate set of coordinates (normally used for projection correction) + PaddedVector xPrime2_; + //! Velocities + PaddedVector v_; + //! Backup for velocities (for reset) + PaddedVector v0_; + + //! Constraints data (type1-i1-j1-type2-i2-j2-...) + std::vector constraints_; + //! Target lengths for all constraint types + std::vector constraintsR0_; + + /*! \brief + * Constructor for the object with all parameters and variables needed by constraints algorithms. + * + * This constructor assembles stubs for all the data structures, required to initialize + * and apply LINCS and SHAKE constraints. The coordinates and velocities before constraining + * are saved to allow for reset. The constraints data are stored for testing after constraints + * were applied. + * + * \param[in] title Human-friendly name of the system. + * \param[in] numAtoms Number of atoms in the system. + * \param[in] masses Atom masses. Size of this vector should be equal to numAtoms. + * \param[in] constraints List of constraints, organized in triples of integers. + * First integer is the index of type for a constraint, second + * and third are the indices of constrained atoms. The types + * of constraints should be sequential but not necessarily + * start from zero (which is the way they normally are in + * GROMACS). + * \param[in] constraintsR0 Target values for bond lengths for bonds of each type. The + * size of this vector should be equal to the total number of + * unique types in constraints vector. + * \param[in] computeVirial Whether the virial should be computed. + * \param[in] virialScaledRef Reference values for scaled virial tensor. + * \param[in] compute_dHdLambda Whether free energy should be computed. + * \param[in] dHdLambdaRef Reference value for dHdLambda. + * \param[in] initialTime Initial time. + * \param[in] timestep Timestep. + * \param[in] x Coordinates before integration step. + * \param[in] xPrime Coordinates after integration step, but before constraining. + * \param[in] v Velocities before constraining. + * \param[in] shakeTolerance Target tolerance for SHAKE. + * \param[in] shakeUseSOR Use successive over-relaxation method for SHAKE iterations. + * The general formula is: + * x_n+1 = (1-omega)*x_n + omega*f(x_n), + * where omega = 1 if SOR is off and may be < 1 if SOR is on. + * \param[in] lincsNumIterations Number of iterations used to compute the inverse matrix. + * \param[in] lincsExpansionOrder The order for algorithm that adjusts the direction of the + * bond after constraints are applied. + * \param[in] lincsWarnAngle The threshold value for the change in bond angle. When + * exceeded the program will issue a warning. + * + */ + ConstraintsTestData(const std::string &title, + int numAtoms, std::vector masses, + std::vector constraints, std::vector constraintsR0, + bool computeVirial, tensor virialScaledRef, + bool compute_dHdLambda, float dHdLambdaRef, + real initialTime, real timestep, + const std::vector &x, const std::vector &xPrime, const std::vector &v, + real shakeTolerance, gmx_bool shakeUseSOR, + int lincsNumIterations, int lincsExpansionOrder, real lincsWarnAngle); + + /*! \brief + * Reset the data structure so it can be reused. + * + * Set the coordinates and velocities back to their values before + * constraining. The scaled virial tensor and dHdLambda are zeroed. + * + */ + void reset(); + + /*! \brief + * Cleaning up the memory. + */ + ~ConstraintsTestData(); +}; + +/*! \brief Apply SHAKE constraints to the test data. + */ +void applyShake(ConstraintsTestData *testData, t_pbc pbc); +/*! \brief Apply LINCS constraints to the test data. + */ +void applyLincs(ConstraintsTestData *testData, t_pbc pbc); +/*! \brief Apply CUDA version of LINCS constraints to the test data. + * + * All the data is copied to the GPU device, then LINCS is applied and + * the resulting coordinates are copied back. + */ +void applyLincsCuda(ConstraintsTestData *testData, t_pbc pbc); + +} // namespace test +} // namespace gmx + +#endif // GMX_MDLIB_TESTS_CONSTR_IMPL_H diff --git a/src/gromacs/mdlib/update_constrain_cuda_impl.cu b/src/gromacs/mdlib/update_constrain_cuda_impl.cu index e1ec5e4022..4f66f1260e 100644 --- a/src/gromacs/mdlib/update_constrain_cuda_impl.cu +++ b/src/gromacs/mdlib/update_constrain_cuda_impl.cu @@ -64,12 +64,12 @@ #include "gromacs/gpu_utils/gputraits.cuh" #include "gromacs/gpu_utils/vectype_ops.cuh" #include "gromacs/math/vec.h" +#include "gromacs/mdlib/lincs_cuda.cuh" #include "gromacs/mdlib/update_constrain_cuda.h" #include "gromacs/pbcutil/pbc.h" #include "gromacs/pbcutil/pbc_aiuc_cuda.cuh" #include "leapfrog_cuda_impl.h" -#include "lincs_cuda_impl.h" #include "settle_cuda_impl.h" namespace gmx @@ -121,7 +121,7 @@ UpdateConstrainCuda::Impl::Impl(int numAtoms, GMX_RELEASE_ASSERT(numAtoms == mtop.natoms, "State and topology number of atoms should be the same."); integrator_ = std::make_unique(); - lincsCuda_ = std::make_unique(ir.nLincsIter, ir.nProjOrder); + lincsCuda_ = std::make_unique(ir.nLincsIter, ir.nProjOrder); settleCuda_ = std::make_unique(mtop); } diff --git a/src/gromacs/mdlib/update_constrain_cuda_impl.h b/src/gromacs/mdlib/update_constrain_cuda_impl.h index cd360b901f..1e6869c21b 100644 --- a/src/gromacs/mdlib/update_constrain_cuda_impl.h +++ b/src/gromacs/mdlib/update_constrain_cuda_impl.h @@ -46,6 +46,7 @@ #ifndef GMX_MDLIB_UPDATE_CONSTRAIN_CUDA_IMPL_H #define GMX_MDLIB_UPDATE_CONSTRAIN_CUDA_IMPL_H +#include "gromacs/mdlib/lincs_cuda.cuh" #include "gromacs/mdlib/update_constrain_cuda.h" #include "gromacs/mdtypes/inputrec.h" #include "gromacs/pbcutil/pbc.h" @@ -53,7 +54,6 @@ #include "gromacs/topology/idef.h" #include "leapfrog_cuda_impl.h" -#include "lincs_cuda_impl.h" #include "settle_cuda_impl.h" namespace gmx @@ -205,7 +205,7 @@ class UpdateConstrainCuda::Impl //! Leap-Frog integrator std::unique_ptr integrator_; //! LINCS CUDA object to use for non-water constraints - std::unique_ptr lincsCuda_; + std::unique_ptr lincsCuda_; //! SETTLE CUDA object for water constrains std::unique_ptr settleCuda_; diff --git a/src/gromacs/pbcutil/pbc_aiuc.h b/src/gromacs/pbcutil/pbc_aiuc.h index 32003071f9..72811bcdf9 100644 --- a/src/gromacs/pbcutil/pbc_aiuc.h +++ b/src/gromacs/pbcutil/pbc_aiuc.h @@ -96,9 +96,10 @@ struct PbcAiuc * \param[out] pbcAiuc Pointer to PbcAiuc structure to be initialized. * */ -static void setPbcAiuc(int numPbcDim, - const matrix box, - PbcAiuc *pbcAiuc) +static inline +void setPbcAiuc(int numPbcDim, + const matrix box, + PbcAiuc *pbcAiuc) { pbcAiuc->invBoxDiagZ = 0.0f; -- 2.11.4.GIT