c7613dd1843b18af27af427a7ce59052793a8677
[alexxy/gromacs.git] / src / gromacs / gpu_utils / tests / hostallocator.cpp
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 2017,2018,2019, 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.
8  *
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.
13  *
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.
18  *
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.
23  *
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.
31  *
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.
34  */
35 /*! \internal \file
36  * \brief
37  * Tests for GPU host allocator.
38  *
39  * \author Mark Abraham <mark.j.abraham@gmail.com>
40  */
41 #include "gmxpre.h"
42
43 #include "gromacs/gpu_utils/hostallocator.h"
44
45 #include "config.h"
46
47 #include <type_traits>
48 #include <vector>
49
50 #include <gtest/gtest.h>
51
52 #include "gromacs/gpu_utils/gpu_utils.h"
53 #include "gromacs/math/vectypes.h"
54 #include "gromacs/utility/arrayref.h"
55 #include "gromacs/utility/real.h"
56
57 #include "gromacs/math/tests/testarrayrefs.h"
58
59 #include "devicetransfers.h"
60 #include "gputest.h"
61
62 namespace gmx
63 {
64
65 namespace test
66 {
67
68 /*! \internal \brief Typed test fixture for infrastructure for
69  * host-side memory used for GPU transfers. */
70 template <typename T>
71 class HostMemoryTest : public test::GpuTest
72 {
73     public:
74         //! Convenience type
75         using ValueType = T;
76 };
77
78 /*! \brief Convenience function to transform a view into one with base
79  * type of (non-const) char.
80  *
81  * This transformation is useful for using containers with C APIs
82  * where the function signature is not declared const even where the
83  * semantics of the usage actually are const.
84  *
85  * \param[in]    data   The data pointer.
86  * \param[in]    size   The size of the data pointer (in T).
87  * \tparam       T      The base type of the container
88  * */
89 template <typename T>
90 ArrayRef<char> charArrayRefFromArray(T *data, size_t size)
91 {
92     // Make a type like T, but without its possible const qualifier.
93     using NonConstT = typename std::remove_const<T>::type;
94     return arrayRefFromArray<char>(reinterpret_cast<char *>(const_cast<NonConstT *>(data)), size * sizeof(T));
95 }
96
97 //! Does a device transfer of \c input to the device in \c gpuInfo, and back to \c output.
98 template <typename T>
99 void runTest(const gmx_gpu_info_t &gpuInfo,
100              ArrayRef<T>           input,
101              ArrayRef<T>           output)
102 {
103     // Convert the views of input and output to flat non-const chars,
104     // so that there's no templating when we call doDeviceTransfers.
105     auto inputRef  = charArrayRefFromArray(input.data(), input.size());
106     auto outputRef = charArrayRefFromArray(output.data(), output.size());
107
108     ASSERT_EQ(inputRef.size(), outputRef.size());
109     doDeviceTransfers(gpuInfo, inputRef, outputRef);
110     compareViews(input, output);
111 }
112
113 struct MoveOnly {
114     MoveOnly(real x = 0) : x(x) {}
115     MoveOnly(const MoveOnly &)            = delete;
116     MoveOnly(MoveOnly &&)                 = default;
117     MoveOnly &operator=(const MoveOnly &) = delete;
118     MoveOnly &operator=(MoveOnly &&)      = default;
119     bool operator==(const MoveOnly &o) const { return x == o.x; }
120     real operator*=(int scaleFactor) { return x *= scaleFactor; }
121     real x;
122 };
123
124 }   // namespace test
125
126 namespace detail
127 {
128
129 template<>
130 struct PaddingTraits<test::MoveOnly>
131 {
132     using SimdBaseType = real;
133     static constexpr int maxSimdWidthOfBaseType = GMX_REAL_MAX_SIMD_WIDTH;
134 };
135
136 }   // namespace detail
137
138 namespace test
139 {
140
141 //! The types used in testing of all operations.
142 typedef ::testing::Types<int32_t, real, RVec, test::MoveOnly> TestTypes;
143
144 //! Typed test fixture
145 template <typename T>
146 struct HostAllocatorTest : HostMemoryTest<T>
147 {
148     using VectorType = PaddedHostVector<T>; //!< PaddedHostVector of type tested
149 };
150 TYPED_TEST_CASE(HostAllocatorTest, TestTypes);
151
152 //! Typed test fixture (no mem/gpu initializtion - much faster)
153 template <typename T>
154 struct HostAllocatorTestNoMem : ::testing::Test
155 {
156     using VectorType = PaddedHostVector<T>; //!< PaddedHostVector of type tested
157 };
158 TYPED_TEST_CASE(HostAllocatorTestNoMem, TestTypes);
159
160 //! Typed test fixture for tests requiring a copyable type
161 template <typename T>
162 struct HostAllocatorTestNoMemCopyable : HostAllocatorTestNoMem<T> {};
163 //! The types used in testing minus move only types
164 using TestTypesCopyable = ::testing::Types<int32_t, real, RVec>;
165
166 TYPED_TEST_CASE(HostAllocatorTestNoMemCopyable, TestTypesCopyable);
167
168 //! Typed test fixture for tests requiring a copyable type
169 template <typename T>
170 using HostAllocatorTestCopyable = HostAllocatorTest<T>;
171 TYPED_TEST_CASE(HostAllocatorTestCopyable, TestTypesCopyable);
172
173 // Note that in GoogleTest typed tests, the use of TestFixture:: and
174 // this-> is sometimes required to get access to things in the fixture
175 // class (or its base classes).
176
177 // Note also that aspects of this code can be tested even when a GPU
178 // device is not available.
179
180 TYPED_TEST(HostAllocatorTest, EmptyMemoryAlwaysWorks)
181 {
182     typename TestFixture::VectorType v;
183 }
184
185 TYPED_TEST(HostAllocatorTestCopyable, VectorsWithDefaultHostAllocatorAlwaysWorks)
186 {
187     typename TestFixture::VectorType input(3), output;
188     output.resizeWithPadding(input.size());
189 }
190
191 // Several tests actually do CUDA transfers. This is not necessary
192 // because the state of page alignment or pinning is not currently
193 // relevant to the success of a CUDA transfer. CUDA checks happen only
194 // during cudaHostRegister and cudaHostUnregister. Such tests are of
195 // value only when this behaviour changes, if ever.
196
197 TYPED_TEST(HostAllocatorTestCopyable, TransfersWithoutPinningWork)
198 {
199     typename TestFixture::VectorType input;
200     fillInput(&input, 1);
201     typename TestFixture::VectorType output;
202     output.resizeWithPadding(input.size());
203
204     runTest(*this->gpuInfo_, makeArrayRef(input), makeArrayRef(output));
205 }
206
207 TYPED_TEST(HostAllocatorTestCopyable, FillInputAlsoWorksAfterCallingReserve)
208 {
209     typename TestFixture::VectorType input;
210     input.reserveWithPadding(3);
211     fillInput(&input, 1);
212 }
213
214 TYPED_TEST(HostAllocatorTestNoMem, CreateVector)
215 {
216     typename TestFixture::VectorType input1;
217     EXPECT_FALSE(input1.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
218     typename TestFixture::VectorType input2({PinningPolicy::PinnedIfSupported});
219     EXPECT_TRUE (input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
220 }
221
222 TYPED_TEST(HostAllocatorTestNoMem, MoveAssignment)
223 {
224     typename TestFixture::VectorType input1({PinningPolicy::PinnedIfSupported});
225     input1 = typename TestFixture::VectorType();
226     EXPECT_FALSE(input1.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
227
228     typename TestFixture::VectorType input2;
229     input2 = typename TestFixture::VectorType({PinningPolicy::PinnedIfSupported});
230     EXPECT_TRUE (input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
231 }
232
233 TYPED_TEST(HostAllocatorTestNoMem, MoveConstruction)
234 {
235     typename TestFixture::VectorType input1;
236     typename TestFixture::VectorType input2(std::move(input1));
237     EXPECT_FALSE(input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
238
239     typename TestFixture::VectorType input3({PinningPolicy::PinnedIfSupported});
240     typename TestFixture::VectorType input4(std::move(input3));
241     EXPECT_TRUE(input4.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
242 }
243
244 TYPED_TEST(HostAllocatorTestNoMemCopyable, CopyAssignment)
245 {
246     typename TestFixture::VectorType input1;
247     typename TestFixture::VectorType input2({PinningPolicy::PinnedIfSupported});
248     input1 = input2;
249     EXPECT_FALSE(input1.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
250     EXPECT_TRUE (input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
251     input2 = input1;
252     EXPECT_FALSE(input1.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
253     EXPECT_TRUE (input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
254 }
255
256 TYPED_TEST(HostAllocatorTestNoMemCopyable, CopyConstruction)
257 {
258     typename TestFixture::VectorType input1;
259     typename TestFixture::VectorType input2(input1); //NOLINT(performance-unnecessary-copy-initialization)
260     EXPECT_FALSE(input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
261
262     typename TestFixture::VectorType input3({PinningPolicy::PinnedIfSupported});
263     typename TestFixture::VectorType input4(input3); //NOLINT(performance-unnecessary-copy-initialization)
264     EXPECT_FALSE(input4.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
265 }
266
267 TYPED_TEST(HostAllocatorTestNoMem, Swap)
268 {
269     typename TestFixture::VectorType input1;
270     typename TestFixture::VectorType input2({PinningPolicy::PinnedIfSupported});
271     std::swap(input1, input2);
272     EXPECT_TRUE (input1.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
273     EXPECT_FALSE(input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
274     std::swap(input2, input1);
275     EXPECT_FALSE(input1.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
276     EXPECT_TRUE (input2.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
277 }
278
279 TYPED_TEST(HostAllocatorTestNoMem, Comparison)
280 {
281     using AllocatorType = typename TestFixture::VectorType::allocator_type;
282     EXPECT_EQ(AllocatorType {}, AllocatorType {});
283     //Should be false for different pinning policy
284     EXPECT_NE(AllocatorType {}, AllocatorType {PinningPolicy::PinnedIfSupported});
285 }
286
287 #if GMX_GPU == GMX_GPU_CUDA
288
289 // Policy suitable for pinning is only supported for a CUDA build
290
291 TYPED_TEST(HostAllocatorTestCopyable, TransfersWithPinningWorkWithCuda)
292 {
293     if (!this->haveValidGpus())
294     {
295         return;
296     }
297
298     typename TestFixture::VectorType input;
299     changePinningPolicy(&input, PinningPolicy::PinnedIfSupported);
300     fillInput(&input, 1);
301     typename TestFixture::VectorType output;
302     changePinningPolicy(&output, PinningPolicy::PinnedIfSupported);
303     output.resizeWithPadding(input.size());
304
305     runTest(*this->gpuInfo_, makeArrayRef(input), makeArrayRef(output));
306 }
307
308 //! Helper function for wrapping a call to isHostMemoryPinned.
309 template <typename VectorType>
310 bool isPinned(const VectorType &v)
311 {
312     void *data = const_cast<void *>(static_cast<const void *>(v.data()));
313     return isHostMemoryPinned(data);
314 }
315
316 TYPED_TEST(HostAllocatorTestCopyable, ManualPinningOperationsWorkWithCuda)
317 {
318     if (!this->haveValidGpus())
319     {
320         return;
321     }
322
323     typename TestFixture::VectorType input;
324     changePinningPolicy(&input, PinningPolicy::PinnedIfSupported);
325     EXPECT_TRUE(input.get_allocator().pinningPolicy() == PinningPolicy::PinnedIfSupported);
326     EXPECT_EQ(0, input.size());
327     EXPECT_EQ(0, input.paddedSize());
328     EXPECT_TRUE(input.empty());
329     EXPECT_FALSE(isPinned(input));
330
331     // Fill some contents, which will be pinned because of the policy.
332     fillInput(&input, 1);
333     EXPECT_TRUE(isPinned(input));
334
335     // Switching policy to CannotBePinned must unpin the buffer (via
336     // realloc and copy).
337     auto oldInputData = input.data();
338     changePinningPolicy(&input, PinningPolicy::CannotBePinned);
339     EXPECT_FALSE(isPinned(input));
340     // These cannot be equal as both had to be allocated at the same
341     // time for the contents to be able to be copied.
342     EXPECT_NE(oldInputData, input.data());
343
344     // Switching policy to PinnedIfSupported must pin the buffer (via
345     // realloc and copy).
346     oldInputData = input.data();
347     changePinningPolicy(&input, PinningPolicy::PinnedIfSupported);
348     EXPECT_TRUE(isPinned(input));
349     // These cannot be equal as both had to be allocated at the same
350     // time for the contents to be able to be copied.
351     EXPECT_NE(oldInputData, input.data());
352 }
353
354 #endif
355
356 TYPED_TEST(HostAllocatorTest, StatefulAllocatorUsesMemory)
357 {
358     // The HostAllocator has state, so a container using it will be
359     // larger than a normal vector, whose default allocator is
360     // stateless.
361     EXPECT_LT(sizeof(std::vector<typename TestFixture::VectorType::value_type>),
362               sizeof(typename TestFixture::VectorType));
363 }
364
365 TEST(HostAllocatorUntypedTest, Comparison)
366 {
367     //Should always be true for the same policy, indpendent of value_type
368     EXPECT_EQ(HostAllocator<float>{}, HostAllocator<double>{});
369 }
370
371 //! Declare allocator types to test.
372 using AllocatorTypesToTest = ::testing::Types<HostAllocator<real>,
373                                               HostAllocator<int32_t>,
374                                               HostAllocator<RVec>,
375                                               HostAllocator<MoveOnly>
376                                               >;
377
378 TYPED_TEST_CASE(AllocatorTest, AllocatorTypesToTest);
379
380 } // namespace test
381 } // namespace gmx
382
383 // Includes tests common to all allocation policies.
384 #include "gromacs/utility/tests/alignedallocator_impl.h"