PME force gathering - CUDA kernel + unit tests
[gromacs.git] / src / gromacs / ewald / tests / pmegathertest.cpp
blob404e6c57c8c439313c65d4c18ca49449971a1eb0
1 /*
2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2016,2017, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
35 /*! \internal \file
36 * \brief
37 * Implements PME force gathering tests.
39 * \author Aleksei Iupinov <a.yupinov@gmail.com>
40 * \ingroup module_ewald
43 #include "gmxpre.h"
45 #include <string>
47 #include <gmock/gmock.h>
49 #include "gromacs/mdtypes/inputrec.h"
50 #include "gromacs/utility/stringutil.h"
52 #include "testutils/refdata.h"
53 #include "testutils/testasserts.h"
55 #include "pmetestcommon.h"
57 namespace gmx
59 namespace test
61 namespace
64 /* Valid input instances */
66 //! A couple of valid inputs for boxes.
67 static std::vector<Matrix3x3> const c_sampleBoxes
69 // normal box
70 Matrix3x3 {{
71 8.0f, 0.0f, 0.0f,
72 0.0f, 3.4f, 0.0f,
73 0.0f, 0.0f, 2.0f
74 }},
75 // triclinic box
76 Matrix3x3 {{
77 7.0f, 0.0f, 0.0f,
78 0.0f, 4.1f, 0.0f,
79 3.5f, 2.0f, 12.2f
80 }},
83 //! A couple of valid inputs for grid sizes
84 static std::vector<IVec> const c_sampleGridSizes
86 IVec {
87 16, 12, 14
89 IVec {
90 13, 15, 11
93 //! Random charges
94 static std::vector<real> const c_sampleChargesFull
96 4.95f, 3.11f, 3.97f, 1.08f, 2.09f, 1.1f, 4.13f, 3.31f, 2.8f, 5.83f, 5.09f, 6.1f, 2.86f, 0.24f, 5.76f, 5.19f, 0.72f
99 //! All the input atom gridline indices
100 static std::vector<IVec> const c_sampleGridLineIndicesFull
102 IVec {
103 4, 2, 6
105 IVec {
106 1, 4, 10
108 IVec {
109 0, 6, 6
111 IVec {
112 0, 1, 4
114 IVec {
115 6, 3, 0
117 IVec {
118 7, 2, 2
120 IVec {
121 8, 3, 1
123 IVec {
124 4, 0, 3
126 IVec {
127 0, 0, 0
129 IVec {
130 8, 5, 8
132 IVec {
133 4, 4, 2
135 IVec {
136 7, 1, 7
138 IVec {
139 8, 5, 5
141 IVec {
142 2, 6, 5
144 IVec {
145 1, 6, 2
147 IVec {
148 7, 1, 8
150 IVec {
151 3, 5, 1
155 // Spline values/derivatives below are also generated randomly, so they are bogus,
156 // but that should not affect the reproducibility, which we're after
158 //! A lot of bogus input spline values - should have at list (max PME order = 5) * (DIM = 3) * (total unique atom number in all test cases = 16) values
159 static std::vector<real> const c_sampleSplineValuesFull
161 0.12f, 0.81f, 0.29f, 0.22f, 0.13f, 0.19f, 0.12f, 0.8f, 0.44f, 0.38f, 0.32f, 0.36f, 0.27f, 0.11f, 0.17f, 0.94f, 0.07f, 0.9f, 0.98f, 0.96f, 0.07f, 0.94f, 0.77f, 0.24f, 0.84f, 0.16f, 0.77f, 0.57f, 0.52f, 0.27f, 0.39f, 0.45f, 0.6f, 0.59f, 0.44f, 0.91f, 0.97f, 0.43f, 0.24f, 0.52f, 0.73f, 0.55f, 0.99f, 0.39f, 0.97f, 0.35f, 0.1f, 0.68f, 0.19f, 0.1f, 0.77f, 0.2f, 0.43f, 0.69f, 0.76f, 0.32f, 0.31f, 0.94f, 0.53f, 0.6f, 0.93f, 0.57f, 0.94f, 0.88f, 0.75f, 0.77f, 0.91f, 0.72f, 0.07f, 0.78f, 0.09f, 0.02f, 0.48f, 0.97f, 0.89f, 0.39f, 0.48f, 0.19f, 0.02f, 0.92f, 0.8f, 0.41f, 0.53f, 0.32f, 0.38f, 0.58f, 0.36f, 0.46f, 0.92f, 0.91f, 0.01f, 0.86f, 0.54f, 0.86f, 0.94f, 0.37f, 0.35f, 0.81f, 0.89f, 0.48f,
162 0.34f, 0.18f, 0.11f, 0.02f, 0.87f, 0.95f, 0.66f, 0.67f, 0.38f, 0.45f, 0.04f, 0.94f, 0.54f, 0.76f, 0.58f, 0.83f, 0.31f, 0.73f, 0.71f, 0.06f, 0.35f, 0.32f, 0.35f, 0.61f, 0.27f, 0.98f, 0.83f, 0.11f, 0.3f, 0.42f, 0.95f, 0.69f, 0.58f, 0.29f, 0.1f, 0.68f, 0.94f, 0.62f, 0.51f, 0.47f, 0.04f, 0.47f, 0.34f, 0.71f, 0.52f, 0.19f, 0.69f, 0.5f, 0.59f, 0.05f, 0.74f, 0.11f, 0.4f, 0.81f, 0.24f, 0.53f, 0.71f, 0.07f, 0.17f, 0.41f, 0.23f, 0.78f, 0.27f, 0.1f, 0.71f, 0.36f, 0.67f, 0.6f, 0.94f, 0.69f, 0.19f, 0.58f, 0.68f, 0.5f, 0.62f, 0.38f, 0.29f, 0.44f, 0.04f, 0.89f, 0.0f, 0.76f, 0.22f, 0.16f, 0.08f, 0.62f, 0.51f, 0.62f, 0.83f, 0.72f, 0.96f, 0.99f, 0.4f, 0.79f, 0.83f, 0.21f, 0.43f, 0.32f, 0.44f, 0.72f,
163 0.21f, 0.4f, 0.93f, 0.07f, 0.11f, 0.41f, 0.24f, 0.04f, 0.36f, 0.15f, 0.92f, 0.08f, 0.99f, 0.35f, 0.42f, 0.7f, 0.17f, 0.39f, 0.69f, 0.0f, 0.86f, 0.89f, 0.59f, 0.81f, 0.77f, 0.15f, 0.89f, 0.17f, 0.76f, 0.67f, 0.58f, 0.78f, 0.26f, 0.19f, 0.69f, 0.18f, 0.46f, 0.6f, 0.69f, 0.23f, 0.34f, 0.3f, 0.64f, 0.34f, 0.6f, 0.99f, 0.69f, 0.57f, 0.75f, 0.07f, 0.36f, 0.75f, 0.81f, 0.8f, 0.42f, 0.09f, 0.94f, 0.66f, 0.35f, 0.67f, 0.34f, 0.66f, 0.02f, 0.47f, 0.78f, 0.21f, 0.02f, 0.18f, 0.42f, 0.2f, 0.46f, 0.34f, 0.4f, 0.46f, 0.96f, 0.86f, 0.25f, 0.25f, 0.22f, 0.37f, 0.59f, 0.19f, 0.45f, 0.61f, 0.04f, 0.71f, 0.77f, 0.51f, 0.77f, 0.15f, 0.78f, 0.36f, 0.62f, 0.24f, 0.86f, 0.2f, 0.77f, 0.08f, 0.09f, 0.3f,
164 0.0f, 0.6f, 0.99f, 0.69f,
167 //! A lot of bogus input spline derivatives - should have at list (max PME order = 5) * (DIM = 3) * (total unique atom number in all test cases = 16) values
168 static std::vector<real> const c_sampleSplineDerivativesFull
170 0.82f, 0.88f, 0.83f, 0.11f, 0.93f, 0.32f, 0.71f, 0.37f, 0.69f, 0.88f, 0.11f, 0.38f, 0.25f, 0.5f, 0.36f, 0.81f, 0.78f, 0.31f, 0.66f, 0.32f, 0.27f, 0.35f, 0.53f, 0.83f, 0.08f, 0.08f, 0.94f, 0.71f, 0.65f, 0.24f, 0.13f, 0.01f, 0.33f, 0.65f, 0.24f, 0.53f, 0.45f, 0.84f, 0.33f, 0.97f, 0.31f, 0.7f, 0.03f, 0.31f, 0.41f, 0.76f, 0.12f, 0.3f, 0.57f, 0.65f, 0.87f, 0.99f, 0.42f, 0.97f, 0.32f, 0.39f, 0.73f, 0.23f, 0.03f, 0.67f, 0.97f, 0.57f, 0.42f, 0.38f, 0.54f, 0.17f, 0.53f, 0.54f, 0.18f, 0.8f, 0.76f, 0.13f, 0.29f, 0.83f, 0.77f, 0.56f, 0.4f, 0.87f, 0.36f, 0.18f, 0.59f, 0.04f, 0.05f, 0.61f, 0.26f, 0.91f, 0.62f, 0.16f, 0.89f, 0.23f, 0.26f, 0.59f, 0.33f, 0.2f, 0.49f, 0.41f, 0.25f, 0.4f, 0.16f, 0.83f,
171 0.44f, 0.82f, 0.21f, 0.95f, 0.14f, 0.8f, 0.37f, 0.31f, 0.41f, 0.53f, 0.15f, 0.85f, 0.78f, 0.17f, 0.92f, 0.03f, 0.13f, 0.2f, 0.03f, 0.33f, 0.87f, 0.38f, 0, 0.08f, 0.79f, 0.36f, 0.53f, 0.05f, 0.07f, 0.94f, 0.23f, 0.85f, 0.13f, 0.27f, 0.23f, 0.22f, 0.26f, 0.38f, 0.15f, 0.48f, 0.18f, 0.33f, 0.23f, 0.62f, 0.1f, 0.36f, 0.99f, 0.07f, 0.02f, 0.04f, 0.09f, 0.29f, 0.52f, 0.29f, 0.83f, 0.97f, 0.61f, 0.81f, 0.49f, 0.56f, 0.08f, 0.09f, 0.03f, 0.65f, 0.46f, 0.1f, 0.06f, 0.06f, 0.39f, 0.29f, 0.04f, 0.03f, 0.1f, 0.83f, 0.94f, 0.59f, 0.97f, 0.82f, 0.2f, 0.66f, 0.23f, 0.11f, 0.03f, 0.16f, 0.27f, 0.53f, 0.94f, 0.46f, 0.43f, 0.29f, 0.97f, 0.64f, 0.46f, 0.37f, 0.43f, 0.48f, 0.37f, 0.93f, 0.5f, 0.2f,
172 0.92f, 0.09f, 0.74f, 0.55f, 0.44f, 0.05f, 0.13f, 0.17f, 0.79f, 0.44f, 0.11f, 0.6f, 0.64f, 0.05f, 0.96f, 0.3f, 0.45f, 0.47f, 0.42f, 0.74f, 0.91f, 0.06f, 0.89f, 0.24f, 0.26f, 0.68f, 0.4f, 0.88f, 0.5f, 0.65f, 0.48f, 0.15f, 0.0f, 0.41f, 0.67f, 0.4f, 0.31f, 0.73f, 0.77f, 0.36f, 0.26f, 0.74f, 0.46f, 0.56f, 0.78f, 0.92f, 0.32f, 0.9f, 0.06f, 0.55f, 0.6f, 0.13f, 0.38f, 0.93f, 0.5f, 0.92f, 0.96f, 0.82f, 0.0f, 0.04f, 0.9f, 0.55f, 0.97f, 1.0f, 0.23f, 0.46f, 0.52f, 0.49f, 0.0f, 0.32f, 0.16f, 0.4f, 0.62f, 0.36f, 0.03f, 0.63f, 0.16f, 0.58f, 0.97f, 0.03f, 0.44f, 0.07f, 0.22f, 0.75f, 0.32f, 0.61f, 0.94f, 0.33f, 0.7f, 0.57f, 0.5f, 0.84f, 0.7f, 0.47f, 0.18f, 0.09f, 0.25f, 0.77f, 0.94f, 0.85f,
173 0.09f, 0.83f, 0.02f, 0.91f,
176 //! 2 c_sample grids - only non-zero values have to be listed
177 static std::vector<SparseRealGridValuesInput> const c_sampleGrids
179 SparseRealGridValuesInput {{
180 IVec {
181 0, 0, 0
182 }, 3.5f
183 }, {
184 IVec {
185 7, 0, 0
186 }, -2.5f
189 IVec {
190 3, 5, 7
191 }, -0.006f
192 }, {
193 IVec {
194 1, 6, 7
195 }, -2.76f
196 }, {
197 IVec {
198 3, 1, 2
199 }, 0.6f
200 }, {
201 IVec {
202 6, 2, 4
203 }, 7.1f
204 }, {
205 IVec {
206 4, 5, 6
207 }, 4.1f
208 }, {
209 IVec {
210 4, 4, 6
211 }, -3.7f
212 }, },
213 SparseRealGridValuesInput {{
214 IVec {
215 0, 4, 0
216 }, 6.f
217 }, {
218 IVec {
219 4, 2, 7
220 }, 13.76f
221 }, {
222 IVec {
223 0, 6, 7
224 }, 3.6f
225 }, {
226 IVec {
227 3, 2, 8
228 }, 0.61f
231 IVec {
232 5, 4, 3
233 }, 2.1f
236 IVec {
237 2, 5, 10
238 }, 3.6f
239 }, {
240 IVec {
241 5, 3, 6
242 }, 2.1f
243 }, {
244 IVec {
245 6, 6, 6
246 }, 6.1f
247 }, }
250 //! Input forces for reduction
251 static std::vector<RVec> const c_sampleForcesFull {
252 RVec {
253 0.02f, 0.87f, 0.95f
254 }, RVec {
255 0.66f, 0.67f, 0.38f
256 }, RVec {
257 0.45f, 0.04f, 0.94f
258 }, RVec {
259 0.54f, 0.76f, 0.58f
261 RVec {
262 0.83f, 0.31f, 0.73f
263 }, RVec {
264 0.71f, 0.06f, 0.35f
265 }, RVec {
266 0.32f, 0.35f, 0.61f
267 }, RVec {
268 0.27f, 0.98f, 0.83f
270 RVec {
271 0.11f, 0.3f, 0.42f
272 }, RVec {
273 0.95f, 0.69f, 0.58f
274 }, RVec {
275 0.29f, 0.1f, 0.68f
276 }, RVec {
277 0.94f, 0.62f, 0.51f
279 RVec {
280 0.47f, 0.04f, 0.47f
281 }, RVec {
282 0.34f, 0.71f, 0.52f
286 //! PME orders to test
287 static std::vector<int> const pmeOrders {
288 3, 4, 5
290 //! Atom counts to test
291 static std::vector<size_t> const atomCounts {
292 1, 2, 13
295 /* Helper structures for test input */
297 //! A structure for all the spline data which depends in size both on the PME order and atom count
298 struct AtomAndPmeOrderSizedData
300 //! Spline values
301 SplineParamsVector splineValues;
302 //! Spline derivatives
303 SplineParamsVector splineDerivatives;
306 //! A structure for all the input atom data which depends in size on atom count - including a range of spline data for different PME orders
307 struct AtomSizedData
309 //! Gridline indices
310 GridLineIndicesVector gridLineIndices;
311 //! Charges
312 ChargesVector charges;
313 //! Coordinates
314 CoordinatesVector coordinates;
315 //! Spline data for different orders
316 std::map<int, AtomAndPmeOrderSizedData> splineDataByPmeOrder;
319 //! A range of the test input data sets, uniquely identified by the atom count
320 typedef std::map<size_t, AtomSizedData> InputDataByAtomCount;
322 /*! \brief Convenience typedef of the test input parameters - unit cell box, PME interpolation order, grid dimensions,
323 * grid values, overwriting/reducing the input forces, atom count.
325 * The rest of the atom-related input data - gridline indices, spline theta values, spline dtheta values, atom charges -
326 * is looked up in the inputAtomDataSets_ test fixture variable.
328 typedef std::tuple<Matrix3x3, int, IVec, SparseRealGridValuesInput, PmeGatherInputHandling, size_t> GatherInputParameters;
330 //! Test fixture
331 class PmeGatherTest : public ::testing::TestWithParam<GatherInputParameters>
333 private:
334 //! Storage of all the input atom datasets
335 static InputDataByAtomCount s_inputAtomDataSets_;
337 public:
338 //! Default constructor
339 PmeGatherTest() = default;
340 //! Sets the input atom data references once
341 static void SetUpTestCase()
343 auto gridLineIndicesIt = c_sampleGridLineIndicesFull.begin();
344 auto chargesIt = c_sampleChargesFull.begin();
345 for (auto atomCount : atomCounts)
347 AtomSizedData atomData;
348 atomData.gridLineIndices = GridLineIndicesVector::fromVector(gridLineIndicesIt, gridLineIndicesIt + atomCount);
349 gridLineIndicesIt += atomCount;
350 atomData.charges = ChargesVector::fromVector(chargesIt, chargesIt + atomCount);
351 chargesIt += atomCount;
352 atomData.coordinates.resize(atomCount, RVec {1e6, 1e7, -1e8});
353 /* The coordinates are intentionally bogus in this test - only the size matters; the gridline indices are fed directly as inputs */
354 for (auto pmeOrder : pmeOrders)
356 AtomAndPmeOrderSizedData splineData;
357 const size_t dimSize = atomCount * pmeOrder;
358 for (int dimIndex = 0; dimIndex < DIM; dimIndex++)
360 splineData.splineValues[dimIndex] = SplineParamsDimVector::fromVector(c_sampleSplineValuesFull.begin() + dimIndex * dimSize,
361 c_sampleSplineValuesFull.begin() + (dimIndex + 1) * dimSize);
362 splineData.splineDerivatives[dimIndex] = SplineParamsDimVector::fromVector(c_sampleSplineDerivativesFull.begin() + dimIndex * dimSize,
363 c_sampleSplineDerivativesFull.begin() + (dimIndex + 1) * dimSize);
365 atomData.splineDataByPmeOrder[pmeOrder] = splineData;
367 s_inputAtomDataSets_[atomCount] = atomData;
370 //! The test
371 void runTest()
373 /* Getting the input */
374 Matrix3x3 box;
375 int pmeOrder;
376 IVec gridSize;
377 size_t atomCount;
378 SparseRealGridValuesInput nonZeroGridValues;
379 PmeGatherInputHandling inputForceTreatment;
380 std::tie(box, pmeOrder, gridSize, nonZeroGridValues, inputForceTreatment, atomCount) = GetParam();
381 auto inputAtomData = s_inputAtomDataSets_[atomCount];
382 auto inputAtomSplineData = inputAtomData.splineDataByPmeOrder[pmeOrder];
384 /* Storing the input where it's needed, running the test */
385 t_inputrec inputRec;
386 inputRec.nkx = gridSize[XX];
387 inputRec.nky = gridSize[YY];
388 inputRec.nkz = gridSize[ZZ];
389 inputRec.pme_order = pmeOrder;
390 inputRec.coulombtype = eelPME;
391 inputRec.epsilon_r = 1.0;
393 TestReferenceData refData;
394 const std::map<CodePath, std::string> modesToTest = {{CodePath::CPU, "CPU"},
395 {CodePath::CUDA, "CUDA"}};
396 for (const auto &mode : modesToTest)
398 const bool supportedInput = pmeSupportsInputForMode(&inputRec, mode.first);
399 if (!supportedInput)
401 /* Testing the failure for the unsupported input */
402 EXPECT_THROW(pmeInitAtoms(&inputRec, mode.first, nullptr, inputAtomData.coordinates, inputAtomData.charges, box), NotImplementedError);
403 continue;
406 const auto contextsToTest = pmeEnv->getHardwareContexts(mode.first);
407 for (const auto &context : contextsToTest)
409 /* Describing the test uniquely */
410 SCOPED_TRACE(formatString("Testing force gathering with %s %sfor PME grid size %d %d %d"
411 ", order %d, %zu atoms, %s",
412 mode.second.c_str(), context.getDescription().c_str(),
413 gridSize[XX], gridSize[YY], gridSize[ZZ],
414 pmeOrder,
415 atomCount,
416 (inputForceTreatment == PmeGatherInputHandling::ReduceWith) ? "with reduction" : "without reduction"
419 PmeSafePointer pmeSafe = pmeInitAtoms(&inputRec, mode.first, context.getDeviceInfo(), inputAtomData.coordinates, inputAtomData.charges, box);
421 /* Setting some more inputs */
422 pmeSetRealGrid(pmeSafe.get(), mode.first, nonZeroGridValues);
423 pmeSetGridLineIndices(pmeSafe.get(), mode.first, inputAtomData.gridLineIndices);
424 for (int dimIndex = 0; dimIndex < DIM; dimIndex++)
426 pmeSetSplineData(pmeSafe.get(), mode.first, inputAtomSplineData.splineValues[dimIndex], PmeSplineDataType::Values, dimIndex);
427 pmeSetSplineData(pmeSafe.get(), mode.first, inputAtomSplineData.splineDerivatives[dimIndex], PmeSplineDataType::Derivatives, dimIndex);
430 /* Explicitly copying the sample forces to be able to modify them */
431 auto inputForcesFull(c_sampleForcesFull);
432 GMX_RELEASE_ASSERT(inputForcesFull.size() >= atomCount, "Bad input forces size");
433 auto forces = ForcesVector::fromVector(inputForcesFull.begin(), inputForcesFull.begin() + atomCount);
435 /* Running the force gathering itself */
436 pmePerformGather(pmeSafe.get(), mode.first, inputForceTreatment, forces);
437 pmeFinalizeTest(pmeSafe.get(), mode.first);
439 /* Check the output forces correctness */
440 TestReferenceChecker forceChecker(refData.rootChecker());
441 const auto ulpTolerance = 3 * pmeOrder;
442 forceChecker.setDefaultTolerance(relativeToleranceAsUlp(1.0, ulpTolerance));
443 forceChecker.checkSequence(forces.begin(), forces.end(), "Forces");
449 // An instance of static atom data
450 InputDataByAtomCount PmeGatherTest::s_inputAtomDataSets_;
452 //! Test for PME force gathering
453 TEST_P(PmeGatherTest, ReproducesOutputs)
455 EXPECT_NO_THROW(runTest());
458 //! Instantiation of the PME gathering test
459 INSTANTIATE_TEST_CASE_P(SaneInput, PmeGatherTest, ::testing::Combine(::testing::ValuesIn(c_sampleBoxes),
460 ::testing::ValuesIn(pmeOrders),
461 ::testing::ValuesIn(c_sampleGridSizes),
462 ::testing::ValuesIn(c_sampleGrids),
463 ::testing::Values(PmeGatherInputHandling::Overwrite, PmeGatherInputHandling::ReduceWith),
464 ::testing::ValuesIn(atomCounts)));