2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 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.
36 * \brief Tests for the halo exchange
38 * The test sets up the rank topology and performs a coordinate halo
39 * exchange (for both CPU and GPU codepaths) for several 1D and 2D
40 * pulse configirations. Each pulse involves a few non-contiguous
41 * indices. The sending rank, atom number and spatial 3D index are
42 * encoded in the x values, to allow correctness checking following
47 * \author Alan Gray <alang@nvidia.com>
48 * \ingroup module_domdec
58 #include <gtest/gtest.h>
60 #include "gromacs/domdec/atomdistribution.h"
61 #include "gromacs/domdec/domdec_internal.h"
62 #include "gromacs/domdec/gpuhaloexchange.h"
64 # include "gromacs/gpu_utils/device_stream.h"
65 # include "gromacs/gpu_utils/devicebuffer.h"
67 #include "gromacs/gpu_utils/gpueventsynchronizer.h"
68 #include "gromacs/gpu_utils/hostallocator.h"
69 #include "gromacs/mdtypes/inputrec.h"
71 #include "testutils/mpitest.h"
72 #include "testutils/test_hardware_environment.h"
81 /*! \brief Get encoded numerical value for sending rank, atom number and spatial 3D index
83 * \param [in] sendRank MPI rank of sender
84 * \param [in] atomNumber Atom number
85 * \param [in] spatial3dIndex Spatial 3D Index
87 * \returns Encoded value
89 float encodedValue(const int sendRank, const int atomNumber, const int spatial3dIndex)
91 return sendRank * 1000 + atomNumber * 100 + spatial3dIndex;
94 /*! \brief Initialize halo array
96 * \param [in] x Atom coordinate data array
97 * \param [in] numHomeAtoms Number of home atoms
98 * \param [in] numAtomsTotal Total number of atoms, including halo
100 void initHaloData(RVec* x, const int numHomeAtoms, const int numAtomsTotal)
103 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
105 for (int i = 0; i < numAtomsTotal; i++)
107 for (int j = 0; j < DIM; j++)
109 x[i][j] = i < numHomeAtoms ? encodedValue(rank, i, j) : -1;
114 /*! \brief Perform GPU halo exchange, including required setup and data transfers
116 * \param [in] dd Domain decomposition object
117 * \param [in] box Box matrix
118 * \param [in] h_x Atom coordinate data array on host
119 * \param [in] numAtomsTotal Total number of atoms, including halo
121 void gpuHalo(gmx_domdec_t* dd, matrix box, HostVector<RVec>* h_x, int numAtomsTotal)
123 #if (GMX_GPU_CUDA && GMX_THREAD_MPI)
124 // pin memory if possible
125 changePinningPolicy(h_x, PinningPolicy::PinnedIfSupported);
126 // Set up GPU hardware environment and assign this MPI rank to a device
128 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
129 int numDevices = getTestHardwareEnvironment()->getTestDeviceList().size();
130 const auto& testDevice = getTestHardwareEnvironment()->getTestDeviceList()[rank % numDevices];
131 const auto& deviceContext = testDevice->deviceContext();
132 setActiveDevice(testDevice->deviceInfo());
133 DeviceStream deviceStream(deviceContext, DeviceStreamPriority::Normal, false);
135 // Set up GPU buffer and copy input data from host
136 DeviceBuffer<RVec> d_x;
138 int d_x_size_alloc = -1;
139 reallocateDeviceBuffer(&d_x, numAtomsTotal, &d_x_size, &d_x_size_alloc, deviceContext);
141 copyToDeviceBuffer(&d_x, h_x->data(), 0, numAtomsTotal, deviceStream, GpuApiCallBehavior::Sync, nullptr);
143 GpuEventSynchronizer coordinatesReadyOnDeviceEvent;
144 coordinatesReadyOnDeviceEvent.markEvent(deviceStream);
146 std::array<std::vector<GpuHaloExchange>, DIM> gpuHaloExchange;
148 // Create halo exchange objects
149 for (int d = 0; d < dd->ndim; d++)
151 for (int pulse = 0; pulse < dd->comm->cd[d].numPulses(); pulse++)
153 gpuHaloExchange[d].push_back(
154 GpuHaloExchange(dd, d, MPI_COMM_WORLD, deviceContext, pulse, nullptr));
158 // Perform GPU halo exchange
159 for (int d = 0; d < dd->ndim; d++)
161 for (int pulse = 0; pulse < dd->comm->cd[d].numPulses(); pulse++)
163 gpuHaloExchange[d][pulse].reinitHalo(d_x, nullptr);
164 gpuHaloExchange[d][pulse].communicateHaloCoordinates(box, &coordinatesReadyOnDeviceEvent);
167 MPI_Barrier(MPI_COMM_WORLD);
169 GpuEventSynchronizer haloCompletedEvent;
170 haloCompletedEvent.markEvent(deviceStream);
171 haloCompletedEvent.waitForEvent();
173 // Copy results back to host
174 copyFromDeviceBuffer(
175 h_x->data(), &d_x, 0, numAtomsTotal, deviceStream, GpuApiCallBehavior::Sync, nullptr);
177 freeDeviceBuffer(d_x);
179 GMX_UNUSED_VALUE(dd);
180 GMX_UNUSED_VALUE(box);
181 GMX_UNUSED_VALUE(h_x);
182 GMX_UNUSED_VALUE(numAtomsTotal);
186 /*! \brief Define 1D rank topology with 4 MPI tasks
188 * \param [in] dd Domain decomposition object
190 void define1dRankTopology(gmx_domdec_t* dd)
193 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
195 const int numRanks = getNumberOfTestMpiRanks();
196 dd->neighbor[0][0] = (rank + 1) % numRanks;
197 dd->neighbor[0][1] = (rank == 0) ? (numRanks - 1) : rank - 1;
200 /*! \brief Define 2D rank topology with 4 MPI tasks
207 * \param [in] dd Domain decomposition object
209 void define2dRankTopology(gmx_domdec_t* dd)
213 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
218 dd->neighbor[0][0] = 1;
219 dd->neighbor[0][1] = 1;
220 dd->neighbor[1][0] = 2;
221 dd->neighbor[1][1] = 2;
224 dd->neighbor[0][0] = 0;
225 dd->neighbor[0][1] = 0;
226 dd->neighbor[1][0] = 3;
227 dd->neighbor[1][1] = 3;
230 dd->neighbor[0][0] = 3;
231 dd->neighbor[0][1] = 3;
232 dd->neighbor[1][0] = 0;
233 dd->neighbor[1][1] = 0;
236 dd->neighbor[0][0] = 2;
237 dd->neighbor[0][1] = 2;
238 dd->neighbor[1][0] = 1;
239 dd->neighbor[1][1] = 1;
244 /*! \brief Define a 1D halo with 1 pulses
246 * \param [in] dd Domain decomposition object
247 * \param [in] indvec Vector of index vectors
249 void define1dHaloWith1Pulse(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
253 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
255 std::vector<int> indexvec;
256 gmx_domdec_ind_t ind;
262 // Set up indices involved in halo
266 dd->comm->cd[dimIndex].receiveInPlace = true;
267 dd->dim[dimIndex] = 0;
268 dd->ci[dimIndex] = rank;
270 // First pulse involves (arbitrary) indices 1 and 3
271 indexvec.push_back(1);
272 indexvec.push_back(3);
274 ind.index = indexvec;
275 ind.nsend[nzone + 1] = 2;
276 ind.nrecv[nzone + 1] = 2;
277 indvec->push_back(ind);
279 dd->comm->cd[dimIndex].ind = *indvec;
282 /*! \brief Define a 1D halo with 2 pulses
284 * \param [in] dd Domain decomposition object
285 * \param [in] indvec Vector of index vectors
287 void define1dHaloWith2Pulses(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
291 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
293 std::vector<int> indexvec;
294 gmx_domdec_ind_t ind;
300 // Set up indices involved in halo
304 dd->comm->cd[dimIndex].receiveInPlace = true;
305 dd->dim[dimIndex] = 0;
306 dd->ci[dimIndex] = rank;
308 // First pulse involves (arbitrary) indices 1 and 3
309 indexvec.push_back(1);
310 indexvec.push_back(3);
312 ind.index = indexvec;
313 ind.nsend[nzone + 1] = 2;
314 ind.nrecv[nzone + 1] = 2;
315 indvec->push_back(ind);
317 // Add another pulse with (arbitrary) indices 4,5,7
320 indexvec.push_back(4);
321 indexvec.push_back(5);
322 indexvec.push_back(7);
324 ind.index = indexvec;
325 ind.nsend[nzone + 1] = 3;
326 ind.nrecv[nzone + 1] = 3;
327 indvec->push_back(ind);
329 dd->comm->cd[dimIndex].ind = *indvec;
332 /*! \brief Define a 2D halo with 1 pulse in each dimension
334 * \param [in] dd Domain decomposition object
335 * \param [in] indvec Vector of index vectors
337 void define2dHaloWith1PulseInEachDim(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
341 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
343 std::vector<int> indexvec;
344 gmx_domdec_ind_t ind;
348 for (int dimIndex = 0; dimIndex < dd->ndim; dimIndex++)
351 // Set up indices involved in halo
355 dd->comm->cd[dimIndex].receiveInPlace = true;
356 dd->dim[dimIndex] = 0;
357 dd->ci[dimIndex] = rank;
359 // Single pulse involving (arbitrary) indices 1 and 3
360 indexvec.push_back(1);
361 indexvec.push_back(3);
363 ind.index = indexvec;
364 ind.nsend[nzone + 1] = 2;
365 ind.nrecv[nzone + 1] = 2;
366 indvec->push_back(ind);
368 dd->comm->cd[dimIndex].ind = *indvec;
374 /*! \brief Define a 2D halo with 2 pulses in the first dimension
376 * \param [in] dd Domain decomposition object
377 * \param [in] indvec Vector of index vectors
379 void define2dHaloWith2PulsesInDim1(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
383 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
385 std::vector<int> indexvec;
386 gmx_domdec_ind_t ind;
390 for (int dimIndex = 0; dimIndex < dd->ndim; dimIndex++)
393 // Set up indices involved in halo
397 dd->comm->cd[dimIndex].receiveInPlace = true;
398 dd->dim[dimIndex] = 0;
399 dd->ci[dimIndex] = rank;
401 // First pulse involves (arbitrary) indices 1 and 3
402 indexvec.push_back(1);
403 indexvec.push_back(3);
405 ind.index = indexvec;
406 ind.nsend[nzone + 1] = 2;
407 ind.nrecv[nzone + 1] = 2;
408 indvec->push_back(ind);
410 if (dimIndex == 0) // Add another pulse with (arbitrary) indices 4,5,7
414 indexvec.push_back(4);
415 indexvec.push_back(5);
416 indexvec.push_back(7);
418 ind.index = indexvec;
419 ind.nsend[nzone + 1] = 3;
420 ind.nrecv[nzone + 1] = 3;
421 indvec->push_back(ind);
424 dd->comm->cd[dimIndex].ind = *indvec;
430 /*! \brief Check results for above-defined 1D halo with 1 pulse
432 * \param [in] x Atom coordinate data array
433 * \param [in] dd Domain decomposition object
434 * \param [in] numHomeAtoms Number of home atoms
436 void checkResults1dHaloWith1Pulse(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
438 // Check results are expected from values encoded in x data
439 for (int j = 0; j < DIM; j++)
441 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
442 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
443 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
447 /*! \brief Check results for above-defined 1D halo with 2 pulses
449 * \param [in] x Atom coordinate data array
450 * \param [in] dd Domain decomposition object
451 * \param [in] numHomeAtoms Number of home atoms
453 void checkResults1dHaloWith2Pulses(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
455 // Check results are expected from values encoded in x data
456 for (int j = 0; j < DIM; j++)
458 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
459 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
460 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
461 // Second Pulse in first dim: atoms 4,5,7 from forward horizontal neighbour
462 EXPECT_EQ(x[numHomeAtoms + 2][j], encodedValue(dd->neighbor[0][0], 4, j));
463 EXPECT_EQ(x[numHomeAtoms + 3][j], encodedValue(dd->neighbor[0][0], 5, j));
464 EXPECT_EQ(x[numHomeAtoms + 4][j], encodedValue(dd->neighbor[0][0], 7, j));
468 /*! \brief Check results for above-defined 2D halo with 1 pulse in each dimension
470 * \param [in] x Atom coordinate data array
471 * \param [in] dd Domain decomposition object
472 * \param [in] numHomeAtoms Number of home atoms
474 void checkResults2dHaloWith1PulseInEachDim(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
476 // Check results are expected from values encoded in x data
477 for (int j = 0; j < DIM; j++)
479 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
480 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
481 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
482 // First Pulse in second dim: atoms 1 and 3 from forward vertical neighbour
483 EXPECT_EQ(x[numHomeAtoms + 2][j], encodedValue(dd->neighbor[1][0], 1, j));
484 EXPECT_EQ(x[numHomeAtoms + 3][j], encodedValue(dd->neighbor[1][0], 3, j));
488 /*! \brief Check results for above-defined 2D halo with 2 pulses in the first dimension
490 * \param [in] x Atom coordinate data array
491 * \param [in] dd Domain decomposition object
492 * \param [in] numHomeAtoms Number of home atoms
494 void checkResults2dHaloWith2PulsesInDim1(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
496 // Check results are expected from values encoded in x data
497 for (int j = 0; j < DIM; j++)
499 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
500 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
501 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
502 // Second Pulse in first dim: atoms 4,5,7 from forward horizontal neighbour
503 EXPECT_EQ(x[numHomeAtoms + 2][j], encodedValue(dd->neighbor[0][0], 4, j));
504 EXPECT_EQ(x[numHomeAtoms + 3][j], encodedValue(dd->neighbor[0][0], 5, j));
505 EXPECT_EQ(x[numHomeAtoms + 4][j], encodedValue(dd->neighbor[0][0], 7, j));
506 // First Pulse in second dim: atoms 1 and 3 from forward vertical neighbour
507 EXPECT_EQ(x[numHomeAtoms + 5][j], encodedValue(dd->neighbor[1][0], 1, j));
508 EXPECT_EQ(x[numHomeAtoms + 6][j], encodedValue(dd->neighbor[1][0], 3, j));
512 TEST(HaloExchangeTest, Coordinates1dHaloWith1Pulse)
514 GMX_MPI_TEST(RequireRankCount<4>);
517 const int numHomeAtoms = 10;
518 const int numHaloAtoms = 2;
519 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
520 HostVector<RVec> h_x;
521 h_x.resize(numAtomsTotal);
523 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
528 dd.mpi_comm_all = MPI_COMM_WORLD;
529 gmx_domdec_comm_t comm;
531 dd.unitCellInfo.haveScrewPBC = false;
533 DDAtomRanges atomRanges;
534 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
535 dd.comm->atomRanges = atomRanges;
537 define1dRankTopology(&dd);
539 std::vector<gmx_domdec_ind_t> indvec;
540 define1dHaloWith1Pulse(&dd, &indvec);
542 // Perform halo exchange
543 matrix box = { { 0., 0., 0. } };
544 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
547 checkResults1dHaloWith1Pulse(h_x.data(), &dd, numHomeAtoms);
549 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
551 // early return if no devices are available.
552 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
557 // Re-initialize input
558 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
560 // Perform GPU halo exchange
561 gpuHalo(&dd, box, &h_x, numAtomsTotal);
564 checkResults1dHaloWith1Pulse(h_x.data(), &dd, numHomeAtoms);
568 TEST(HaloExchangeTest, Coordinates1dHaloWith2Pulses)
570 GMX_MPI_TEST(RequireRankCount<4>);
573 const int numHomeAtoms = 10;
574 const int numHaloAtoms = 5;
575 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
576 HostVector<RVec> h_x;
577 h_x.resize(numAtomsTotal);
579 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
584 dd.mpi_comm_all = MPI_COMM_WORLD;
585 gmx_domdec_comm_t comm;
587 dd.unitCellInfo.haveScrewPBC = false;
589 DDAtomRanges atomRanges;
590 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
591 dd.comm->atomRanges = atomRanges;
593 define1dRankTopology(&dd);
595 std::vector<gmx_domdec_ind_t> indvec;
596 define1dHaloWith2Pulses(&dd, &indvec);
598 // Perform halo exchange
599 matrix box = { { 0., 0., 0. } };
600 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
603 checkResults1dHaloWith2Pulses(h_x.data(), &dd, numHomeAtoms);
605 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
607 // early return if no devices are available.
608 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
613 // Re-initialize input
614 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
616 // Perform GPU halo exchange
617 gpuHalo(&dd, box, &h_x, numAtomsTotal);
620 checkResults1dHaloWith2Pulses(h_x.data(), &dd, numHomeAtoms);
625 TEST(HaloExchangeTest, Coordinates2dHaloWith1PulseInEachDim)
627 GMX_MPI_TEST(RequireRankCount<4>);
630 const int numHomeAtoms = 10;
631 const int numHaloAtoms = 4;
632 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
633 HostVector<RVec> h_x;
634 h_x.resize(numAtomsTotal);
636 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
641 dd.mpi_comm_all = MPI_COMM_WORLD;
642 gmx_domdec_comm_t comm;
644 dd.unitCellInfo.haveScrewPBC = false;
646 DDAtomRanges atomRanges;
647 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
648 dd.comm->atomRanges = atomRanges;
650 define2dRankTopology(&dd);
652 std::vector<gmx_domdec_ind_t> indvec;
653 define2dHaloWith1PulseInEachDim(&dd, &indvec);
655 // Perform halo exchange
656 matrix box = { { 0., 0., 0. } };
657 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
660 checkResults2dHaloWith1PulseInEachDim(h_x.data(), &dd, numHomeAtoms);
662 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
664 // early return if no devices are available.
665 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
670 // Re-initialize input
671 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
673 // Perform GPU halo exchange
674 gpuHalo(&dd, box, &h_x, numAtomsTotal);
677 checkResults2dHaloWith1PulseInEachDim(h_x.data(), &dd, numHomeAtoms);
681 TEST(HaloExchangeTest, Coordinates2dHaloWith2PulsesInDim1)
683 GMX_MPI_TEST(RequireRankCount<4>);
686 const int numHomeAtoms = 10;
687 const int numHaloAtoms = 7;
688 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
689 HostVector<RVec> h_x;
690 h_x.resize(numAtomsTotal);
692 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
697 dd.mpi_comm_all = MPI_COMM_WORLD;
698 gmx_domdec_comm_t comm;
700 dd.unitCellInfo.haveScrewPBC = false;
702 DDAtomRanges atomRanges;
703 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
704 dd.comm->atomRanges = atomRanges;
706 define2dRankTopology(&dd);
708 std::vector<gmx_domdec_ind_t> indvec;
709 define2dHaloWith2PulsesInDim1(&dd, &indvec);
711 // Perform halo exchange
712 matrix box = { { 0., 0., 0. } };
713 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
716 checkResults2dHaloWith2PulsesInDim1(h_x.data(), &dd, numHomeAtoms);
718 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
720 // early return if no devices are available.
721 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
726 // Re-initialize input
727 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
729 // Perform GPU halo exchange
730 gpuHalo(&dd, box, &h_x, numAtomsTotal);
733 checkResults2dHaloWith2PulsesInDim1(h_x.data(), &dd, numHomeAtoms);