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 dd->neighbor[0][0] = (rank + 1) % 4;
196 dd->neighbor[0][1] = (rank == 0) ? 3 : rank - 1;
199 /*! \brief Define 2D rank topology with 4 MPI tasks
206 * \param [in] dd Domain decomposition object
208 void define2dRankTopology(gmx_domdec_t* dd)
212 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
217 dd->neighbor[0][0] = 1;
218 dd->neighbor[0][1] = 1;
219 dd->neighbor[1][0] = 2;
220 dd->neighbor[1][1] = 2;
223 dd->neighbor[0][0] = 0;
224 dd->neighbor[0][1] = 0;
225 dd->neighbor[1][0] = 3;
226 dd->neighbor[1][1] = 3;
229 dd->neighbor[0][0] = 3;
230 dd->neighbor[0][1] = 3;
231 dd->neighbor[1][0] = 0;
232 dd->neighbor[1][1] = 0;
235 dd->neighbor[0][0] = 2;
236 dd->neighbor[0][1] = 2;
237 dd->neighbor[1][0] = 1;
238 dd->neighbor[1][1] = 1;
243 /*! \brief Define a 1D halo with 1 pulses
245 * \param [in] dd Domain decomposition object
246 * \param [in] indvec Vector of index vectors
248 void define1dHaloWith1Pulse(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
252 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
254 std::vector<int> indexvec;
255 gmx_domdec_ind_t ind;
261 // Set up indices involved in halo
265 dd->comm->cd[dimIndex].receiveInPlace = true;
266 dd->dim[dimIndex] = 0;
267 dd->ci[dimIndex] = rank;
269 // First pulse involves (arbitrary) indices 1 and 3
270 indexvec.push_back(1);
271 indexvec.push_back(3);
273 ind.index = indexvec;
274 ind.nsend[nzone + 1] = 2;
275 ind.nrecv[nzone + 1] = 2;
276 indvec->push_back(ind);
278 dd->comm->cd[dimIndex].ind = *indvec;
281 /*! \brief Define a 1D halo with 2 pulses
283 * \param [in] dd Domain decomposition object
284 * \param [in] indvec Vector of index vectors
286 void define1dHaloWith2Pulses(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
290 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
292 std::vector<int> indexvec;
293 gmx_domdec_ind_t ind;
299 // Set up indices involved in halo
303 dd->comm->cd[dimIndex].receiveInPlace = true;
304 dd->dim[dimIndex] = 0;
305 dd->ci[dimIndex] = rank;
307 // First pulse involves (arbitrary) indices 1 and 3
308 indexvec.push_back(1);
309 indexvec.push_back(3);
311 ind.index = indexvec;
312 ind.nsend[nzone + 1] = 2;
313 ind.nrecv[nzone + 1] = 2;
314 indvec->push_back(ind);
316 // Add another pulse with (arbitrary) indices 4,5,7
319 indexvec.push_back(4);
320 indexvec.push_back(5);
321 indexvec.push_back(7);
323 ind.index = indexvec;
324 ind.nsend[nzone + 1] = 3;
325 ind.nrecv[nzone + 1] = 3;
326 indvec->push_back(ind);
328 dd->comm->cd[dimIndex].ind = *indvec;
331 /*! \brief Define a 2D halo with 1 pulse in each dimension
333 * \param [in] dd Domain decomposition object
334 * \param [in] indvec Vector of index vectors
336 void define2dHaloWith1PulseInEachDim(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
340 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
342 std::vector<int> indexvec;
343 gmx_domdec_ind_t ind;
347 for (int dimIndex = 0; dimIndex < dd->ndim; dimIndex++)
350 // Set up indices involved in halo
354 dd->comm->cd[dimIndex].receiveInPlace = true;
355 dd->dim[dimIndex] = 0;
356 dd->ci[dimIndex] = rank;
358 // Single pulse involving (arbitrary) indices 1 and 3
359 indexvec.push_back(1);
360 indexvec.push_back(3);
362 ind.index = indexvec;
363 ind.nsend[nzone + 1] = 2;
364 ind.nrecv[nzone + 1] = 2;
365 indvec->push_back(ind);
367 dd->comm->cd[dimIndex].ind = *indvec;
373 /*! \brief Define a 2D halo with 2 pulses in the first dimension
375 * \param [in] dd Domain decomposition object
376 * \param [in] indvec Vector of index vectors
378 void define2dHaloWith2PulsesInDim1(gmx_domdec_t* dd, std::vector<gmx_domdec_ind_t>* indvec)
382 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
384 std::vector<int> indexvec;
385 gmx_domdec_ind_t ind;
389 for (int dimIndex = 0; dimIndex < dd->ndim; dimIndex++)
392 // Set up indices involved in halo
396 dd->comm->cd[dimIndex].receiveInPlace = true;
397 dd->dim[dimIndex] = 0;
398 dd->ci[dimIndex] = rank;
400 // First pulse involves (arbitrary) indices 1 and 3
401 indexvec.push_back(1);
402 indexvec.push_back(3);
404 ind.index = indexvec;
405 ind.nsend[nzone + 1] = 2;
406 ind.nrecv[nzone + 1] = 2;
407 indvec->push_back(ind);
409 if (dimIndex == 0) // Add another pulse with (arbitrary) indices 4,5,7
413 indexvec.push_back(4);
414 indexvec.push_back(5);
415 indexvec.push_back(7);
417 ind.index = indexvec;
418 ind.nsend[nzone + 1] = 3;
419 ind.nrecv[nzone + 1] = 3;
420 indvec->push_back(ind);
423 dd->comm->cd[dimIndex].ind = *indvec;
429 /*! \brief Check results for above-defined 1D halo with 1 pulse
431 * \param [in] x Atom coordinate data array
432 * \param [in] dd Domain decomposition object
433 * \param [in] numHomeAtoms Number of home atoms
435 void checkResults1dHaloWith1Pulse(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
437 // Check results are expected from values encoded in x data
438 for (int j = 0; j < DIM; j++)
440 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
441 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
442 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
446 /*! \brief Check results for above-defined 1D halo with 2 pulses
448 * \param [in] x Atom coordinate data array
449 * \param [in] dd Domain decomposition object
450 * \param [in] numHomeAtoms Number of home atoms
452 void checkResults1dHaloWith2Pulses(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
454 // Check results are expected from values encoded in x data
455 for (int j = 0; j < DIM; j++)
457 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
458 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
459 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
460 // Second Pulse in first dim: atoms 4,5,7 from forward horizontal neighbour
461 EXPECT_EQ(x[numHomeAtoms + 2][j], encodedValue(dd->neighbor[0][0], 4, j));
462 EXPECT_EQ(x[numHomeAtoms + 3][j], encodedValue(dd->neighbor[0][0], 5, j));
463 EXPECT_EQ(x[numHomeAtoms + 4][j], encodedValue(dd->neighbor[0][0], 7, j));
467 /*! \brief Check results for above-defined 2D halo with 1 pulse in each dimension
469 * \param [in] x Atom coordinate data array
470 * \param [in] dd Domain decomposition object
471 * \param [in] numHomeAtoms Number of home atoms
473 void checkResults2dHaloWith1PulseInEachDim(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
475 // Check results are expected from values encoded in x data
476 for (int j = 0; j < DIM; j++)
478 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
479 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
480 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
481 // First Pulse in second dim: atoms 1 and 3 from forward vertical neighbour
482 EXPECT_EQ(x[numHomeAtoms + 2][j], encodedValue(dd->neighbor[1][0], 1, j));
483 EXPECT_EQ(x[numHomeAtoms + 3][j], encodedValue(dd->neighbor[1][0], 3, j));
487 /*! \brief Check results for above-defined 2D halo with 2 pulses in the first dimension
489 * \param [in] x Atom coordinate data array
490 * \param [in] dd Domain decomposition object
491 * \param [in] numHomeAtoms Number of home atoms
493 void checkResults2dHaloWith2PulsesInDim1(const RVec* x, const gmx_domdec_t* dd, const int numHomeAtoms)
495 // Check results are expected from values encoded in x data
496 for (int j = 0; j < DIM; j++)
498 // First Pulse in first dim: atoms 1 and 3 from forward horizontal neighbour
499 EXPECT_EQ(x[numHomeAtoms][j], encodedValue(dd->neighbor[0][0], 1, j));
500 EXPECT_EQ(x[numHomeAtoms + 1][j], encodedValue(dd->neighbor[0][0], 3, j));
501 // Second Pulse in first dim: atoms 4,5,7 from forward horizontal neighbour
502 EXPECT_EQ(x[numHomeAtoms + 2][j], encodedValue(dd->neighbor[0][0], 4, j));
503 EXPECT_EQ(x[numHomeAtoms + 3][j], encodedValue(dd->neighbor[0][0], 5, j));
504 EXPECT_EQ(x[numHomeAtoms + 4][j], encodedValue(dd->neighbor[0][0], 7, j));
505 // First Pulse in second dim: atoms 1 and 3 from forward vertical neighbour
506 EXPECT_EQ(x[numHomeAtoms + 5][j], encodedValue(dd->neighbor[1][0], 1, j));
507 EXPECT_EQ(x[numHomeAtoms + 6][j], encodedValue(dd->neighbor[1][0], 3, j));
511 TEST(HaloExchangeTest, Coordinates1dHaloWith1Pulse)
516 const int numHomeAtoms = 10;
517 const int numHaloAtoms = 2;
518 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
519 HostVector<RVec> h_x;
520 h_x.resize(numAtomsTotal);
522 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
527 dd.mpi_comm_all = MPI_COMM_WORLD;
528 gmx_domdec_comm_t comm;
530 dd.unitCellInfo.haveScrewPBC = false;
532 DDAtomRanges atomRanges;
533 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
534 dd.comm->atomRanges = atomRanges;
536 define1dRankTopology(&dd);
538 std::vector<gmx_domdec_ind_t> indvec;
539 define1dHaloWith1Pulse(&dd, &indvec);
541 // Perform halo exchange
542 matrix box = { { 0., 0., 0. } };
543 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
546 checkResults1dHaloWith1Pulse(h_x.data(), &dd, numHomeAtoms);
548 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
550 // early return if no devices are available.
551 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
556 // Re-initialize input
557 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
559 // Perform GPU halo exchange
560 gpuHalo(&dd, box, &h_x, numAtomsTotal);
563 checkResults1dHaloWith1Pulse(h_x.data(), &dd, numHomeAtoms);
567 TEST(HaloExchangeTest, Coordinates1dHaloWith2Pulses)
572 const int numHomeAtoms = 10;
573 const int numHaloAtoms = 5;
574 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
575 HostVector<RVec> h_x;
576 h_x.resize(numAtomsTotal);
578 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
583 dd.mpi_comm_all = MPI_COMM_WORLD;
584 gmx_domdec_comm_t comm;
586 dd.unitCellInfo.haveScrewPBC = false;
588 DDAtomRanges atomRanges;
589 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
590 dd.comm->atomRanges = atomRanges;
592 define1dRankTopology(&dd);
594 std::vector<gmx_domdec_ind_t> indvec;
595 define1dHaloWith2Pulses(&dd, &indvec);
597 // Perform halo exchange
598 matrix box = { { 0., 0., 0. } };
599 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
602 checkResults1dHaloWith2Pulses(h_x.data(), &dd, numHomeAtoms);
604 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
606 // early return if no devices are available.
607 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
612 // Re-initialize input
613 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
615 // Perform GPU halo exchange
616 gpuHalo(&dd, box, &h_x, numAtomsTotal);
619 checkResults1dHaloWith2Pulses(h_x.data(), &dd, numHomeAtoms);
624 TEST(HaloExchangeTest, Coordinates2dHaloWith1PulseInEachDim)
629 const int numHomeAtoms = 10;
630 const int numHaloAtoms = 4;
631 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
632 HostVector<RVec> h_x;
633 h_x.resize(numAtomsTotal);
635 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
640 dd.mpi_comm_all = MPI_COMM_WORLD;
641 gmx_domdec_comm_t comm;
643 dd.unitCellInfo.haveScrewPBC = false;
645 DDAtomRanges atomRanges;
646 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
647 dd.comm->atomRanges = atomRanges;
649 define2dRankTopology(&dd);
651 std::vector<gmx_domdec_ind_t> indvec;
652 define2dHaloWith1PulseInEachDim(&dd, &indvec);
654 // Perform halo exchange
655 matrix box = { { 0., 0., 0. } };
656 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
659 checkResults2dHaloWith1PulseInEachDim(h_x.data(), &dd, numHomeAtoms);
661 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
663 // early return if no devices are available.
664 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
669 // Re-initialize input
670 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
672 // Perform GPU halo exchange
673 gpuHalo(&dd, box, &h_x, numAtomsTotal);
676 checkResults2dHaloWith1PulseInEachDim(h_x.data(), &dd, numHomeAtoms);
680 TEST(HaloExchangeTest, Coordinates2dHaloWith2PulsesInDim1)
685 const int numHomeAtoms = 10;
686 const int numHaloAtoms = 7;
687 const int numAtomsTotal = numHomeAtoms + numHaloAtoms;
688 HostVector<RVec> h_x;
689 h_x.resize(numAtomsTotal);
691 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
696 dd.mpi_comm_all = MPI_COMM_WORLD;
697 gmx_domdec_comm_t comm;
699 dd.unitCellInfo.haveScrewPBC = false;
701 DDAtomRanges atomRanges;
702 atomRanges.setEnd(DDAtomRanges::Type::Home, numHomeAtoms);
703 dd.comm->atomRanges = atomRanges;
705 define2dRankTopology(&dd);
707 std::vector<gmx_domdec_ind_t> indvec;
708 define2dHaloWith2PulsesInDim1(&dd, &indvec);
710 // Perform halo exchange
711 matrix box = { { 0., 0., 0. } };
712 dd_move_x(&dd, box, static_cast<ArrayRef<RVec>>(h_x), nullptr);
715 checkResults2dHaloWith2PulsesInDim1(h_x.data(), &dd, numHomeAtoms);
717 if (GMX_GPU_CUDA && GMX_THREAD_MPI) // repeat with GPU halo codepath
719 // early return if no devices are available.
720 if (getTestHardwareEnvironment()->getTestDeviceList().empty())
725 // Re-initialize input
726 initHaloData(h_x.data(), numHomeAtoms, numAtomsTotal);
728 // Perform GPU halo exchange
729 gpuHalo(&dd, box, &h_x, numAtomsTotal);
732 checkResults2dHaloWith2PulsesInDim1(h_x.data(), &dd, numHomeAtoms);