2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2020, 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
57 #include <gtest/gtest.h>
59 #include "gromacs/domdec/atomdistribution.h"
60 #include "gromacs/domdec/domdec_internal.h"
61 #include "gromacs/domdec/gpuhaloexchange.h"
63 # include "gromacs/gpu_utils/device_stream.h"
64 # include "gromacs/gpu_utils/devicebuffer.h"
65 # include "gromacs/gpu_utils/gpueventsynchronizer.cuh"
67 #include "gromacs/gpu_utils/hostallocator.h"
68 #include "gromacs/mdtypes/inputrec.h"
70 #include "testutils/mpitest.h"
71 #include "testutils/test_hardware_environment.h"
80 /*! \brief Get encoded numerical value for sending rank, atom number and spatial 3D index
82 * \param [in] sendRank MPI rank of sender
83 * \param [in] atomNumber Atom number
84 * \param [in] spatial3dIndex Spatial 3D Index
86 * \returns Encoded value
88 float encodedValue(const int sendRank
, const int atomNumber
, const int spatial3dIndex
)
90 return sendRank
* 1000 + atomNumber
* 100 + spatial3dIndex
;
93 /*! \brief Initialize halo array
95 * \param [in] x Atom coordinate data array
96 * \param [in] numHomeAtoms Number of home atoms
97 * \param [in] numAtomsTotal Total number of atoms, including halo
99 void initHaloData(RVec
* x
, const int numHomeAtoms
, const int numAtomsTotal
)
102 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
104 for (int i
= 0; i
< numAtomsTotal
; i
++)
106 for (int j
= 0; j
< DIM
; j
++)
108 x
[i
][j
] = i
< numHomeAtoms
? encodedValue(rank
, i
, j
) : -1;
113 /*! \brief Perform GPU halo exchange, including required setup and data transfers
115 * \param [in] dd Domain decomposition object
116 * \param [in] box Box matrix
117 * \param [in] h_x Atom coordinate data array on host
118 * \param [in] numAtomsTotal Total number of atoms, including halo
120 void gpuHalo(gmx_domdec_t
* dd
, matrix box
, HostVector
<RVec
>* h_x
, int numAtomsTotal
)
122 #if (GMX_GPU_CUDA && GMX_THREAD_MPI)
123 // pin memory if possible
124 changePinningPolicy(h_x
, PinningPolicy::PinnedIfSupported
);
125 // Set up GPU hardware environment and assign this MPI rank to a device
127 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
128 int numDevices
= getTestHardwareEnvironment()->getTestDeviceList().size();
129 const auto& testDevice
= getTestHardwareEnvironment()->getTestDeviceList()[rank
% numDevices
];
130 const auto& deviceContext
= testDevice
->deviceContext();
131 setActiveDevice(testDevice
->deviceInfo());
132 DeviceStream
deviceStream(deviceContext
, DeviceStreamPriority::Normal
, false);
134 // Set up GPU buffer and copy input data from host
135 DeviceBuffer
<RVec
> d_x
;
137 int d_x_size_alloc
= -1;
138 reallocateDeviceBuffer(&d_x
, numAtomsTotal
, &d_x_size
, &d_x_size_alloc
, deviceContext
);
140 copyToDeviceBuffer(&d_x
, h_x
->data(), 0, numAtomsTotal
, deviceStream
, GpuApiCallBehavior::Sync
, nullptr);
142 GpuEventSynchronizer coordinatesReadyOnDeviceEvent
;
143 coordinatesReadyOnDeviceEvent
.markEvent(deviceStream
);
145 // Perform GPU halo exchange
146 for (int d
= 0; d
< dd
->ndim
; d
++)
148 for (int pulse
= 0; pulse
< dd
->comm
->cd
[d
].numPulses(); pulse
++)
150 GpuHaloExchange
gpuHaloExchange(dd
, d
, MPI_COMM_WORLD
, deviceContext
, deviceStream
,
151 deviceStream
, pulse
, nullptr);
152 gpuHaloExchange
.reinitHalo(d_x
, nullptr);
153 gpuHaloExchange
.communicateHaloCoordinates(box
, &coordinatesReadyOnDeviceEvent
);
157 GpuEventSynchronizer haloCompletedEvent
;
158 haloCompletedEvent
.markEvent(deviceStream
);
159 haloCompletedEvent
.waitForEvent();
161 // Copy results back to host
162 copyFromDeviceBuffer(h_x
->data(), &d_x
, 0, numAtomsTotal
, deviceStream
,
163 GpuApiCallBehavior::Sync
, nullptr);
165 freeDeviceBuffer(d_x
);
167 GMX_UNUSED_VALUE(dd
);
168 GMX_UNUSED_VALUE(box
);
169 GMX_UNUSED_VALUE(h_x
);
170 GMX_UNUSED_VALUE(numAtomsTotal
);
174 /*! \brief Define 1D rank topology with 4 MPI tasks
176 * \param [in] dd Domain decomposition object
178 void define1dRankTopology(gmx_domdec_t
* dd
)
181 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
183 dd
->neighbor
[0][0] = (rank
+ 1) % 4;
184 dd
->neighbor
[0][1] = (rank
== 0) ? 3 : rank
- 1;
187 /*! \brief Define 2D rank topology with 4 MPI tasks
194 * \param [in] dd Domain decomposition object
196 void define2dRankTopology(gmx_domdec_t
* dd
)
200 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
205 dd
->neighbor
[0][0] = 1;
206 dd
->neighbor
[0][1] = 1;
207 dd
->neighbor
[1][0] = 2;
208 dd
->neighbor
[1][1] = 2;
211 dd
->neighbor
[0][0] = 0;
212 dd
->neighbor
[0][1] = 0;
213 dd
->neighbor
[1][0] = 3;
214 dd
->neighbor
[1][1] = 3;
217 dd
->neighbor
[0][0] = 3;
218 dd
->neighbor
[0][1] = 3;
219 dd
->neighbor
[1][0] = 0;
220 dd
->neighbor
[1][1] = 0;
223 dd
->neighbor
[0][0] = 2;
224 dd
->neighbor
[0][1] = 2;
225 dd
->neighbor
[1][0] = 1;
226 dd
->neighbor
[1][1] = 1;
231 /*! \brief Define a 1D halo with 1 pulses
233 * \param [in] dd Domain decomposition object
234 * \param [in] indvec Vector of index vectors
236 void define1dHaloWith1Pulse(gmx_domdec_t
* dd
, std::vector
<gmx_domdec_ind_t
>* indvec
)
240 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
242 std::vector
<int> indexvec
;
243 gmx_domdec_ind_t ind
;
249 // Set up indices involved in halo
253 dd
->comm
->cd
[dimIndex
].receiveInPlace
= true;
254 dd
->dim
[dimIndex
] = 0;
255 dd
->ci
[dimIndex
] = rank
;
257 // First pulse involves (arbitrary) indices 1 and 3
258 indexvec
.push_back(1);
259 indexvec
.push_back(3);
261 ind
.index
= indexvec
;
262 ind
.nsend
[nzone
+ 1] = 2;
263 ind
.nrecv
[nzone
+ 1] = 2;
264 indvec
->push_back(ind
);
266 dd
->comm
->cd
[dimIndex
].ind
= *indvec
;
269 /*! \brief Define a 1D halo with 2 pulses
271 * \param [in] dd Domain decomposition object
272 * \param [in] indvec Vector of index vectors
274 void define1dHaloWith2Pulses(gmx_domdec_t
* dd
, std::vector
<gmx_domdec_ind_t
>* indvec
)
278 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
280 std::vector
<int> indexvec
;
281 gmx_domdec_ind_t ind
;
287 // Set up indices involved in halo
291 dd
->comm
->cd
[dimIndex
].receiveInPlace
= true;
292 dd
->dim
[dimIndex
] = 0;
293 dd
->ci
[dimIndex
] = rank
;
295 // First pulse involves (arbitrary) indices 1 and 3
296 indexvec
.push_back(1);
297 indexvec
.push_back(3);
299 ind
.index
= indexvec
;
300 ind
.nsend
[nzone
+ 1] = 2;
301 ind
.nrecv
[nzone
+ 1] = 2;
302 indvec
->push_back(ind
);
304 // Add another pulse with (arbitrary) indices 4,5,7
307 indexvec
.push_back(4);
308 indexvec
.push_back(5);
309 indexvec
.push_back(7);
311 ind
.index
= indexvec
;
312 ind
.nsend
[nzone
+ 1] = 3;
313 ind
.nrecv
[nzone
+ 1] = 3;
314 indvec
->push_back(ind
);
316 dd
->comm
->cd
[dimIndex
].ind
= *indvec
;
319 /*! \brief Define a 2D halo with 1 pulse in each dimension
321 * \param [in] dd Domain decomposition object
322 * \param [in] indvec Vector of index vectors
324 void define2dHaloWith1PulseInEachDim(gmx_domdec_t
* dd
, std::vector
<gmx_domdec_ind_t
>* indvec
)
328 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
330 std::vector
<int> indexvec
;
331 gmx_domdec_ind_t ind
;
335 for (int dimIndex
= 0; dimIndex
< dd
->ndim
; dimIndex
++)
338 // Set up indices involved in halo
342 dd
->comm
->cd
[dimIndex
].receiveInPlace
= true;
343 dd
->dim
[dimIndex
] = 0;
344 dd
->ci
[dimIndex
] = rank
;
346 // Single pulse involving (arbitrary) indices 1 and 3
347 indexvec
.push_back(1);
348 indexvec
.push_back(3);
350 ind
.index
= indexvec
;
351 ind
.nsend
[nzone
+ 1] = 2;
352 ind
.nrecv
[nzone
+ 1] = 2;
353 indvec
->push_back(ind
);
355 dd
->comm
->cd
[dimIndex
].ind
= *indvec
;
361 /*! \brief Define a 2D halo with 2 pulses in the first dimension
363 * \param [in] dd Domain decomposition object
364 * \param [in] indvec Vector of index vectors
366 void define2dHaloWith2PulsesInDim1(gmx_domdec_t
* dd
, std::vector
<gmx_domdec_ind_t
>* indvec
)
370 MPI_Comm_rank(MPI_COMM_WORLD
, &rank
);
372 std::vector
<int> indexvec
;
373 gmx_domdec_ind_t ind
;
377 for (int dimIndex
= 0; dimIndex
< dd
->ndim
; dimIndex
++)
380 // Set up indices involved in halo
384 dd
->comm
->cd
[dimIndex
].receiveInPlace
= true;
385 dd
->dim
[dimIndex
] = 0;
386 dd
->ci
[dimIndex
] = rank
;
388 // First pulse involves (arbitrary) indices 1 and 3
389 indexvec
.push_back(1);
390 indexvec
.push_back(3);
392 ind
.index
= indexvec
;
393 ind
.nsend
[nzone
+ 1] = 2;
394 ind
.nrecv
[nzone
+ 1] = 2;
395 indvec
->push_back(ind
);
397 if (dimIndex
== 0) // Add another pulse with (arbitrary) indices 4,5,7
401 indexvec
.push_back(4);
402 indexvec
.push_back(5);
403 indexvec
.push_back(7);
405 ind
.index
= indexvec
;
406 ind
.nsend
[nzone
+ 1] = 3;
407 ind
.nrecv
[nzone
+ 1] = 3;
408 indvec
->push_back(ind
);
411 dd
->comm
->cd
[dimIndex
].ind
= *indvec
;
417 /*! \brief Check results for above-defined 1D halo with 1 pulse
419 * \param [in] x Atom coordinate data array
420 * \param [in] dd Domain decomposition object
421 * \param [in] numHomeAtoms Number of home atoms
423 void checkResults1dHaloWith1Pulse(const RVec
* x
, const gmx_domdec_t
* dd
, const int numHomeAtoms
)
425 // Check results are expected from values encoded in x data
426 for (int j
= 0; j
< DIM
; j
++)
428 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
429 EXPECT_EQ(x
[numHomeAtoms
][j
], encodedValue(dd
->neighbor
[0][0], 1, j
));
430 EXPECT_EQ(x
[numHomeAtoms
+ 1][j
], encodedValue(dd
->neighbor
[0][0], 3, j
));
434 /*! \brief Check results for above-defined 1D halo with 2 pulses
436 * \param [in] x Atom coordinate data array
437 * \param [in] dd Domain decomposition object
438 * \param [in] numHomeAtoms Number of home atoms
440 void checkResults1dHaloWith2Pulses(const RVec
* x
, const gmx_domdec_t
* dd
, const int numHomeAtoms
)
442 // Check results are expected from values encoded in x data
443 for (int j
= 0; j
< DIM
; j
++)
445 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
446 EXPECT_EQ(x
[numHomeAtoms
][j
], encodedValue(dd
->neighbor
[0][0], 1, j
));
447 EXPECT_EQ(x
[numHomeAtoms
+ 1][j
], encodedValue(dd
->neighbor
[0][0], 3, j
));
448 // Second Pulse in first dim: atoms 4,5,7 from forward horizontal neighbour
449 EXPECT_EQ(x
[numHomeAtoms
+ 2][j
], encodedValue(dd
->neighbor
[0][0], 4, j
));
450 EXPECT_EQ(x
[numHomeAtoms
+ 3][j
], encodedValue(dd
->neighbor
[0][0], 5, j
));
451 EXPECT_EQ(x
[numHomeAtoms
+ 4][j
], encodedValue(dd
->neighbor
[0][0], 7, j
));
455 /*! \brief Check results for above-defined 2D halo with 1 pulse in each dimension
457 * \param [in] x Atom coordinate data array
458 * \param [in] dd Domain decomposition object
459 * \param [in] numHomeAtoms Number of home atoms
461 void checkResults2dHaloWith1PulseInEachDim(const RVec
* x
, const gmx_domdec_t
* dd
, const int numHomeAtoms
)
463 // Check results are expected from values encoded in x data
464 for (int j
= 0; j
< DIM
; j
++)
466 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
467 EXPECT_EQ(x
[numHomeAtoms
][j
], encodedValue(dd
->neighbor
[0][0], 1, j
));
468 EXPECT_EQ(x
[numHomeAtoms
+ 1][j
], encodedValue(dd
->neighbor
[0][0], 3, j
));
469 // First Pulse in second dim: atoms 1 and 3 from forward vertical neighbour
470 EXPECT_EQ(x
[numHomeAtoms
+ 2][j
], encodedValue(dd
->neighbor
[1][0], 1, j
));
471 EXPECT_EQ(x
[numHomeAtoms
+ 3][j
], encodedValue(dd
->neighbor
[1][0], 3, j
));
475 /*! \brief Check results for above-defined 2D halo with 2 pulses in the first dimension
477 * \param [in] x Atom coordinate data array
478 * \param [in] dd Domain decomposition object
479 * \param [in] numHomeAtoms Number of home atoms
481 void checkResults2dHaloWith2PulsesInDim1(const RVec
* x
, const gmx_domdec_t
* dd
, const int numHomeAtoms
)
483 // Check results are expected from values encoded in x data
484 for (int j
= 0; j
< DIM
; j
++)
486 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
487 EXPECT_EQ(x
[numHomeAtoms
][j
], encodedValue(dd
->neighbor
[0][0], 1, j
));
488 EXPECT_EQ(x
[numHomeAtoms
+ 1][j
], encodedValue(dd
->neighbor
[0][0], 3, j
));
489 // Second Pulse in first dim: atoms 4,5,7 from forward horizontal neighbour
490 EXPECT_EQ(x
[numHomeAtoms
+ 2][j
], encodedValue(dd
->neighbor
[0][0], 4, j
));
491 EXPECT_EQ(x
[numHomeAtoms
+ 3][j
], encodedValue(dd
->neighbor
[0][0], 5, j
));
492 EXPECT_EQ(x
[numHomeAtoms
+ 4][j
], encodedValue(dd
->neighbor
[0][0], 7, j
));
493 // First Pulse in second dim: atoms 1 and 3 from forward vertical neighbour
494 EXPECT_EQ(x
[numHomeAtoms
+ 5][j
], encodedValue(dd
->neighbor
[1][0], 1, j
));
495 EXPECT_EQ(x
[numHomeAtoms
+ 6][j
], encodedValue(dd
->neighbor
[1][0], 3, j
));
499 TEST(HaloExchangeTest
, Coordinates1dHaloWith1Pulse
)
504 const int numHomeAtoms
= 10;
505 const int numHaloAtoms
= 2;
506 const int numAtomsTotal
= numHomeAtoms
+ numHaloAtoms
;
507 HostVector
<RVec
> h_x
;
508 h_x
.resize(numAtomsTotal
);
510 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
515 dd
.mpi_comm_all
= MPI_COMM_WORLD
;
516 gmx_domdec_comm_t comm
;
518 dd
.unitCellInfo
.haveScrewPBC
= false;
520 DDAtomRanges atomRanges
;
521 atomRanges
.setEnd(DDAtomRanges::Type::Home
, numHomeAtoms
);
522 dd
.comm
->atomRanges
= atomRanges
;
524 define1dRankTopology(&dd
);
526 std::vector
<gmx_domdec_ind_t
> indvec
;
527 define1dHaloWith1Pulse(&dd
, &indvec
);
529 // Perform halo exchange
530 matrix box
= { { 0., 0., 0. } };
531 dd_move_x(&dd
, box
, static_cast<ArrayRef
<RVec
>>(h_x
), nullptr);
534 checkResults1dHaloWith1Pulse(h_x
.data(), &dd
, numHomeAtoms
);
536 if (GMX_GPU_CUDA
&& GMX_THREAD_MPI
) // repeat with GPU halo codepath
538 // early return if no devices are available.
539 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
544 // Re-initialize input
545 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
547 // Perform GPU halo exchange
548 gpuHalo(&dd
, box
, &h_x
, numAtomsTotal
);
551 checkResults1dHaloWith1Pulse(h_x
.data(), &dd
, numHomeAtoms
);
555 TEST(HaloExchangeTest
, Coordinates1dHaloWith2Pulses
)
560 const int numHomeAtoms
= 10;
561 const int numHaloAtoms
= 5;
562 const int numAtomsTotal
= numHomeAtoms
+ numHaloAtoms
;
563 HostVector
<RVec
> h_x
;
564 h_x
.resize(numAtomsTotal
);
566 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
571 dd
.mpi_comm_all
= MPI_COMM_WORLD
;
572 gmx_domdec_comm_t comm
;
574 dd
.unitCellInfo
.haveScrewPBC
= false;
576 DDAtomRanges atomRanges
;
577 atomRanges
.setEnd(DDAtomRanges::Type::Home
, numHomeAtoms
);
578 dd
.comm
->atomRanges
= atomRanges
;
580 define1dRankTopology(&dd
);
582 std::vector
<gmx_domdec_ind_t
> indvec
;
583 define1dHaloWith2Pulses(&dd
, &indvec
);
585 // Perform halo exchange
586 matrix box
= { { 0., 0., 0. } };
587 dd_move_x(&dd
, box
, static_cast<ArrayRef
<RVec
>>(h_x
), nullptr);
590 checkResults1dHaloWith2Pulses(h_x
.data(), &dd
, numHomeAtoms
);
592 if (GMX_GPU_CUDA
&& GMX_THREAD_MPI
) // repeat with GPU halo codepath
594 // early return if no devices are available.
595 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
600 // Re-initialize input
601 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
603 // Perform GPU halo exchange
604 gpuHalo(&dd
, box
, &h_x
, numAtomsTotal
);
607 checkResults1dHaloWith2Pulses(h_x
.data(), &dd
, numHomeAtoms
);
612 TEST(HaloExchangeTest
, Coordinates2dHaloWith1PulseInEachDim
)
617 const int numHomeAtoms
= 10;
618 const int numHaloAtoms
= 4;
619 const int numAtomsTotal
= numHomeAtoms
+ numHaloAtoms
;
620 HostVector
<RVec
> h_x
;
621 h_x
.resize(numAtomsTotal
);
623 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
628 dd
.mpi_comm_all
= MPI_COMM_WORLD
;
629 gmx_domdec_comm_t comm
;
631 dd
.unitCellInfo
.haveScrewPBC
= false;
633 DDAtomRanges atomRanges
;
634 atomRanges
.setEnd(DDAtomRanges::Type::Home
, numHomeAtoms
);
635 dd
.comm
->atomRanges
= atomRanges
;
637 define2dRankTopology(&dd
);
639 std::vector
<gmx_domdec_ind_t
> indvec
;
640 define2dHaloWith1PulseInEachDim(&dd
, &indvec
);
642 // Perform halo exchange
643 matrix box
= { { 0., 0., 0. } };
644 dd_move_x(&dd
, box
, static_cast<ArrayRef
<RVec
>>(h_x
), nullptr);
647 checkResults2dHaloWith1PulseInEachDim(h_x
.data(), &dd
, numHomeAtoms
);
649 if (GMX_GPU_CUDA
&& GMX_THREAD_MPI
) // repeat with GPU halo codepath
651 // early return if no devices are available.
652 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
657 // Re-initialize input
658 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
660 // Perform GPU halo exchange
661 gpuHalo(&dd
, box
, &h_x
, numAtomsTotal
);
664 checkResults2dHaloWith1PulseInEachDim(h_x
.data(), &dd
, numHomeAtoms
);
668 TEST(HaloExchangeTest
, Coordinates2dHaloWith2PulsesInDim1
)
673 const int numHomeAtoms
= 10;
674 const int numHaloAtoms
= 7;
675 const int numAtomsTotal
= numHomeAtoms
+ numHaloAtoms
;
676 HostVector
<RVec
> h_x
;
677 h_x
.resize(numAtomsTotal
);
679 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
684 dd
.mpi_comm_all
= MPI_COMM_WORLD
;
685 gmx_domdec_comm_t comm
;
687 dd
.unitCellInfo
.haveScrewPBC
= false;
689 DDAtomRanges atomRanges
;
690 atomRanges
.setEnd(DDAtomRanges::Type::Home
, numHomeAtoms
);
691 dd
.comm
->atomRanges
= atomRanges
;
693 define2dRankTopology(&dd
);
695 std::vector
<gmx_domdec_ind_t
> indvec
;
696 define2dHaloWith2PulsesInDim1(&dd
, &indvec
);
698 // Perform halo exchange
699 matrix box
= { { 0., 0., 0. } };
700 dd_move_x(&dd
, box
, static_cast<ArrayRef
<RVec
>>(h_x
), nullptr);
703 checkResults2dHaloWith2PulsesInDim1(h_x
.data(), &dd
, numHomeAtoms
);
705 if (GMX_GPU_CUDA
&& GMX_THREAD_MPI
) // repeat with GPU halo codepath
707 // early return if no devices are available.
708 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
713 // Re-initialize input
714 initHaloData(h_x
.data(), numHomeAtoms
, numAtomsTotal
);
716 // Perform GPU halo exchange
717 gpuHalo(&dd
, box
, &h_x
, numAtomsTotal
);
720 checkResults2dHaloWith2PulsesInDim1(h_x
.data(), &dd
, numHomeAtoms
);