Refactor Leap-Frog tests and connect them to CPU version
authorArtem Zhmurov <zhmurov@gmail.com>
Tue, 30 Jul 2019 14:30:45 +0000 (16:30 +0200)
committerMark Abraham <mark.j.abraham@gmail.com>
Thu, 8 Aug 2019 09:13:34 +0000 (11:13 +0200)
This introduces test data object and runners to the Leap-Frog
tests, which are now connected to the CPU version of Leap-Frog.
This also makes possible to include tests based on the reference
values, which are needed to make sure that the temperature and(or)
pressure control works fine in new implementations.

Refs. #2816, #2888.

Change-Id: Id2d934c43138889ad178a94126cab4da2895bb5a

src/gromacs/mdlib/tests/CMakeLists.txt
src/gromacs/mdlib/tests/leapfrog.cpp
src/gromacs/mdlib/tests/leapfrogtestdata.cpp [new file with mode: 0644]
src/gromacs/mdlib/tests/leapfrogtestdata.h [new file with mode: 0644]
src/gromacs/mdlib/tests/leapfrogtestrunners.cpp [new file with mode: 0644]
src/gromacs/mdlib/tests/leapfrogtestrunners.cu [moved from src/gromacs/mdlib/tests/leapfrog.cu with 63% similarity]
src/gromacs/mdlib/tests/leapfrogtestrunners.h [moved from src/gromacs/mdlib/tests/leapfrogtest.h with 68% similarity]

index 8491d15bb25eb868698482ec87cecedb3acbefe9..c7cf5620a5b2121cb5579f70d9a49bd07701ffdf 100644 (file)
@@ -40,6 +40,8 @@ gmx_add_unit_test(MdlibUnitTest mdlib-test
                   ebin.cpp
                   energyoutput.cpp
                   leapfrog.cpp
+                  leapfrogtestdata.cpp
+                  leapfrogtestrunners.cpp
                   settle.cpp
                   shake.cpp
                   simulationsignal.cpp
@@ -49,6 +51,6 @@ gmx_add_unit_test(MdlibUnitTest mdlib-test
 # TODO: Make CUDA source to compile inside the testing framework
 if(GMX_USE_CUDA)
     gmx_add_libgromacs_sources(constrtestrunners.cu
-                               leapfrog.cu
+                               leapfrogtestrunners.cu
                                settle_runners.cu)
 endif()
index b4492b8ecefe1f55c843c9e6f85b38ddba6dfe0e..8838ac492f4866cffd08f455c8be3524887978b2 100644 (file)
  * the research papers on the package. Check out http://www.gromacs.org.
  */
 /*! \internal \file
- * \brief Integrator tests
+ * \brief Tests for the Leap-Frog integrator
  *
- * Test for CUDA implementation of the Leap-Frog integrator.
- *
- * \todo Connect to the CPU-based version.
- * \todo Prepare for temperature and pressure controlled integrators.
+ * \todo Add tests for temperature and pressure controlled integrators.
  * \todo Add PBC handling test.
- * \todo Take over coordinates/velocities/forces handling - this infrastructure
- *       will be removed from integrator.
+ * \todo Reference values tests.
  *
  * \author Artem Zhmurov <zhmurov@gmail.com>
  * \ingroup module_mdlib
 
 #include "gromacs/gpu_utils/gpu_testutils.h"
 #include "gromacs/math/vec.h"
+#include "gromacs/math/vectypes.h"
 #include "gromacs/mdtypes/mdatom.h"
 #include "gromacs/utility/stringutil.h"
 
 #include "testutils/refdata.h"
 #include "testutils/testasserts.h"
 
-#include "leapfrogtest.h"
+#include "leapfrogtestdata.h"
+#include "leapfrogtestrunners.h"
 
 namespace gmx
 {
 namespace test
 {
 
-#if GMX_GPU == GMX_GPU_CUDA
-
 /*! \brief The parameter space for the test.
  *
  * The test will run for all possible combinations of accessible
@@ -88,84 +84,32 @@ namespace test
  * 6-8. Force components
  * 9. Number of steps
  */
-typedef std::tuple<int, real, real, real, real, real, real, real, int> IntegratorTestParameters;
+typedef std::tuple<int, real, real, real, real, real, real, real, int> LeapFrogTestParameters;
 
-/*! \brief Test fixture for integrator.
+/*! \brief Test fixture for LeapFrog integrator.
  *
  *  Creates a system of independent particles exerting constant external forces,
  *  makes several numerical integration timesteps and compares the result
  *  with analytical solution.
  *
  */
-class IntegratorTest : public ::testing::TestWithParam<IntegratorTestParameters>
+class LeapFrogTest : public ::testing::TestWithParam<LeapFrogTestParameters>
 {
     public:
-        //! Number of atoms in the system
-        int  numAtoms_;
-        //! Integration timestep
-        real timestep_;
-
-        //! Initial coordinates
-        std::vector<RVec> x0_;
-        //! Current coordinates
-        std::vector<RVec> x_;
-        //! Coordinates after integrator update
-        std::vector<RVec> xPrime_;
-        //! Initial velocities
-        std::vector<RVec> v0_;
-        //! Current velocities
-        std::vector<RVec> v_;
-        //! External forces
-        std::vector<RVec> f_;
-        //! Inverse masses of the particles
-        std::vector<real> inverseMasses_;
-        //! MD atoms structure in which inverse masses will be passed to the integrator
-        t_mdatoms         mdAtoms_;
+        //! Availiable runners (CPU and GPU versions of the Leap-Frog)
+        std::unordered_map <std::string, void(*)(LeapFrogTestData *testData,
+                                                 const int         numSteps)> runners_;
 
-        /*! \brief Test initialization function.
-         *
-         * \param[in]  numAtoms  Number of atoms in the system
-         * \param[in]  timestep  Integration timestep
-         * \param[in]  v0        Initial velocity (same for all particles)
-         * \param[in]  f0        External constant force, acting on all particles
-         */
-        void init(int numAtoms, real timestep, rvec v0, rvec f0)
+        LeapFrogTest()
         {
-            numAtoms_    = numAtoms;
-            timestep_    = timestep;
-
-            x0_.resize(numAtoms_);
-            x_.resize(numAtoms_);
-            xPrime_.resize(numAtoms_);
-            v0_.resize(numAtoms_);
-            v_.resize(numAtoms_);
-            f_.resize(numAtoms_);
-            inverseMasses_.resize(numAtoms_);
-
-            mdAtoms_.nr = numAtoms_;
-
-            for (unsigned i = 0; i < x_.size(); i++)
+            //
+            // All runners should be registered here under appropriate conditions
+            //
+            runners_["LeapFrogSimple"]  = integrateLeapFrogSimple;
+            if (GMX_GPU == GMX_GPU_CUDA && s_hasCompatibleGpus)
             {
-                // Typical PBC box size is tens of nanometers
-                x_[i][XX] = (i%21)*1.0;
-                x_[i][YY] = 6.5 + (i%13)*(-1.0);
-                x_[i][ZZ] = (i%32)*(0.0);
-
-                for (int d = 0; d < DIM; d++)
-                {
-                    xPrime_[i][d] = 0.0;
-                    // Thermal velocity is ~1 nm/ps (|v0| = 1-2 nm/ps)
-                    v_[i][d] = v0[d];
-                    // TODO Check what value typical MD forces have (now ~ 1 kJ/mol/nm)
-                    f_[i][d] = f0[d];
-
-                    x0_[i][d] = x_[i][d];
-                    v0_[i][d] = v_[i][d];
-                }
-                // Atom masses are ~1-100 g/mol
-                inverseMasses_[i] = 1.0/(1.0 + i%100);
+                runners_["LeapFrogGpu"] = integrateLeapFrogGpu;
             }
-            mdAtoms_.invmass = inverseMasses_.data();
         }
 
         //! Store whether any compatible GPUs exist.
@@ -177,71 +121,61 @@ class IntegratorTest : public ::testing::TestWithParam<IntegratorTestParameters>
         }
 };
 
-bool IntegratorTest::s_hasCompatibleGpus = false;
+bool LeapFrogTest::s_hasCompatibleGpus = false;
 
-// The test will run only if:
-// 1. The code was compiled with CUDA
-// 2. There is a CUDA-capable GPU in a system
-// 3. This GPU is detectable
-// 4. GPU detection was not disabled by GMX_DISABLE_GPU_DETECTION environment variable
-TEST_P(IntegratorTest, SimpleIntegration)
+TEST_P(LeapFrogTest, SimpleIntegration)
 {
-    if (!s_hasCompatibleGpus)
-    {
-        return;
-    }
-
-    int  numAtoms;    // 1. Number of atoms
-    real timestep;    // 2. Timestep
-    rvec v0;          // 3. Velocity
-    rvec f0;          // 4. Force
-    int  numSteps;    // 5. Number of steps
-    std::tie(numAtoms, timestep, v0[XX], v0[YY], v0[ZZ], f0[XX], f0[YY], f0[ZZ], numSteps) = GetParam();
-
-    std::string testDescription = formatString("while testing %d atoms for %d timestep (dt = %f, v0=(%f, %f, %f), f0=(%f, %f, %f))",
-                                               numAtoms, numSteps, timestep,
-                                               v0[XX], v0[YY], v0[ZZ],
-                                               f0[XX], f0[YY], f0[ZZ]);
-
-    init(numAtoms, timestep, v0, f0);
-
-    integrateLeapFrogCuda(numAtoms,
-                          (rvec*)(x_.data()),
-                          (rvec*)(xPrime_.data()),
-                          (rvec*)(v_.data()),
-                          (rvec*)(f_.data()),
-                          timestep_,
-                          numSteps,
-                          mdAtoms_);
-
-    real totalTime = numSteps*timestep;
-    // TODO For the case of constant force, the numerical scheme is exact and
-    //      the only source of errors is floating point arithmetic. Hence,
-    //      the tolerance can be calculated.
-    FloatingPointTolerance tolerance = absoluteTolerance(numSteps*0.000005);
-
-    for (unsigned i = 0; i < x_.size(); i++)
+    // Cycle through all available runners
+    for (const auto &runner : runners_)
     {
-        rvec xAnalytical;
-        rvec vAnalytical;
-        for (int d = 0; d < DIM; d++)
+        std::string runnerName = runner.first;
+
+        int         numAtoms; // 1. Number of atoms
+        real        timestep; // 2. Timestep
+        rvec        v0;       // 3. Velocity
+        rvec        f0;       // 4. Force
+        int         numSteps; // 5. Number of steps
+        std::tie(numAtoms, timestep, v0[XX], v0[YY], v0[ZZ], f0[XX], f0[YY], f0[ZZ], numSteps) = GetParam();
+
+        std::string testDescription = formatString("Testing %s with %d atoms for %d timestep (dt = %f, v0=(%f, %f, %f), f0=(%f, %f, %f))",
+                                                   runnerName.c_str(),
+                                                   numAtoms, numSteps, timestep,
+                                                   v0[XX], v0[YY], v0[ZZ],
+                                                   f0[XX], f0[YY], f0[ZZ]);
+        SCOPED_TRACE(testDescription);
+
+        std::unique_ptr<LeapFrogTestData> testData =
+            std::make_unique<LeapFrogTestData>(numAtoms, timestep, v0, f0);
+
+        runner.second(testData.get(), numSteps);
+
+        real totalTime = numSteps*timestep;
+        // TODO For the case of constant force, the numerical scheme is exact and
+        //      the only source of errors is floating point arithmetic. Hence,
+        //      the tolerance can be calculated.
+        FloatingPointTolerance tolerance = absoluteTolerance(numSteps*0.000005);
+
+        for (int i = 0; i < testData->numAtoms_; i++)
         {
-            // Analytical solution for constant-force particle movement
-            xAnalytical[d] = x0_[i][d] + v0_[i][d]*totalTime + 0.5*f_[i][d]*totalTime*totalTime*inverseMasses_[i];
-            vAnalytical[d] = v0_[i][d] + f_[i][d]*totalTime*inverseMasses_[i];
+            rvec xAnalytical;
+            rvec vAnalytical;
+            for (int d = 0; d < DIM; d++)
+            {
+                // Analytical solution for constant-force particle movement
+                xAnalytical[d] = testData->x0_[i][d] + testData->v0_[i][d]*totalTime + 0.5*testData->f_[i][d]*totalTime*totalTime*testData->inverseMasses_[i];
+                vAnalytical[d] = testData->v0_[i][d] + testData->f_[i][d]*totalTime*testData->inverseMasses_[i];
 
-            EXPECT_REAL_EQ_TOL(xPrime_[i][d], xAnalytical[d], tolerance)
-            << gmx::formatString("Coordinate %d of atom %d is different from analytical solution", d, i)
-            << testDescription;
+                EXPECT_REAL_EQ_TOL(xAnalytical[d], testData->xPrime_[i][d], tolerance)
+                << gmx::formatString("Coordinate %d of atom %d is different from analytical solution.", d, i);
 
-            EXPECT_REAL_EQ_TOL(v_[i][d], vAnalytical[d], tolerance)
-            << gmx::formatString("Velocity component %d of atom %d is different from analytical solution", d, i)
-            << testDescription;
+                EXPECT_REAL_EQ_TOL(vAnalytical[d], testData->v_[i][d], tolerance)
+                << gmx::formatString("Velocity component %d of atom %d is different from analytical solution.", d, i);
+            }
         }
     }
 }
 
-INSTANTIATE_TEST_CASE_P(WithParameters, IntegratorTest,
+INSTANTIATE_TEST_CASE_P(WithParameters, LeapFrogTest,
                             ::testing::Combine(
                                     ::testing::Values(1, 10, 300),    // Number of atoms
                                     ::testing::Values(0.001, 0.0005), // Timestep
@@ -253,7 +187,6 @@ INSTANTIATE_TEST_CASE_P(WithParameters, IntegratorTest,
                                     ::testing::Values( 2.0),          // fz
                                     ::testing::Values(1, 10)          // Number of steps
                                 ));
-#endif                                                                // GMX_GPU == GMX_GPU_CUDA
 
 }                                                                     // namespace test
 }                                                                     // namespace gmx
diff --git a/src/gromacs/mdlib/tests/leapfrogtestdata.cpp b/src/gromacs/mdlib/tests/leapfrogtestdata.cpp
new file mode 100644 (file)
index 0000000..38ae157
--- /dev/null
@@ -0,0 +1,174 @@
+/*
+ * 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 Defines the class to accumulate the data needed for the Leap-Frog integrator tests
+ *
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ * \ingroup module_mdlib
+ */
+#include "gmxpre.h"
+
+#include "leapfrogtestdata.h"
+
+#include <assert.h>
+
+#include <cmath>
+
+#include <algorithm>
+#include <unordered_map>
+#include <vector>
+
+#include <gtest/gtest.h>
+
+#include "gromacs/gpu_utils/gpu_utils.h"
+#include "gromacs/math/vec.h"
+#include "gromacs/math/vectypes.h"
+#include "gromacs/mdtypes/mdatom.h"
+#include "gromacs/utility/smalloc.h"
+#include "gromacs/utility/stringutil.h"
+
+#include "testutils/refdata.h"
+#include "testutils/testasserts.h"
+
+namespace gmx
+{
+namespace test
+{
+
+LeapFrogTestData::LeapFrogTestData(int numAtoms, real timestep, const rvec v0, const rvec f0) :
+    numAtoms_(numAtoms),
+    timestep_(timestep),
+    x0_(numAtoms_),
+    x_(numAtoms_),
+    xPrime_(numAtoms_),
+    v0_(numAtoms_),
+    v_(numAtoms_),
+    f_(numAtoms_),
+    inverseMasses_(numAtoms_),
+    inverseMassesPerDim_(numAtoms_)
+{
+    mdAtoms_.nr = numAtoms_;
+
+    for (int i = 0; i < numAtoms_; i++)
+    {
+        // Typical PBC box size is tens of nanometers
+        x_[i][XX] = (i%21)*1.0;
+        x_[i][YY] = 6.5 + (i%13)*(-1.0);
+        x_[i][ZZ] = (i%32)*(0.0);
+
+        for (int d = 0; d < DIM; d++)
+        {
+            xPrime_[i][d] = 0.0;
+            // Thermal velocity is ~1 nm/ps (|v0| = 1-2 nm/ps)
+            v_[i][d] = v0[d];
+            // TODO Check what value typical MD forces have (now ~ 1 kJ/mol/nm)
+            f_[i][d] = f0[d];
+
+            x0_[i][d] = x_[i][d];
+            v0_[i][d] = v_[i][d];
+        }
+        // Atom masses are ~1-100 g/mol
+        inverseMasses_[i] = 1.0/(1.0 + i%100);
+        for (int d = 0; d < DIM; d++)
+        {
+            inverseMassesPerDim_[i][d] = inverseMasses_[i];
+        }
+    }
+    mdAtoms_.invmass       = inverseMasses_.data();
+    mdAtoms_.invMassPerDim = as_rvec_array(inverseMassesPerDim_.data());
+
+    // Data needed for current CPU-based implementation
+    inputRecord_.eI      = eiMD;
+    inputRecord_.delta_t = timestep_;
+    inputRecord_.etc     = etcNO;
+    inputRecord_.epc     = epcNO;
+
+    state_.flags = 0;
+
+    state_.box[XX][XX] = 10.0;
+    state_.box[XX][YY] = 0.0;
+    state_.box[XX][ZZ] = 0.0;
+
+    state_.box[YY][XX] = 0.0;
+    state_.box[YY][YY] = 10.0;
+    state_.box[YY][ZZ] = 0.0;
+
+    state_.box[ZZ][XX] = 0.0;
+    state_.box[ZZ][YY] = 0.0;
+    state_.box[ZZ][ZZ] = 10.0;
+
+    kineticEnergyData_.bNEMD            = false;
+    kineticEnergyData_.cosacc.cos_accel = 0.0;
+    t_grp_tcstat temperatureCouplingGroupData;
+    temperatureCouplingGroupData.lambda = 1;
+    kineticEnergyData_.tcstat.emplace_back(temperatureCouplingGroupData);
+
+    kineticEnergyData_.nthreads = 1;
+    snew(kineticEnergyData_.ekin_work_alloc, kineticEnergyData_.nthreads);
+    snew(kineticEnergyData_.ekin_work, kineticEnergyData_.nthreads);
+    snew(kineticEnergyData_.dekindl_work, kineticEnergyData_.nthreads);
+
+    mdAtoms_.homenr                   = numAtoms_;
+    mdAtoms_.haveVsites               = false;
+    mdAtoms_.havePartiallyFrozenAtoms = false;
+    snew(mdAtoms_.cTC, numAtoms_);
+    for (int i = 0; i < numAtoms_; i++)
+    {
+        mdAtoms_.cTC[i] = 0;
+    }
+
+    prVScalingMatrix_[XX][XX] = 1.0;
+    prVScalingMatrix_[XX][YY] = 0.0;
+    prVScalingMatrix_[XX][ZZ] = 0.0;
+
+    prVScalingMatrix_[YY][XX] = 0.0;
+    prVScalingMatrix_[YY][YY] = 1.0;
+    prVScalingMatrix_[YY][ZZ] = 0.0;
+
+    prVScalingMatrix_[ZZ][XX] = 0.0;
+    prVScalingMatrix_[ZZ][YY] = 0.0;
+    prVScalingMatrix_[ZZ][ZZ] = 1.0;
+
+    update_ = std::make_unique<Update>(&inputRecord_, nullptr);
+    update_->setNumAtoms(numAtoms);
+}
+
+LeapFrogTestData::~LeapFrogTestData()
+{
+    sfree(mdAtoms_.cTC);
+}
+
+}  // namespace test
+}  // namespace gmx
diff --git a/src/gromacs/mdlib/tests/leapfrogtestdata.h b/src/gromacs/mdlib/tests/leapfrogtestdata.h
new file mode 100644 (file)
index 0000000..79bc3df
--- /dev/null
@@ -0,0 +1,125 @@
+/*
+ * 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 Tests for the Leap-Frog integrator
+ *
+ * \todo Prepare for temperature and pressure controlled integrators.
+ * \todo Add PBC handling test.
+ * \todo Reference values tests.
+ *
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ * \ingroup module_mdlib
+ */
+
+#ifndef GMX_MDLIB_TESTS_LEAPFROGTESTDATA_H
+#define GMX_MDLIB_TESTS_LEAPFROGTESTDATA_H
+
+#include <vector>
+
+#include "gromacs/gpu_utils/gpu_utils.h"
+#include "gromacs/math/vec.h"
+#include "gromacs/math/vectypes.h"
+#include "gromacs/mdlib/update.h"
+#include "gromacs/mdtypes/fcdata.h"
+#include "gromacs/mdtypes/group.h"
+#include "gromacs/mdtypes/inputrec.h"
+#include "gromacs/mdtypes/mdatom.h"
+#include "gromacs/mdtypes/state.h"
+#include "gromacs/utility/stringutil.h"
+
+namespace gmx
+{
+namespace test
+{
+
+/* \brief Declares the class to accumulate the data needed for the Leap-Frog integrator tests
+ *
+ */
+class LeapFrogTestData
+{
+    public:
+        //! Number of atoms in the system
+        int                     numAtoms_;
+        //! Integration timestep
+        real                    timestep_;
+
+        //! Initial coordinates
+        PaddedVector<RVec>      x0_;
+        //! Current coordinates
+        PaddedVector<RVec>      x_;
+        //! Coordinates after integrator update
+        PaddedVector<RVec>      xPrime_;
+        //! Initial velocities
+        PaddedVector<RVec>      v0_;
+        //! Current velocities
+        PaddedVector<RVec>      v_;
+        //! External forces
+        PaddedVector<RVec>      f_;
+        //! Inverse masses of the particles
+        PaddedVector<real>      inverseMasses_;
+        //! Inverse masses of the particles per dimension
+        PaddedVector<RVec>      inverseMassesPerDim_;
+
+        //! MD atoms structure in which inverse masses will be passed to the integrator
+        t_mdatoms               mdAtoms_;
+        //! Input record (to get integrator type, temperature and pressure coupling)
+        t_inputrec              inputRecord_;
+        //! System state
+        t_state                 state_;
+        //! Force calculation data
+        t_fcdata                forceCalculationData_;
+        //! Kinetic energy data (to disable non-equilibrium MD integration)
+        gmx_ekindata_t          kineticEnergyData_;
+        //! Parrinnello-Rahman velocity rescalling matrix
+        matrix                  prVScalingMatrix_;
+        //! Update data
+        std::unique_ptr<Update> update_;
+
+        /*! \brief Constructor.
+         *
+         * \param[in]  numAtoms  Number of atoms in the system
+         * \param[in]  timestep  Integration timestep
+         * \param[in]  v0        Initial velocity (same for all particles)
+         * \param[in]  f0        External constant force, acting on all particles
+         */
+        LeapFrogTestData(int numAtoms, real timestep, const rvec v0, const rvec f0);
+
+        ~LeapFrogTestData();
+};
+
+}      // namespace test
+}      // namespace gmx
+
+#endif // GMX_MDLIB_TESTS_LEAPFROGTESTDATA_H
diff --git a/src/gromacs/mdlib/tests/leapfrogtestrunners.cpp b/src/gromacs/mdlib/tests/leapfrogtestrunners.cpp
new file mode 100644 (file)
index 0000000..6fcdff5
--- /dev/null
@@ -0,0 +1,132 @@
+/*
+ * 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 Runner for CPU-based implementation of the integrator.
+ *
+ * \author Artem Zhmurov <zhmurov@gmail.com>
+ * \ingroup module_mdlib
+ */
+#include "gmxpre.h"
+
+#include "leapfrogtestrunners.h"
+
+#include "config.h"
+
+#include <assert.h>
+
+#include <cmath>
+
+#include <algorithm>
+#include <unordered_map>
+#include <vector>
+
+#include "gromacs/math/vec.h"
+#include "gromacs/math/vectypes.h"
+#include "gromacs/mdlib/gmx_omp_nthreads.h"
+#include "gromacs/mdlib/update.h"
+#include "gromacs/mdtypes/group.h"
+#include "gromacs/mdtypes/mdatom.h"
+#include "gromacs/utility/arrayref.h"
+
+#include "testutils/testasserts.h"
+
+#include "leapfrogtestdata.h"
+
+
+namespace gmx
+{
+namespace test
+{
+
+void integrateLeapFrogSimple(LeapFrogTestData *testData,
+                             int               numSteps)
+{
+    testData->state_.x.reserveWithPadding(testData->numAtoms_);
+    testData->state_.v.reserveWithPadding(testData->numAtoms_);
+    for (int i = 0; i < testData->numAtoms_; i++)
+    {
+        testData->state_.x[i] = testData->x_[i];
+        testData->state_.v[i] = testData->v_[i];
+    }
+
+    gmx_omp_nthreads_set(emntUpdate, 1);
+
+    for (int step = 0; step < numSteps; step++)
+    {
+        update_coords(step,
+                      &testData->inputRecord_,
+                      &testData->mdAtoms_,
+                      &testData->state_,
+                      testData->f_,
+                      &testData->forceCalculationData_,
+                      &testData->kineticEnergyData_,
+                      testData->prVScalingMatrix_,
+                      testData->update_.get(),
+                      0,
+                      nullptr,
+                      nullptr);
+        finish_update(&testData->inputRecord_,
+                      &testData->mdAtoms_,
+                      &testData->state_,
+                      nullptr,
+                      nullptr,
+                      nullptr,
+                      testData->update_.get(),
+                      nullptr);
+    }
+    auto xp = makeArrayRef(*testData->update_->xp()).subArray(0, testData->numAtoms_);
+    for (int i = 0; i < testData->numAtoms_; i++)
+    {
+        for (int d = 0; d < DIM; d++)
+        {
+            testData->x_[i][d]      = testData->state_.x[i][d];
+            testData->v_[i][d]      = testData->state_.v[i][d];
+            testData->xPrime_[i][d] = xp[i][d];
+        }
+    }
+}
+
+#if GMX_GPU != GMX_GPU_CUDA
+
+void integrateLeapFrogGpu(gmx_unused  LeapFrogTestData *testData,
+                          gmx_unused  int               numSteps)
+{
+    FAIL() << "Dummy Leap-Frog CUDA function was called instead of the real one.";
+}
+
+#endif // GMX_GPU != GMX_GPU_CUDA
+
+}      // namespace test
+}      // namespace gmx
similarity index 63%
rename from src/gromacs/mdlib/tests/leapfrog.cu
rename to src/gromacs/mdlib/tests/leapfrogtestrunners.cu
index 2281e4945690dea817fe7346ed7d3f284f384a36..7de6892e45c0e930e291e5049e18f4a2efeba7ae 100644 (file)
@@ -33,7 +33,7 @@
  * the research papers on the package. Check out http://www.gromacs.org.
  */
 /*! \internal \file
- * \brief Run integration using CUDA
+ * \brief Runner for CUDA version of the integrator
  *
  * Handles GPU data management and actual numerical integration.
  *
  */
 #include "gmxpre.h"
 
+#include "leapfrogtestrunners.h"
+
+#include "config.h"
+
 #include <assert.h>
 
 #include <cmath>
 #include "gromacs/math/vec.h"
 #include "gromacs/mdlib/leapfrog_cuda.cuh"
 
-#include "leapfrogtest.h"
-
 namespace gmx
 {
 namespace test
 {
 
-void integrateLeapFrogCuda(const int        numAtoms,
-                           const rvec      *h_x,
-                           rvec            *h_xp,
-                           rvec            *h_v,
-                           const rvec      *h_f,
-                           const real       dt,
-                           const int        numSteps,
-                           t_mdatoms        mdAtoms)
+#if GMX_GPU == GMX_GPU_CUDA
+
+void integrateLeapFrogGpu(LeapFrogTestData *testData,
+                          int               numSteps)
 {
+    int     numAtoms = testData->numAtoms_;
+
+    float3 *h_x  = reinterpret_cast<float3*>(testData->x_.data());
+    float3 *h_xp = reinterpret_cast<float3*>(testData->xPrime_.data());
+    float3 *h_v  = reinterpret_cast<float3*>(testData->v_.data());
+    float3 *h_f  = reinterpret_cast<float3*>(testData->f_.data());
+
     float3 *d_x, *d_xp, *d_v, *d_f;
 
     allocateDeviceBuffer(&d_x,  numAtoms, nullptr);
@@ -78,25 +83,25 @@ void integrateLeapFrogCuda(const int        numAtoms,
     allocateDeviceBuffer(&d_v,  numAtoms, nullptr);
     allocateDeviceBuffer(&d_f,  numAtoms, nullptr);
 
-    copyToDeviceBuffer(&d_x,  (float3*)h_x,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
-    copyToDeviceBuffer(&d_xp, (float3*)h_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
-    copyToDeviceBuffer(&d_v,  (float3*)h_v,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
-    copyToDeviceBuffer(&d_f,  (float3*)h_f,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+    copyToDeviceBuffer(&d_x,  h_x,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+    copyToDeviceBuffer(&d_xp, h_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+    copyToDeviceBuffer(&d_v,  h_v,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+    copyToDeviceBuffer(&d_f,  h_f,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
 
     auto integrator = std::make_unique<LeapFrogCuda>();
 
-    integrator->set(mdAtoms);
+    integrator->set(testData->mdAtoms_);
 
     for (int step = 0; step < numSteps; step++)
     {
-        integrator->integrate(d_x, d_xp, d_v, d_f, dt);
+        integrator->integrate(d_x, d_xp, d_v, d_f, testData->timestep_);
 
-        copyFromDeviceBuffer((float3*)h_xp, &d_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
-        copyToDeviceBuffer(&d_x,    (float3*)h_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+        copyFromDeviceBuffer(h_xp, &d_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+        copyToDeviceBuffer(&d_x,    h_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
     }
 
-    copyFromDeviceBuffer((float3*)h_xp, &d_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
-    copyFromDeviceBuffer((float3*)h_v,  &d_v,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+    copyFromDeviceBuffer(h_xp, &d_xp, 0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
+    copyFromDeviceBuffer(h_v,  &d_v,  0, numAtoms, nullptr, GpuApiCallBehavior::Sync, nullptr);
 
     freeDeviceBuffer(&d_x);
     freeDeviceBuffer(&d_xp);
@@ -104,5 +109,7 @@ void integrateLeapFrogCuda(const int        numAtoms,
     freeDeviceBuffer(&d_f);
 }
 
-}   // namespace test
-}   // namespace gmx
+#endif // GMX_GPU == GMX_GPU_CUDA
+
+}      // namespace test
+}      // namespace gmx
similarity index 68%
rename from src/gromacs/mdlib/tests/leapfrogtest.h
rename to src/gromacs/mdlib/tests/leapfrogtestrunners.h
index ca2124137f604d2d2e74261d58b2443e03e3c9a2..9ef313640fa459e2a7d97ec05dab623d0fbecf3d 100644 (file)
  * the research papers on the package. Check out http://www.gromacs.org.
  */
 /*! \internal \file
- * \brief Declaration of interfaces to run integrator on various devices
+ * \brief Declaration of interfaces to run various implementations of integrator (runners)
  *
  * \author Artem Zhmurov <zhmurov@gmail.com>
  * \ingroup module_mdlib
  */
 
-#include "config.h"
+#ifndef GMX_MDLIB_TESTS_LEAPFROGTESTRUNNERS_H
+#define GMX_MDLIB_TESTS_LEAPFROGTESTRUNNERS_H
 
 #include "gromacs/math/vec.h"
 
-struct t_mdatoms;
+#include "leapfrogtestdata.h"
 
 namespace gmx
 {
 namespace test
 {
 
-#if GMX_GPU == GMX_GPU_CUDA
+/*! \brief Integrate using CPU version of Leap-Frog
+ *
+ * \param[in]     testData  Data needed for the integrator
+ * \param[in]     numSteps  Total number of steps to run integration for.
+ */
+void integrateLeapFrogSimple(LeapFrogTestData *testData,
+                             int               numSteps);
 
 /*! \brief Integrate using CUDA version of Leap-Frog
  *
@@ -58,25 +65,13 @@ namespace test
  * for requested number of steps using Leap-Frog algorithm, copies
  * the result back.
  *
- * \param[in]     numAtoms  Number of atoms.
- * \param[in]     h_x       Initial coordinates.
- * \param[out]    h_xp      Place to save the resulting coordinates to.
- * \param[in,out] h_v       Velocities (will be updated).
- * \param[in]     h_f       Forces.
- * \param[in]     dt        Timestep.
- * \param[in]     numSteps  Total number of steps to runtegration for.
- * \param[in]     mdAtoms   Infromation on atoms (atom masses are used here).
+ * \param[in]     testData  Data needed for the integrator
+ * \param[in]     numSteps  Total number of steps to run integration for.
  */
-void integrateLeapFrogCuda(int         numAtoms,
-                           const rvec *h_x,
-                           rvec       *h_xp,
-                           rvec       *h_v,
-                           const rvec *h_f,
-                           real        dt,
-                           int         numSteps,
-                           t_mdatoms   mdAtoms);
-
-#endif // GMX_GPU == GMX_GPU_CUDA
+void integrateLeapFrogGpu(LeapFrogTestData *testData,
+                          int               numSteps);
 
 }      // namespace test
 }      // namespace gmx
+
+#endif // GMX_MDLIB_TESTS_LEAPFROGTESTRUNNERS_H