Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-02-23 09:16:25

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008 
0009 #include <boost/test/unit_test.hpp>
0010 
0011 #include <Acts/Plugins/ExaTrkX/detail/ConnectedComponents.cuh>
0012 
0013 #include <filesystem>
0014 #include <fstream>
0015 #include <random>
0016 #include <set>
0017 #include <vector>
0018 
0019 #include <boost/graph/adjacency_list.hpp>
0020 #include <boost/graph/connected_components.hpp>
0021 
0022 using namespace Acts::detail;
0023 
0024 using BoostGraph =
0025     boost::adjacency_list<boost::vecS, boost::vecS, boost::undirectedS>;
0026 
0027 using Vi = std::vector<int>;
0028 
0029 Vi checkLabeling(const std::vector<int> &src, const std::vector<int> &tgt) {
0030   std::size_t numNodes = std::max(*std::max_element(src.begin(), src.end()),
0031                                   *std::max_element(tgt.begin(), tgt.end())) +
0032                          1;
0033 
0034   int *cudaSrc, *cudaTgt;
0035   BOOST_REQUIRE_EQUAL(cudaMalloc(&cudaSrc, src.size() * sizeof(int)),
0036                       cudaSuccess);
0037   BOOST_REQUIRE_EQUAL(cudaMalloc(&cudaTgt, tgt.size() * sizeof(int)),
0038                       cudaSuccess);
0039   BOOST_REQUIRE_EQUAL(cudaMemcpy(cudaSrc, src.data(), src.size() * sizeof(int),
0040                                  cudaMemcpyHostToDevice),
0041                       cudaSuccess);
0042   BOOST_REQUIRE_EQUAL(cudaMemcpy(cudaTgt, tgt.data(), src.size() * sizeof(int),
0043                                  cudaMemcpyHostToDevice),
0044                       cudaSuccess);
0045 
0046   int *cudaLabels;
0047   BOOST_REQUIRE_EQUAL(cudaMalloc(&cudaLabels, numNodes * sizeof(int)),
0048                       cudaSuccess);
0049   int *cudaLabelsNext;
0050   BOOST_REQUIRE_EQUAL(cudaMalloc(&cudaLabelsNext, numNodes * sizeof(int)),
0051                       cudaSuccess);
0052 
0053   labelConnectedComponents<<<1, 1024>>>(src.size(), cudaSrc, cudaTgt, numNodes,
0054                                         cudaLabels, cudaLabelsNext);
0055 
0056   std::vector<int> labelsFromCuda(numNodes);
0057   BOOST_REQUIRE_EQUAL(
0058       cudaMemcpy(labelsFromCuda.data(), cudaLabels, numNodes * sizeof(int),
0059                  cudaMemcpyDeviceToHost),
0060       cudaSuccess);
0061 
0062   BoostGraph G(numNodes);
0063 
0064   for (auto i = 0ul; i < src.size(); ++i) {
0065     boost::add_edge(src[i], tgt[i], G);
0066   }
0067 
0068   std::vector<std::size_t> cpuLabels(numNodes);
0069   boost::connected_components(G, &cpuLabels[0]);
0070 
0071   // print
0072   std::cout << "cpu labels:     ";
0073   for (auto i = 0ul; i < numNodes; ++i) {
0074     std::cout << cpuLabels[i] << " ";
0075   }
0076   std::cout << std::endl;
0077 
0078   std::cout << "my CUDA labels: ";
0079   for (auto i = 0ul; i < numNodes; ++i) {
0080     std::cout << labelsFromCuda[i] << " ";
0081   }
0082   std::cout << std::endl;
0083 
0084   // check systematically
0085   std::map<int, int> boostToCuda;
0086   for (auto i = 0ul; i < numNodes; ++i) {
0087     if (boostToCuda.contains(cpuLabels[i])) {
0088       BOOST_CHECK_EQUAL(labelsFromCuda[i], boostToCuda.at(cpuLabels[i]));
0089     } else {
0090       auto [it, success] =
0091           boostToCuda.insert({cpuLabels[i], labelsFromCuda[i]});
0092       BOOST_CHECK(success);
0093     }
0094   }
0095 
0096   return labelsFromCuda;
0097 }
0098 
0099 BOOST_AUTO_TEST_CASE(simple_test_1) {
0100   Vi src{0, 1, 2, 3};
0101   Vi tgt{1, 2, 3, 4};
0102   checkLabeling(src, tgt);
0103 }
0104 
0105 BOOST_AUTO_TEST_CASE(simple_test_2) {
0106   Vi src{0, 1, 2, 4, 5, 6};
0107   Vi tgt{1, 2, 3, 5, 6, 7};
0108   checkLabeling(src, tgt);
0109 }
0110 
0111 BOOST_AUTO_TEST_CASE(simple_test_3) {
0112   Vi src{4, 3, 2, 1};
0113   Vi tgt{3, 2, 1, 0};
0114   checkLabeling(src, tgt);
0115 }
0116 
0117 void testRelabeling(const Vi &labels, const Vi &refLabelMask,
0118                     const Vi &refPrefixSum, const Vi &refLabels) {
0119   dim3 blockDim = 32;
0120   dim3 gridDim = (labels.size() + blockDim.x - 1) / blockDim.x;
0121 
0122   // Copy labels to device
0123   int *cudaLabels;
0124   BOOST_REQUIRE_EQUAL(cudaMalloc(&cudaLabels, labels.size() * sizeof(int)),
0125                       cudaSuccess);
0126   BOOST_REQUIRE_EQUAL(
0127       cudaMemcpy(cudaLabels, labels.data(), labels.size() * sizeof(int),
0128                  cudaMemcpyHostToDevice),
0129       cudaSuccess);
0130 
0131   // Init label mask
0132   int *cudaLabelMask;
0133   BOOST_REQUIRE_EQUAL(cudaMalloc(&cudaLabelMask, labels.size() * sizeof(int)),
0134                       cudaSuccess);
0135   BOOST_REQUIRE_EQUAL(cudaMemset(cudaLabelMask, 0, labels.size() * sizeof(int)),
0136                       cudaSuccess);
0137 
0138   makeLabelMask<<<1, 256>>>(labels.size(), cudaLabels, cudaLabelMask);
0139   BOOST_REQUIRE_EQUAL(cudaDeviceSynchronize(), cudaSuccess);
0140 
0141   std::vector<int> labelMask(labels.size());
0142   BOOST_REQUIRE_EQUAL(
0143       cudaMemcpy(labelMask.data(), cudaLabelMask,
0144                  labelMask.size() * sizeof(int), cudaMemcpyDeviceToHost),
0145       cudaSuccess);
0146 
0147   BOOST_CHECK_EQUAL_COLLECTIONS(labelMask.begin(), labelMask.end(),
0148                                 refLabelMask.begin(), refLabelMask.end());
0149 
0150   // Prefix sum
0151   int *cudaPrefixSum;
0152   BOOST_REQUIRE_EQUAL(cudaMalloc(&cudaPrefixSum, labels.size() * sizeof(int)),
0153                       cudaSuccess);
0154   thrust::exclusive_scan(thrust::device.on(0), cudaLabelMask,
0155                          cudaLabelMask + labels.size(), cudaPrefixSum);
0156 
0157   Vi prefixSum(labels.size());
0158   BOOST_REQUIRE_EQUAL(
0159       cudaMemcpy(prefixSum.data(), cudaPrefixSum, labels.size() * sizeof(int),
0160                  cudaMemcpyDeviceToHost),
0161       cudaSuccess);
0162   BOOST_CHECK_EQUAL_COLLECTIONS(prefixSum.begin(), prefixSum.end(),
0163                                 refPrefixSum.begin(), refPrefixSum.end());
0164 
0165   // Relabel
0166   mapEdgeLabels<<<1, 256>>>(labels.size(), cudaLabels, cudaPrefixSum);
0167   BOOST_REQUIRE_EQUAL(cudaDeviceSynchronize(), cudaSuccess);
0168 
0169   std::vector<int> labelsFromCuda(labels.size());
0170   BOOST_REQUIRE_EQUAL(
0171       cudaMemcpy(labelsFromCuda.data(), cudaLabels, labels.size() * sizeof(int),
0172                  cudaMemcpyDeviceToHost),
0173       cudaSuccess);
0174 
0175   BOOST_CHECK_EQUAL_COLLECTIONS(labelsFromCuda.begin(), labelsFromCuda.end(),
0176                                 refLabels.begin(), refLabels.end());
0177 }
0178 
0179 BOOST_AUTO_TEST_CASE(test_relabeling) {
0180   // clang-format off
0181   Vi labels      {0, 3, 5, 3, 0, 0};
0182   Vi refLabelMask{1, 0, 0, 1, 0, 1};
0183   Vi refPrefixSum{0, 1, 1, 1, 2, 2};
0184   Vi refLabels   {0, 1, 2, 1, 0, 0};
0185   // clang-format on
0186 
0187   testRelabeling(labels, refLabelMask, refPrefixSum, refLabels);
0188 }
0189 
0190 BOOST_AUTO_TEST_CASE(test_relabeling_2) {
0191   // clang-format off
0192   Vi labels      {1, 3, 5, 3, 1, 1};
0193   Vi refLabelMask{0, 1, 0, 1, 0, 1};
0194   Vi refPrefixSum{0, 0, 1, 1, 2, 2};
0195   Vi refLabels   {0, 1, 2, 1, 0, 0};
0196   // clang-format on
0197 
0198   testRelabeling(labels, refLabelMask, refPrefixSum, refLabels);
0199 }
0200 
0201 auto makeRandomGraph(std::size_t nodes, std::size_t edges) {
0202   std::default_random_engine rng(2345);
0203   std::uniform_int_distribution<> dist(0, nodes);
0204   std::set<std::pair<int, int>> set;
0205   Vi src(edges), tgt(edges);
0206   for (auto n = 0ul; n < edges; ++n) {
0207     auto a = dist(rng);
0208     auto b = dist(rng);
0209     if (a == b) {
0210       continue;
0211     }
0212     auto s = std::min(a, b);
0213     auto t = std::max(a, b);
0214     auto [it, success] = set.insert({s, t});
0215     if (success) {
0216       src.at(n) = s;
0217       tgt.at(n) = t;
0218     }
0219   }
0220 
0221   return std::make_pair(src, tgt);
0222 }
0223 
0224 BOOST_AUTO_TEST_CASE(test_random_graph) {
0225   auto [src, tgt] = makeRandomGraph(5, 10);
0226   checkLabeling(src, tgt);
0227 }
0228 
0229 void testFullConnectedComponents(const Vi &src, const Vi &tgt) {
0230   const auto nNodes = std::max(*std::max_element(src.begin(), src.end()),
0231                                *std::max_element(tgt.begin(), tgt.end())) +
0232                       1;
0233 
0234   // print src and tgt
0235   /*
0236     std::cout << "src: ";
0237     for (int i = 0; i < src.size(); ++i) {
0238       std::cout << src[i] << " ";
0239     }
0240     std::cout << std::endl;
0241     std::cout << "tgt: ";
0242     for (int i = 0; i < tgt.size(); ++i) {
0243       std::cout << tgt[i] << " ";
0244     }
0245     std::cout << std::endl;
0246   */
0247   cudaStream_t stream;
0248   BOOST_REQUIRE_EQUAL(cudaStreamCreate(&stream), cudaSuccess);
0249 
0250   // copy src and tgt to device
0251   int *cudaSrc, *cudaTgt;
0252   BOOST_REQUIRE_EQUAL(
0253       cudaMallocAsync(&cudaSrc, src.size() * sizeof(int), stream), cudaSuccess);
0254   BOOST_REQUIRE_EQUAL(
0255       cudaMallocAsync(&cudaTgt, tgt.size() * sizeof(int), stream), cudaSuccess);
0256   BOOST_REQUIRE_EQUAL(
0257       cudaMemcpyAsync(cudaSrc, src.data(), src.size() * sizeof(int),
0258                       cudaMemcpyHostToDevice, stream),
0259       cudaSuccess);
0260   BOOST_REQUIRE_EQUAL(
0261       cudaMemcpyAsync(cudaTgt, tgt.data(), src.size() * sizeof(int),
0262                       cudaMemcpyHostToDevice, stream),
0263       cudaSuccess);
0264 
0265   // init label array
0266   int *cudaLabels;
0267   BOOST_REQUIRE_EQUAL(
0268       cudaMallocAsync(&cudaLabels, nNodes * sizeof(int), stream), cudaSuccess);
0269 
0270   // run connected components
0271   int cudaNumLabels = connectedComponentsCuda(src.size(), cudaSrc, cudaTgt,
0272                                               nNodes, cudaLabels, stream);
0273   BOOST_REQUIRE_EQUAL(cudaStreamSynchronize(stream), cudaSuccess);
0274 
0275   // print message from last cuda error code
0276   std::cout << "CUDA Error msg: " << cudaGetErrorString(cudaPeekAtLastError())
0277             << std::endl;
0278   BOOST_REQUIRE_EQUAL(cudaGetLastError(), cudaSuccess);
0279 
0280   // copy labels back
0281   std::vector<int> labelsFromCuda(nNodes);
0282   BOOST_REQUIRE_EQUAL(
0283       cudaMemcpyAsync(labelsFromCuda.data(), cudaLabels, nNodes * sizeof(int),
0284                       cudaMemcpyDeviceToHost, stream),
0285       cudaSuccess);
0286 
0287   BOOST_REQUIRE_EQUAL(cudaFreeAsync(cudaSrc, stream), cudaSuccess);
0288   BOOST_REQUIRE_EQUAL(cudaFreeAsync(cudaTgt, stream), cudaSuccess);
0289   BOOST_REQUIRE_EQUAL(cudaFreeAsync(cudaLabels, stream), cudaSuccess);
0290 
0291   // sync
0292   BOOST_REQUIRE_EQUAL(cudaStreamSynchronize(stream), cudaSuccess);
0293   BOOST_REQUIRE_EQUAL(cudaStreamDestroy(stream), cudaSuccess);
0294 
0295   // print labelsFromCuda
0296   /*
0297       std::cout << "CUDA labels: ";
0298       for (int i = 0; i < nNodes; ++i) {
0299         std::cout << labelsFromCuda[i] << " ";
0300       }
0301       std::cout << std::endl;
0302   */
0303   // run boost graph for comparison
0304 
0305   BoostGraph G(nNodes);
0306 
0307   for (auto i = 0ul; i < src.size(); ++i) {
0308     boost::add_edge(src[i], tgt[i], G);
0309   }
0310 
0311   std::vector<std::size_t> cpuLabels(boost::num_vertices(G));
0312   int cpuNumLabels = boost::connected_components(G, &cpuLabels[0]);
0313 
0314   // check
0315   BOOST_CHECK_EQUAL(cudaNumLabels, cpuNumLabels);
0316   BOOST_CHECK_EQUAL_COLLECTIONS(labelsFromCuda.begin(), labelsFromCuda.end(),
0317                                 cpuLabels.begin(), cpuLabels.end());
0318 }
0319 
0320 BOOST_AUTO_TEST_CASE(full_test_tiny_graph) {
0321   auto [src, tgt] = makeRandomGraph(5, 10);
0322   testFullConnectedComponents(src, tgt);
0323 }
0324 
0325 BOOST_AUTO_TEST_CASE(full_test_small_graph) {
0326   auto [src, tgt] = makeRandomGraph(100, 500);
0327   testFullConnectedComponents(src, tgt);
0328 }
0329 
0330 BOOST_AUTO_TEST_CASE(full_test_big_graph) {
0331   for (int i = 0; i < 3; ++i) {
0332     std::cout << "Test graph " << i << std::endl;
0333     auto [src, tgt] = makeRandomGraph(100'000, 500'000);
0334     testFullConnectedComponents(src, tgt);
0335   }
0336 }
0337 
0338 BOOST_AUTO_TEST_CASE(test_from_file) {
0339   if (!std::filesystem::exists("edges_cuda_trackbuilding.txt")) {
0340     std::cout << "File edges_cuda_trackbuilding.txt not found" << std::endl;
0341     return;
0342   }
0343 
0344   std::ifstream file("edges_cuda_trackbuilding.txt");
0345   std::vector<int> src, tgt;
0346   int a, b;
0347   while (file >> a >> b) {
0348     src.push_back(a);
0349     tgt.push_back(b);
0350   }
0351 
0352   testFullConnectedComponents(src, tgt);
0353 }
0354 
0355 // try this pathologic case
0356 BOOST_AUTO_TEST_CASE(special_1) {
0357   testFullConnectedComponents({1, 2}, {4, 7});
0358 }
0359 
0360 BOOST_AUTO_TEST_CASE(special_2) {
0361   Vi src{1, 2};
0362   Vi tgt{4, 7};
0363   checkLabeling(src, tgt);
0364 }
0365 
0366 BOOST_AUTO_TEST_CASE(special_3) {
0367   // clang-format off
0368   Vi labels      {0, 1, 2, 3, 1, 5, 6, 2};
0369   Vi refLabelMask{1, 1, 1, 1, 0, 1, 1, 0};
0370   Vi refPrefixSum{0, 1, 2, 3, 4, 4, 5, 6};
0371   Vi refLabels   {0, 1, 2, 3, 1, 4, 5, 2};
0372   // clang-format on
0373 
0374   testRelabeling(labels, refLabelMask, refPrefixSum, refLabels);
0375 }
0376 
0377 BOOST_AUTO_TEST_CASE(special_4) {
0378   Vi src{1, 2};
0379   Vi tgt{4, 7};
0380 
0381   auto labelsFromCuda = checkLabeling(src, tgt);
0382 
0383   Vi refLabelMask{1, 1, 1, 1, 0, 1, 1, 0};
0384   Vi refPrefixSum{0, 1, 2, 3, 4, 4, 5, 6};
0385   Vi refLabels{0, 1, 2, 3, 1, 4, 5, 2};
0386 
0387   testRelabeling(labelsFromCuda, refLabelMask, refPrefixSum, refLabels);
0388 }