2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2016,2017,2018,2019,2020,2021, 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.
37 * Implements PME force gathering tests.
39 * \author Aleksei Iupinov <a.yupinov@gmail.com>
40 * \ingroup module_ewald
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/test_hardware_environment.h"
54 #include "testutils/testasserts.h"
56 #include "pmetestcommon.h"
65 /* Valid input instances */
67 //! A couple of valid inputs for boxes.
68 std::vector<Matrix3x3> const c_sampleBoxes{
70 Matrix3x3{ { 8.0F, 0.0F, 0.0F, 0.0F, 3.4F, 0.0F, 0.0F, 0.0F, 2.0F } },
72 Matrix3x3{ { 7.0F, 0.0F, 0.0F, 0.0F, 4.1F, 0.0F, 3.5F, 2.0F, 12.2F } },
75 //! A couple of valid inputs for grid sizes
76 std::vector<IVec> const c_sampleGridSizes{ IVec{ 16, 12, 14 }, IVec{ 13, 15, 11 } };
78 std::vector<real> const c_sampleChargesFull{ 4.95F, 3.11F, 3.97F, 1.08F, 2.09F, 1.1F,
79 4.13F, 3.31F, 2.8F, 5.83F, 5.09F, 6.1F,
80 2.86F, 0.24F, 5.76F, 5.19F, 0.72F };
82 //! All the input atom gridline indices
83 std::vector<IVec> const c_sampleGridLineIndicesFull{
84 IVec{ 4, 2, 6 }, IVec{ 1, 4, 10 }, IVec{ 0, 6, 6 }, IVec{ 0, 1, 4 }, IVec{ 6, 3, 0 },
85 IVec{ 7, 2, 2 }, IVec{ 8, 3, 1 }, IVec{ 4, 0, 3 }, IVec{ 0, 0, 0 }, IVec{ 8, 5, 8 },
86 IVec{ 4, 4, 2 }, IVec{ 7, 1, 7 }, IVec{ 8, 5, 5 }, IVec{ 2, 6, 5 }, IVec{ 1, 6, 2 },
87 IVec{ 7, 1, 8 }, IVec{ 3, 5, 1 },
90 // Spline values/derivatives below are also generated randomly, so they are bogus,
91 // but that should not affect the reproducibility, which we're after
93 //! 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
94 std::vector<real> const c_sampleSplineValuesFull{
95 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,
96 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,
97 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,
98 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,
99 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,
100 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,
101 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,
102 0.86F, 0.54F, 0.86F, 0.94F, 0.37F, 0.35F, 0.81F, 0.89F, 0.48F, 0.34F, 0.18F, 0.11F, 0.02F,
103 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,
104 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,
105 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,
106 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,
107 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,
108 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,
109 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,
110 0.21F, 0.43F, 0.32F, 0.44F, 0.72F, 0.21F, 0.4F, 0.93F, 0.07F, 0.11F, 0.41F, 0.24F, 0.04F,
111 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,
112 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,
113 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,
114 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,
115 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,
116 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,
117 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,
118 0.3F, 0.0F, 0.6F, 0.99F, 0.69F,
121 //! 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
122 std::vector<real> const c_sampleSplineDerivativesFull{
123 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,
124 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,
125 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,
126 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,
127 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,
128 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,
129 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,
130 0.59F, 0.33F, 0.2F, 0.49F, 0.41F, 0.25F, 0.4F, 0.16F, 0.83F, 0.44F, 0.82F, 0.21F, 0.95F,
131 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,
132 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,
133 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,
134 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,
135 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,
136 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,
137 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,
138 0.48F, 0.37F, 0.93F, 0.5F, 0.2F, 0.92F, 0.09F, 0.74F, 0.55F, 0.44F, 0.05F, 0.13F, 0.17F,
139 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,
140 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,
141 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,
142 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,
143 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,
144 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,
145 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,
146 0.85F, 0.09F, 0.83F, 0.02F, 0.91F,
149 //! 2 c_sample grids - only non-zero values have to be listed
150 std::vector<SparseRealGridValuesInput> const c_sampleGrids{ SparseRealGridValuesInput{
151 { IVec{ 0, 0, 0 }, 3.5F },
152 { IVec{ 7, 0, 0 }, -2.5F },
153 { IVec{ 3, 5, 7 }, -0.006F },
154 { IVec{ 1, 6, 7 }, -2.76F },
155 { IVec{ 3, 1, 2 }, 0.6F },
156 { IVec{ 6, 2, 4 }, 7.1F },
157 { IVec{ 4, 5, 6 }, 4.1F },
158 { IVec{ 4, 4, 6 }, -3.7F },
160 SparseRealGridValuesInput{
161 { IVec{ 0, 4, 0 }, 6.F },
162 { IVec{ 4, 2, 7 }, 13.76F },
163 { IVec{ 0, 6, 7 }, 3.6F },
164 { IVec{ 3, 2, 8 }, 0.61F },
165 { IVec{ 5, 4, 3 }, 2.1F },
166 { IVec{ 2, 5, 10 }, 3.6F },
167 { IVec{ 5, 3, 6 }, 2.1F },
168 { IVec{ 6, 6, 6 }, 6.1F },
171 //! Input forces for reduction
172 std::vector<RVec> const c_sampleForcesFull{
173 RVec{ 0.02F, 0.87F, 0.95F }, RVec{ 0.66F, 0.67F, 0.38F }, RVec{ 0.45F, 0.04F, 0.94F },
174 RVec{ 0.54F, 0.76F, 0.58F }, RVec{ 0.83F, 0.31F, 0.73F }, RVec{ 0.71F, 0.06F, 0.35F },
175 RVec{ 0.32F, 0.35F, 0.61F }, RVec{ 0.27F, 0.98F, 0.83F }, RVec{ 0.11F, 0.3F, 0.42F },
176 RVec{ 0.95F, 0.69F, 0.58F }, RVec{ 0.29F, 0.1F, 0.68F }, RVec{ 0.94F, 0.62F, 0.51F },
177 RVec{ 0.47F, 0.04F, 0.47F }, RVec{ 0.34F, 0.71F, 0.52F }
180 //! PME orders to test
181 std::vector<int> const pmeOrders{ 3, 4, 5 };
182 //! Atom counts to test
183 std::vector<size_t> const atomCounts{ 1, 2, 13 };
185 /* Helper structures for test input */
187 //! A structure for all the spline data which depends in size both on the PME order and atom count
188 struct AtomAndPmeOrderSizedData
191 SplineParamsVector splineValues;
192 //! Spline derivatives
193 SplineParamsVector splineDerivatives;
196 //! 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
200 GridLineIndicesVector gridLineIndices;
202 ChargesVector charges;
204 CoordinatesVector coordinates;
205 //! Spline data for different orders
206 std::map<int, AtomAndPmeOrderSizedData> splineDataByPmeOrder;
209 //! A range of the test input data sets, uniquely identified by the atom count
210 typedef std::map<size_t, AtomSizedData> InputDataByAtomCount;
212 /*! \brief Convenience typedef of the test input parameters - unit cell box, PME interpolation
213 * order, grid dimensions, grid values, overwriting/reducing the input forces, atom count.
215 * The rest of the atom-related input data - gridline indices, spline theta values, spline dtheta
216 * values, atom charges - is looked up in the inputAtomDataSets_ test fixture variable.
218 typedef std::tuple<Matrix3x3, int, IVec, SparseRealGridValuesInput, size_t> GatherInputParameters;
221 class PmeGatherTest : public ::testing::TestWithParam<GatherInputParameters>
224 //! Storage of all the input atom datasets
225 static InputDataByAtomCount s_inputAtomDataSets_;
228 PmeGatherTest() = default;
229 //! Sets the input atom data references and programs once
230 static void SetUpTestSuite()
233 for (auto atomCount : atomCounts)
235 AtomSizedData atomData;
236 atomData.gridLineIndices =
237 GridLineIndicesVector(c_sampleGridLineIndicesFull).subArray(start, atomCount);
238 atomData.charges = ChargesVector(c_sampleChargesFull).subArray(start, atomCount);
240 atomData.coordinates.resize(atomCount, RVec{ 1e6, 1e7, -1e8 });
241 /* The coordinates are intentionally bogus in this test - only the size matters; the gridline indices are fed directly as inputs */
242 for (auto pmeOrder : pmeOrders)
244 AtomAndPmeOrderSizedData splineData;
245 const size_t dimSize = atomCount * pmeOrder;
246 for (int dimIndex = 0; dimIndex < DIM; dimIndex++)
248 splineData.splineValues[dimIndex] =
249 SplineParamsDimVector(c_sampleSplineValuesFull).subArray(dimIndex * dimSize, dimSize);
250 splineData.splineDerivatives[dimIndex] =
251 SplineParamsDimVector(c_sampleSplineDerivativesFull).subArray(dimIndex * dimSize, dimSize);
253 atomData.splineDataByPmeOrder[pmeOrder] = splineData;
255 s_inputAtomDataSets_[atomCount] = atomData;
257 s_pmeTestHardwareContexts = createPmeTestHardwareContextList();
258 g_allowPmeWithSyclForTesting = true; // We support PmeGather with SYCL
261 static void TearDownTestSuite()
263 // Revert the value back.
264 g_allowPmeWithSyclForTesting = false;
268 static void runTest()
270 /* Getting the input */
275 SparseRealGridValuesInput nonZeroGridValues;
276 std::tie(box, pmeOrder, gridSize, nonZeroGridValues, atomCount) = GetParam();
277 auto inputAtomData = s_inputAtomDataSets_[atomCount];
278 auto inputAtomSplineData = inputAtomData.splineDataByPmeOrder[pmeOrder];
280 /* Storing the input where it's needed, running the test */
282 inputRec.nkx = gridSize[XX];
283 inputRec.nky = gridSize[YY];
284 inputRec.nkz = gridSize[ZZ];
285 inputRec.pme_order = pmeOrder;
286 inputRec.coulombtype = CoulombInteractionType::Pme;
287 inputRec.epsilon_r = 1.0;
289 TestReferenceData refData;
290 for (const auto& pmeTestHardwareContext : s_pmeTestHardwareContexts)
292 pmeTestHardwareContext->activate();
293 CodePath codePath = pmeTestHardwareContext->codePath();
294 const bool supportedInput = pmeSupportsInputForMode(
295 *getTestHardwareEnvironment()->hwinfo(), &inputRec, codePath);
298 /* Testing the failure for the unsupported input */
299 EXPECT_THROW_GMX(pmeInitWrapper(&inputRec, codePath, nullptr, nullptr, nullptr, box),
300 NotImplementedError);
304 /* Describing the test uniquely */
306 formatString("Testing force gathering on %s for PME grid size %d %d %d"
307 ", order %d, %zu atoms",
308 pmeTestHardwareContext->description().c_str(),
315 PmeSafePointer pmeSafe = pmeInitWrapper(&inputRec,
317 pmeTestHardwareContext->deviceContext(),
318 pmeTestHardwareContext->deviceStream(),
319 pmeTestHardwareContext->pmeGpuProgram(),
321 std::unique_ptr<StatePropagatorDataGpu> stateGpu =
322 (codePath == CodePath::GPU)
323 ? makeStatePropagatorDataGpu(*pmeSafe.get(),
324 pmeTestHardwareContext->deviceContext(),
325 pmeTestHardwareContext->deviceStream())
328 pmeInitAtoms(pmeSafe.get(),
331 inputAtomData.coordinates,
332 inputAtomData.charges);
334 /* Setting some more inputs */
335 pmeSetRealGrid(pmeSafe.get(), codePath, nonZeroGridValues);
336 pmeSetGridLineIndices(pmeSafe.get(), codePath, inputAtomData.gridLineIndices);
337 for (int dimIndex = 0; dimIndex < DIM; dimIndex++)
339 pmeSetSplineData(pmeSafe.get(),
341 inputAtomSplineData.splineValues[dimIndex],
342 PmeSplineDataType::Values,
344 pmeSetSplineData(pmeSafe.get(),
346 inputAtomSplineData.splineDerivatives[dimIndex],
347 PmeSplineDataType::Derivatives,
351 /* Explicitly copying the sample forces to be able to modify them */
352 auto inputForcesFull(c_sampleForcesFull);
353 GMX_RELEASE_ASSERT(inputForcesFull.size() >= atomCount, "Bad input forces size");
354 auto forces = ForcesVector(inputForcesFull).subArray(0, atomCount);
356 /* Running the force gathering itself */
357 pmePerformGather(pmeSafe.get(), codePath, forces);
358 pmeFinalizeTest(pmeSafe.get(), codePath);
360 /* Check the output forces correctness */
361 TestReferenceChecker forceChecker(refData.rootChecker());
362 const auto ulpTolerance = 3 * pmeOrder;
363 forceChecker.setDefaultTolerance(relativeToleranceAsUlp(1.0, ulpTolerance));
364 forceChecker.checkSequence(forces.begin(), forces.end(), "Forces");
368 static std::vector<std::unique_ptr<PmeTestHardwareContext>> s_pmeTestHardwareContexts;
371 std::vector<std::unique_ptr<PmeTestHardwareContext>> PmeGatherTest::s_pmeTestHardwareContexts;
373 // An instance of static atom data
374 InputDataByAtomCount PmeGatherTest::s_inputAtomDataSets_;
376 //! Test for PME force gathering
377 TEST_P(PmeGatherTest, ReproducesOutputs)
379 EXPECT_NO_THROW_GMX(runTest());
382 //! Instantiation of the PME gathering test
383 INSTANTIATE_TEST_SUITE_P(SaneInput,
385 ::testing::Combine(::testing::ValuesIn(c_sampleBoxes),
386 ::testing::ValuesIn(pmeOrders),
387 ::testing::ValuesIn(c_sampleGridSizes),
388 ::testing::ValuesIn(c_sampleGrids),
389 ::testing::ValuesIn(atomCounts)));