From 4e9c28ad69e93689f5f1b43af27803f3c1b83a19 Mon Sep 17 00:00:00 2001 From: Andrey Prokopenko Date: Mon, 29 Apr 2024 20:15:12 -0400 Subject: [PATCH 1/5] Switch triangulated benchmark to use Kokkos::Array --- benchmarks/CMakeLists.txt | 5 +- .../triangulated_surface_distance.cpp | 78 +++++++++---------- 2 files changed, 43 insertions(+), 40 deletions(-) diff --git a/benchmarks/CMakeLists.txt b/benchmarks/CMakeLists.txt index 37c4edb53..f45a82a55 100644 --- a/benchmarks/CMakeLists.txt +++ b/benchmarks/CMakeLists.txt @@ -40,7 +40,10 @@ if(NOT WIN32) add_subdirectory(develop) add_subdirectory(union_find) endif() -add_subdirectory(triangulated_surface_distance) + +if(Kokkos_VERSION VERSION_GREATER_EQUAL 4.3.99) + add_subdirectory(triangulated_surface_distance) +endif() if (ARBORX_ENABLE_MPI) add_subdirectory(distributed_tree_driver) diff --git a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp index 7145c33d3..2a1df44e6 100644 --- a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp +++ b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp @@ -14,6 +14,7 @@ #include #include +#include #include #include @@ -46,28 +47,28 @@ auto icosahedron() vertices.push_back(Point{b, -a, 0}); vertices.push_back(Point{-b, -a, 0}); - std::vector triangles; - - triangles.insert(triangles.end(), {2, 1, 0}); - triangles.insert(triangles.end(), {1, 2, 3}); - triangles.insert(triangles.end(), {5, 4, 3}); - triangles.insert(triangles.end(), {4, 8, 3}); - triangles.insert(triangles.end(), {7, 6, 0}); - triangles.insert(triangles.end(), {6, 9, 0}); - triangles.insert(triangles.end(), {11, 10, 4}); - triangles.insert(triangles.end(), {10, 11, 6}); - triangles.insert(triangles.end(), {9, 5, 2}); - triangles.insert(triangles.end(), {5, 9, 11}); - triangles.insert(triangles.end(), {8, 7, 1}); - triangles.insert(triangles.end(), {7, 8, 10}); - triangles.insert(triangles.end(), {2, 5, 3}); - triangles.insert(triangles.end(), {8, 1, 3}); - triangles.insert(triangles.end(), {9, 2, 0}); - triangles.insert(triangles.end(), {1, 7, 0}); - triangles.insert(triangles.end(), {11, 9, 6}); - triangles.insert(triangles.end(), {7, 10, 6}); - triangles.insert(triangles.end(), {5, 11, 4}); - triangles.insert(triangles.end(), {10, 8, 4}); + std::vector> triangles; + + triangles.push_back({2, 1, 0}); + triangles.push_back({1, 2, 3}); + triangles.push_back({5, 4, 3}); + triangles.push_back({4, 8, 3}); + triangles.push_back({7, 6, 0}); + triangles.push_back({6, 9, 0}); + triangles.push_back({11, 10, 4}); + triangles.push_back({10, 11, 6}); + triangles.push_back({9, 5, 2}); + triangles.push_back({5, 9, 11}); + triangles.push_back({8, 7, 1}); + triangles.push_back({7, 8, 10}); + triangles.push_back({2, 5, 3}); + triangles.push_back({8, 1, 3}); + triangles.push_back({9, 2, 0}); + triangles.push_back({1, 7, 0}); + triangles.push_back({11, 9, 6}); + triangles.push_back({7, 10, 6}); + triangles.push_back({5, 11, 4}); + triangles.push_back({10, 8, 4}); return std::make_pair(vertices, triangles); } @@ -80,18 +81,18 @@ auto icosahedron() / \ / \ / \ /__________\ /____\/____\ */ -void subdivide(std::vector &vertices, std::vector &triangles) +void subdivide(std::vector &vertices, + std::vector> &triangles) { std::map, int> hash; - int const num_triangles = triangles.size() / 3; + int const num_triangles = triangles.size(); int vindex = vertices.size(); - std::vector new_triangles; + std::vector> new_triangles; for (int i = 0; i < num_triangles; ++i) { - int v[3] = {triangles[3 * i + 0], triangles[3 * i + 1], - triangles[3 * i + 2]}; + int v[3] = {triangles[i][0], triangles[i][1], triangles[i][2]}; int vmid[3]; for (int j = 0; j < 3; ++j) @@ -116,10 +117,10 @@ void subdivide(std::vector &vertices, std::vector &triangles) } } - new_triangles.insert(new_triangles.end(), {v[0], vmid[0], vmid[2]}); - new_triangles.insert(new_triangles.end(), {v[1], vmid[0], vmid[1]}); - new_triangles.insert(new_triangles.end(), {v[2], vmid[2], vmid[1]}); - new_triangles.insert(new_triangles.end(), {vmid[0], vmid[1], vmid[2]}); + new_triangles.push_back({v[0], vmid[0], vmid[2]}); + new_triangles.push_back({v[1], vmid[0], vmid[1]}); + new_triangles.push_back({v[2], vmid[2], vmid[1]}); + new_triangles.push_back({vmid[0], vmid[1], vmid[2]}); } triangles = new_triangles; } @@ -163,7 +164,7 @@ template struct Triangles { Kokkos::View _points; - Kokkos::View _triangle_vertices; + Kokkos::View *, MemorySpace> _triangle_vertices; }; template @@ -176,14 +177,13 @@ class ArborX::AccessTraits, ArborX::PrimitivesTag> static KOKKOS_FUNCTION auto size(Self const &self) { - return self._triangle_vertices.size() / 3; + return self._triangle_vertices.size(); } static KOKKOS_FUNCTION auto get(Self const &self, int i) { auto const &vertices = self._triangle_vertices; - return Triangle{self._points(vertices(3 * i + 0)), - self._points(vertices(3 * i + 1)), - self._points(vertices(3 * i + 2))}; + return Triangle{self._points(vertices(i)[0]), self._points(vertices(i)[1]), + self._points(vertices(i)[2])}; } }; @@ -204,7 +204,7 @@ void writeVtk(std::string const &filename, Points const &vertices, Triangles const &triangles) { int const num_vertices = vertices.size(); - int const num_elements = triangles.size() / 3; + int const num_elements = triangles.size(); constexpr int DIM = ArborX::GeometryTraits::dimension_v; @@ -228,14 +228,14 @@ void writeVtk(std::string const &filename, Points const &vertices, out << '\n'; } - int const num_cell_vertices = 3; + constexpr int num_cell_vertices = 3; out << "\nPOLYGONS " << num_elements << " " << (num_elements * (1 + num_cell_vertices)) << '\n'; for (int i = 0; i < num_elements; ++i) { out << num_cell_vertices; for (int j = 0; j < num_cell_vertices; ++j) - out << " " << triangles(i * num_cell_vertices + j); + out << " " << triangles(i)[j]; out << '\n'; } } From 6e6d550b15db7247b42341b61adb6a5d5f4517f7 Mon Sep 17 00:00:00 2001 From: Andrey Prokopenko Date: Mon, 29 Apr 2024 16:54:16 -0400 Subject: [PATCH 2/5] Rework benchmark generation to build triangles in parallel --- .../triangulated_surface_distance.cpp | 202 +++++++++++++----- 1 file changed, 150 insertions(+), 52 deletions(-) diff --git a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp index 2a1df44e6..b87ca9eb4 100644 --- a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp +++ b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp @@ -70,83 +70,157 @@ auto icosahedron() triangles.push_back({5, 11, 4}); triangles.push_back({10, 8, 4}); - return std::make_pair(vertices, triangles); + return std::make_tuple(vertices, triangles); } -/* Subdivide every triangle into four - /\ /\ - / \ / \ - / \ ---> /____\ - / \ /\ /\ - / \ / \ / \ - /__________\ /____\/____\ -*/ -void subdivide(std::vector &vertices, - std::vector> &triangles) +void convertTriangles2EdgeForm(std::vector> &triangles, + std::vector> &edges) { std::map, int> hash; int const num_triangles = triangles.size(); - - int vindex = vertices.size(); - std::vector> new_triangles; - for (int i = 0; i < num_triangles; ++i) + edges.clear(); + for (int i = 0, eindex = 0; i < num_triangles; ++i) { int v[3] = {triangles[i][0], triangles[i][1], triangles[i][2]}; - int vmid[3]; - + int e[3]; for (int j = 0; j < 3; ++j) { int vmin = std::min(v[j], v[(j + 1) % 3]); int vmax = std::max(v[j], v[(j + 1) % 3]); auto it = hash.find(std::make_pair(vmin, vmax)); - if (it != hash.end()) + if (it == hash.end()) { - vmid[j] = it->second; + edges.push_back({vmin, vmax}); + e[j] = eindex; + hash[std::make_pair(vmin, vmax)] = eindex; + ++eindex; } else { - vmid[j] = vindex; - hash[std::make_pair(vmin, vmax)] = vindex; - ++vindex; - - vertices.push_back(Point{(vertices[vmin][0] + vertices[vmax][0]) / 2, - (vertices[vmin][1] + vertices[vmax][1]) / 2, - (vertices[vmin][2] + vertices[vmax][2]) / 2}); + e[j] = it->second; } } - - new_triangles.push_back({v[0], vmid[0], vmid[2]}); - new_triangles.push_back({v[1], vmid[0], vmid[1]}); - new_triangles.push_back({v[2], vmid[2], vmid[1]}); - new_triangles.push_back({vmid[0], vmid[1], vmid[2]}); + triangles[i] = {e[0], e[1], e[2]}; } - triangles = new_triangles; } -void projectVerticesToSphere(std::vector &vertices, float radius) +template +void convertTriangles2VertexForm( + ExecutionSpace const &space, + Kokkos::View *, MemorySpace> &triangles, + Kokkos::View *, MemorySpace> const &edges) { - for (auto &v : vertices) - { - auto norm = std::sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]); - v[0] *= radius / norm; - v[1] *= radius / norm; - v[2] *= radius / norm; - } + int const num_triangles = triangles.size(); + Kokkos::parallel_for( + "Benchmark::to_vertex_form", + Kokkos::RangePolicy(space, 0, num_triangles), + KOKKOS_LAMBDA(int i) { + auto const &e0 = edges(triangles(i)[0]); + auto const &e1 = edges(triangles(i)[1]); + + KOKKOS_ASSERT(e0[0] == e1[0] || e0[0] == e1[1] || e0[1] == e1[0] || + e0[1] == e1[1]); + if (e0[0] == e1[0] || e0[1] == e1[0]) + triangles(i) = {e0[0], e0[1], e1[1]}; + else + triangles(i) = {e0[0], e0[1], e1[0]}; + }); } -auto buildTriangles(float radius, int num_refinements) +/* Subdivide every triangle into four + /\ /\ + / \ / \ + / \ ---> /____\ + / \ /\ /\ + / \ / \ / \ + /__________\ /____\/____\ +*/ +template +void subdivide(ExecutionSpace const &space, + Kokkos::View &vertices, + Kokkos::View *, MemorySpace> &edges, + Kokkos::View *, MemorySpace> &triangles) { - Kokkos::Profiling::ScopedRegion guard("Benchmark::build_triangles"); + int const num_vertices = vertices.size(); + int const num_edges = edges.size(); + int const num_triangles = triangles.size(); - auto [vertices, triangles] = icosahedron(); - for (int i = 1; i <= num_refinements; ++i) - subdivide(vertices, triangles); + Kokkos::resize(space, vertices, vertices.size() + edges.size()); - projectVerticesToSphere(vertices, radius); + // Each edge is split in two, and each triangle adds three internal edges + Kokkos::View *, MemorySpace> new_edges( + Kokkos::view_alloc(space, Kokkos::WithoutInitializing, + "Benchmark::edges"), + 2 * num_edges + 3 * num_triangles); + Kokkos::parallel_for( + "Benchmark::split_edges", + Kokkos::RangePolicy(space, 0, num_edges), + KOKKOS_LAMBDA(int i) { + int v = edges(i)[0]; + int w = edges(i)[1]; + + int new_vindex = num_vertices + i; + vertices(new_vindex) = Point{(vertices(v)[0] + vertices(w)[0]) / 2, + (vertices(v)[1] + vertices(w)[1]) / 2, + (vertices(v)[2] + vertices(w)[2]) / 2}; + new_edges(2 * i + 0) = {v, new_vindex}; + new_edges(2 * i + 1) = {w, new_vindex}; + }); + + Kokkos::View *, MemorySpace> new_triangles( + Kokkos::view_alloc(space, Kokkos::WithoutInitializing, + "Benchmark::triangles"), + 4 * num_triangles); + Kokkos::parallel_for( + "Benchmark::split_triangles", + Kokkos::RangePolicy(space, 0, num_triangles), + KOKKOS_LAMBDA(int i) { + int e[3] = {triangles(i)[0], triangles(i)[1], triangles(i)[2]}; + + int new_edges_offset = 2 * num_edges + 3 * i; + for (int j = 0; j < 3; ++j) + { + int e0 = 2 * e[j]; + int e1 = 2 * e[(j + 1) % 3]; + if (new_edges(e0)[0] == new_edges(e1 + 1)[0]) + ++e1; + else if (new_edges(e0 + 1)[0] == new_edges(e1)[0]) + ++e0; + else if (new_edges(e0 + 1)[0] == new_edges(e1 + 1)[0]) + { + ++e0; + ++e1; + } + assert(new_edges(e0)[0] == new_edges(e1)[0]); + int enew = new_edges_offset + j; + + new_edges(enew) = {new_edges(e0)[1], new_edges(e1)[1]}; + new_triangles(4 * i + j) = {e0, e1, enew}; + } + new_triangles[4 * i + 3] = {new_edges_offset + 0, new_edges_offset + 1, + new_edges_offset + 2}; + }); + edges = new_edges; + triangles = new_triangles; +} - return std::make_pair(vertices, triangles); +template +void projectVerticesToSphere(ExecutionSpace const &space, + Kokkos::View &points, + float radius) +{ + Kokkos::parallel_for( + "Benchmark::project_to_surface", + Kokkos::RangePolicy(space, 0, points.size()), + KOKKOS_LAMBDA(int i) { + auto &v = points(i); + auto norm = std::sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]); + v[0] *= radius / norm; + v[1] *= radius / norm; + v[2] *= radius / norm; + }); } template @@ -160,6 +234,32 @@ auto vec2view(std::vector const &in, std::string const &label = "") return out; } +template +auto buildTriangles(ExecutionSpace const &space, float radius, + int num_refinements) +{ + Kokkos::Profiling::ScopedRegion guard("Benchmark::build_triangles"); + + auto [vertices_v, triangles_v] = icosahedron(); + + // Convert to edge form + std::vector> edges_v; + convertTriangles2EdgeForm(triangles_v, edges_v); + + auto vertices = vec2view(vertices_v); + auto edges = vec2view(edges_v); + auto triangles = vec2view(triangles_v); + + for (int i = 1; i <= num_refinements; ++i) + subdivide(space, vertices, edges, triangles); + + projectVerticesToSphere(space, vertices, radius); + + convertTriangles2VertexForm(space, triangles, edges); + + return std::make_pair(vertices, triangles); +} + template struct Triangles { @@ -278,16 +378,14 @@ int main(int argc, char *argv[]) return 1; } - auto [vertices_v, triangles_v] = buildTriangles(radius, num_refinements); + ExecutionSpace space; - auto vertices = vec2view(vertices_v); - auto triangles = vec2view(triangles_v); + auto [vertices, triangles] = + buildTriangles(space, radius, num_refinements); if (!vtk_filename.empty()) writeVtk(vtk_filename, vertices, triangles); - ExecutionSpace space; - Kokkos::Profiling::pushRegion("Benchmark::build_points"); Kokkos::View random_points( Kokkos::view_alloc(space, Kokkos::WithoutInitializing, From 3e4d820c412d1087684fd86da356fa13d5286a24 Mon Sep 17 00:00:00 2001 From: Andrey Prokopenko Date: Tue, 30 Apr 2024 11:56:35 -0400 Subject: [PATCH 3/5] Move out generator --- .../generator.hpp | 252 ++++++++++++++++++ .../triangulated_surface_distance.cpp | 236 +--------------- 2 files changed, 254 insertions(+), 234 deletions(-) create mode 100644 benchmarks/triangulated_surface_distance/generator.hpp diff --git a/benchmarks/triangulated_surface_distance/generator.hpp b/benchmarks/triangulated_surface_distance/generator.hpp new file mode 100644 index 000000000..62b4e2cf1 --- /dev/null +++ b/benchmarks/triangulated_surface_distance/generator.hpp @@ -0,0 +1,252 @@ +/**************************************************************************** + * Copyright (c) 2024 by the ArborX authors * + * All rights reserved. * + * * + * This file is part of the ArborX library. ArborX is * + * distributed under a BSD 3-clause license. For the licensing terms see * + * the LICENSE file in the top-level directory. * + * * + * SPDX-License-Identifier: BSD-3-Clause * + ****************************************************************************/ + +#include + +#include +#include + +static auto icosahedron() +{ + auto a = Kokkos::numbers::phi_v; + auto b = 1.f; + + using Point = ArborX::ExperimentalHyperGeometry::Point<3>; + + std::vector vertices; + + vertices.push_back(Point{0, b, -a}); + vertices.push_back(Point{b, a, 0}); + vertices.push_back(Point{-b, a, 0}); + vertices.push_back(Point{0, b, a}); + vertices.push_back(Point{0, -b, a}); + vertices.push_back(Point{-a, 0, b}); + vertices.push_back(Point{0, -b, -a}); + vertices.push_back(Point{a, 0, -b}); + vertices.push_back(Point{a, 0, b}); + vertices.push_back(Point{-a, 0, -b}); + vertices.push_back(Point{b, -a, 0}); + vertices.push_back(Point{-b, -a, 0}); + + std::vector> triangles; + + triangles.push_back({2, 1, 0}); + triangles.push_back({1, 2, 3}); + triangles.push_back({5, 4, 3}); + triangles.push_back({4, 8, 3}); + triangles.push_back({7, 6, 0}); + triangles.push_back({6, 9, 0}); + triangles.push_back({11, 10, 4}); + triangles.push_back({10, 11, 6}); + triangles.push_back({9, 5, 2}); + triangles.push_back({5, 9, 11}); + triangles.push_back({8, 7, 1}); + triangles.push_back({7, 8, 10}); + triangles.push_back({2, 5, 3}); + triangles.push_back({8, 1, 3}); + triangles.push_back({9, 2, 0}); + triangles.push_back({1, 7, 0}); + triangles.push_back({11, 9, 6}); + triangles.push_back({7, 10, 6}); + triangles.push_back({5, 11, 4}); + triangles.push_back({10, 8, 4}); + + return std::make_tuple(vertices, triangles); +} + +void convertTriangles2EdgeForm(std::vector> &triangles, + std::vector> &edges) +{ + std::map, int> hash; + + edges.clear(); + for (int i = 0, eindex = 0; i < (int)triangles.size(); ++i) + { + int e[3]; + for (int j = 0; j < 3; ++j) + { + auto minmax_pair = + std::minmax(triangles[i][j], triangles[i][(j + 1) % 3]); + auto it = hash.find(minmax_pair); + if (it == hash.end()) + { + edges.push_back({minmax_pair.first, minmax_pair.second}); + hash[minmax_pair] = eindex; + e[j] = eindex++; + } + else + { + e[j] = it->second; + } + } + triangles[i] = {e[0], e[1], e[2]}; + } +} + +template +void convertTriangles2VertexForm( + ExecutionSpace const &space, + Kokkos::View *, MemorySpace> &triangles, + Kokkos::View *, MemorySpace> const &edges) +{ + int const num_triangles = triangles.size(); + Kokkos::parallel_for( + "Benchmark::to_vertex_form", + Kokkos::RangePolicy(space, 0, num_triangles), + KOKKOS_LAMBDA(int i) { + auto const &e0 = edges(triangles(i)[0]); + auto const &e1 = edges(triangles(i)[1]); + + KOKKOS_ASSERT(e0[0] == e1[0] || e0[0] == e1[1] || e0[1] == e1[0] || + e0[1] == e1[1]); + if (e0[0] == e1[0] || e0[1] == e1[0]) + triangles(i) = {e0[0], e0[1], e1[1]}; + else + triangles(i) = {e0[0], e0[1], e1[0]}; + }); +} + +/* Subdivide every triangle into four + /\ /\ + / \ / \ + / \ ---> /____\ + / \ /\ /\ + / \ / \ / \ + /__________\ /____\/____\ +*/ +template +void subdivide(ExecutionSpace const &space, + Kokkos::View *, + MemorySpace> &vertices, + Kokkos::View *, MemorySpace> &edges, + Kokkos::View *, MemorySpace> &triangles) +{ + using Point = ArborX::ExperimentalHyperGeometry::Point<3>; + + int const num_vertices = vertices.size(); + int const num_edges = edges.size(); + int const num_triangles = triangles.size(); + + Kokkos::resize(space, vertices, vertices.size() + edges.size()); + + // Each edge is split in two, and each triangle adds three internal edges + Kokkos::View *, MemorySpace> new_edges( + Kokkos::view_alloc(space, Kokkos::WithoutInitializing, + "Benchmark::edges"), + 2 * num_edges + 3 * num_triangles); + Kokkos::parallel_for( + "Benchmark::split_edges", + Kokkos::RangePolicy(space, 0, num_edges), + KOKKOS_LAMBDA(int i) { + int v = edges(i)[0]; + int w = edges(i)[1]; + + int new_vindex = num_vertices + i; + vertices(new_vindex) = Point{(vertices(v)[0] + vertices(w)[0]) / 2, + (vertices(v)[1] + vertices(w)[1]) / 2, + (vertices(v)[2] + vertices(w)[2]) / 2}; + new_edges(2 * i + 0) = {v, new_vindex}; + new_edges(2 * i + 1) = {w, new_vindex}; + }); + + Kokkos::View *, MemorySpace> new_triangles( + Kokkos::view_alloc(space, Kokkos::WithoutInitializing, + "Benchmark::triangles"), + 4 * num_triangles); + Kokkos::parallel_for( + "Benchmark::split_triangles", + Kokkos::RangePolicy(space, 0, num_triangles), + KOKKOS_LAMBDA(int i) { + int e[3] = {triangles(i)[0], triangles(i)[1], triangles(i)[2]}; + + int new_edges_offset = 2 * num_edges + 3 * i; + for (int j = 0; j < 3; ++j) + { + // Find indices of the new edges that share a vertex + int e0 = 2 * e[j]; + int e1 = 2 * e[(j + 1) % 3]; + if (new_edges(e0)[0] == new_edges(e1 + 1)[0]) + ++e1; + else if (new_edges(e0 + 1)[0] == new_edges(e1)[0]) + ++e0; + else if (new_edges(e0 + 1)[0] == new_edges(e1 + 1)[0]) + { + ++e0; + ++e1; + } + assert(new_edges(e0)[0] == new_edges(e1)[0]); + int enew = new_edges_offset + j; + + new_edges(enew) = {new_edges(e0)[1], new_edges(e1)[1]}; + new_triangles(4 * i + j) = {e0, e1, enew}; + } + new_triangles[4 * i + 3] = {new_edges_offset + 0, new_edges_offset + 1, + new_edges_offset + 2}; + }); + edges = new_edges; + triangles = new_triangles; +} + +template +void projectVerticesToSphere( + ExecutionSpace const &space, + Kokkos::View *, MemorySpace> + &points, + float radius) +{ + Kokkos::parallel_for( + "Benchmark::project_to_surface", + Kokkos::RangePolicy(space, 0, points.size()), + KOKKOS_LAMBDA(int i) { + auto &v = points(i); + auto norm = std::sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]); + v[0] *= radius / norm; + v[1] *= radius / norm; + v[2] *= radius / norm; + }); +} + +template +auto vec2view(std::vector const &in, std::string const &label = "") +{ + Kokkos::View out( + Kokkos::view_alloc(label, Kokkos::WithoutInitializing), in.size()); + Kokkos::deep_copy(out, Kokkos::View>{ + in.data(), in.size()}); + return out; +} + +template +auto buildTriangles(ExecutionSpace const &space, float radius, + int num_refinements) +{ + Kokkos::Profiling::ScopedRegion guard("Benchmark::build_triangles"); + + auto [vertices_v, triangles_v] = icosahedron(); + + // Convert to edge form + std::vector> edges_v; + convertTriangles2EdgeForm(triangles_v, edges_v); + + auto vertices = vec2view(vertices_v); + auto edges = vec2view(edges_v); + auto triangles = vec2view(triangles_v); + + for (int i = 1; i <= num_refinements; ++i) + subdivide(space, vertices, edges, triangles); + + projectVerticesToSphere(space, vertices, radius); + + convertTriangles2VertexForm(space, triangles, edges); + + return std::make_pair(vertices, triangles); +} diff --git a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp index b87ca9eb4..c8aecd9dc 100644 --- a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp +++ b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp @@ -15,7 +15,6 @@ #include #include -#include #include #include @@ -24,242 +23,11 @@ #include #include +#include "generator.hpp" + using Point = ArborX::ExperimentalHyperGeometry::Point<3>; using Triangle = ArborX::ExperimentalHyperGeometry::Triangle<3>; -auto icosahedron() -{ - auto a = Kokkos::numbers::phi_v; - auto b = 1.f; - - std::vector vertices; - - vertices.push_back(Point{0, b, -a}); - vertices.push_back(Point{b, a, 0}); - vertices.push_back(Point{-b, a, 0}); - vertices.push_back(Point{0, b, a}); - vertices.push_back(Point{0, -b, a}); - vertices.push_back(Point{-a, 0, b}); - vertices.push_back(Point{0, -b, -a}); - vertices.push_back(Point{a, 0, -b}); - vertices.push_back(Point{a, 0, b}); - vertices.push_back(Point{-a, 0, -b}); - vertices.push_back(Point{b, -a, 0}); - vertices.push_back(Point{-b, -a, 0}); - - std::vector> triangles; - - triangles.push_back({2, 1, 0}); - triangles.push_back({1, 2, 3}); - triangles.push_back({5, 4, 3}); - triangles.push_back({4, 8, 3}); - triangles.push_back({7, 6, 0}); - triangles.push_back({6, 9, 0}); - triangles.push_back({11, 10, 4}); - triangles.push_back({10, 11, 6}); - triangles.push_back({9, 5, 2}); - triangles.push_back({5, 9, 11}); - triangles.push_back({8, 7, 1}); - triangles.push_back({7, 8, 10}); - triangles.push_back({2, 5, 3}); - triangles.push_back({8, 1, 3}); - triangles.push_back({9, 2, 0}); - triangles.push_back({1, 7, 0}); - triangles.push_back({11, 9, 6}); - triangles.push_back({7, 10, 6}); - triangles.push_back({5, 11, 4}); - triangles.push_back({10, 8, 4}); - - return std::make_tuple(vertices, triangles); -} - -void convertTriangles2EdgeForm(std::vector> &triangles, - std::vector> &edges) -{ - std::map, int> hash; - - int const num_triangles = triangles.size(); - edges.clear(); - for (int i = 0, eindex = 0; i < num_triangles; ++i) - { - int v[3] = {triangles[i][0], triangles[i][1], triangles[i][2]}; - int e[3]; - for (int j = 0; j < 3; ++j) - { - int vmin = std::min(v[j], v[(j + 1) % 3]); - int vmax = std::max(v[j], v[(j + 1) % 3]); - - auto it = hash.find(std::make_pair(vmin, vmax)); - if (it == hash.end()) - { - edges.push_back({vmin, vmax}); - e[j] = eindex; - hash[std::make_pair(vmin, vmax)] = eindex; - ++eindex; - } - else - { - e[j] = it->second; - } - } - triangles[i] = {e[0], e[1], e[2]}; - } -} - -template -void convertTriangles2VertexForm( - ExecutionSpace const &space, - Kokkos::View *, MemorySpace> &triangles, - Kokkos::View *, MemorySpace> const &edges) -{ - int const num_triangles = triangles.size(); - Kokkos::parallel_for( - "Benchmark::to_vertex_form", - Kokkos::RangePolicy(space, 0, num_triangles), - KOKKOS_LAMBDA(int i) { - auto const &e0 = edges(triangles(i)[0]); - auto const &e1 = edges(triangles(i)[1]); - - KOKKOS_ASSERT(e0[0] == e1[0] || e0[0] == e1[1] || e0[1] == e1[0] || - e0[1] == e1[1]); - if (e0[0] == e1[0] || e0[1] == e1[0]) - triangles(i) = {e0[0], e0[1], e1[1]}; - else - triangles(i) = {e0[0], e0[1], e1[0]}; - }); -} - -/* Subdivide every triangle into four - /\ /\ - / \ / \ - / \ ---> /____\ - / \ /\ /\ - / \ / \ / \ - /__________\ /____\/____\ -*/ -template -void subdivide(ExecutionSpace const &space, - Kokkos::View &vertices, - Kokkos::View *, MemorySpace> &edges, - Kokkos::View *, MemorySpace> &triangles) -{ - int const num_vertices = vertices.size(); - int const num_edges = edges.size(); - int const num_triangles = triangles.size(); - - Kokkos::resize(space, vertices, vertices.size() + edges.size()); - - // Each edge is split in two, and each triangle adds three internal edges - Kokkos::View *, MemorySpace> new_edges( - Kokkos::view_alloc(space, Kokkos::WithoutInitializing, - "Benchmark::edges"), - 2 * num_edges + 3 * num_triangles); - Kokkos::parallel_for( - "Benchmark::split_edges", - Kokkos::RangePolicy(space, 0, num_edges), - KOKKOS_LAMBDA(int i) { - int v = edges(i)[0]; - int w = edges(i)[1]; - - int new_vindex = num_vertices + i; - vertices(new_vindex) = Point{(vertices(v)[0] + vertices(w)[0]) / 2, - (vertices(v)[1] + vertices(w)[1]) / 2, - (vertices(v)[2] + vertices(w)[2]) / 2}; - new_edges(2 * i + 0) = {v, new_vindex}; - new_edges(2 * i + 1) = {w, new_vindex}; - }); - - Kokkos::View *, MemorySpace> new_triangles( - Kokkos::view_alloc(space, Kokkos::WithoutInitializing, - "Benchmark::triangles"), - 4 * num_triangles); - Kokkos::parallel_for( - "Benchmark::split_triangles", - Kokkos::RangePolicy(space, 0, num_triangles), - KOKKOS_LAMBDA(int i) { - int e[3] = {triangles(i)[0], triangles(i)[1], triangles(i)[2]}; - - int new_edges_offset = 2 * num_edges + 3 * i; - for (int j = 0; j < 3; ++j) - { - int e0 = 2 * e[j]; - int e1 = 2 * e[(j + 1) % 3]; - if (new_edges(e0)[0] == new_edges(e1 + 1)[0]) - ++e1; - else if (new_edges(e0 + 1)[0] == new_edges(e1)[0]) - ++e0; - else if (new_edges(e0 + 1)[0] == new_edges(e1 + 1)[0]) - { - ++e0; - ++e1; - } - assert(new_edges(e0)[0] == new_edges(e1)[0]); - int enew = new_edges_offset + j; - - new_edges(enew) = {new_edges(e0)[1], new_edges(e1)[1]}; - new_triangles(4 * i + j) = {e0, e1, enew}; - } - new_triangles[4 * i + 3] = {new_edges_offset + 0, new_edges_offset + 1, - new_edges_offset + 2}; - }); - edges = new_edges; - triangles = new_triangles; -} - -template -void projectVerticesToSphere(ExecutionSpace const &space, - Kokkos::View &points, - float radius) -{ - Kokkos::parallel_for( - "Benchmark::project_to_surface", - Kokkos::RangePolicy(space, 0, points.size()), - KOKKOS_LAMBDA(int i) { - auto &v = points(i); - auto norm = std::sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]); - v[0] *= radius / norm; - v[1] *= radius / norm; - v[2] *= radius / norm; - }); -} - -template -auto vec2view(std::vector const &in, std::string const &label = "") -{ - Kokkos::View out( - Kokkos::view_alloc(label, Kokkos::WithoutInitializing), in.size()); - Kokkos::deep_copy(out, Kokkos::View>{ - in.data(), in.size()}); - return out; -} - -template -auto buildTriangles(ExecutionSpace const &space, float radius, - int num_refinements) -{ - Kokkos::Profiling::ScopedRegion guard("Benchmark::build_triangles"); - - auto [vertices_v, triangles_v] = icosahedron(); - - // Convert to edge form - std::vector> edges_v; - convertTriangles2EdgeForm(triangles_v, edges_v); - - auto vertices = vec2view(vertices_v); - auto edges = vec2view(edges_v); - auto triangles = vec2view(triangles_v); - - for (int i = 1; i <= num_refinements; ++i) - subdivide(space, vertices, edges, triangles); - - projectVerticesToSphere(space, vertices, radius); - - convertTriangles2VertexForm(space, triangles, edges); - - return std::make_pair(vertices, triangles); -} - template struct Triangles { From a66bace969c03807fef2d146bd69a6676f51ed77 Mon Sep 17 00:00:00 2001 From: Andrey Prokopenko Date: Thu, 2 May 2024 17:08:14 -0400 Subject: [PATCH 4/5] Use poor man's Array --- benchmarks/CMakeLists.txt | 5 +-- .../generator.hpp | 39 +++++++++++++------ .../triangulated_surface_distance.cpp | 3 +- 3 files changed, 30 insertions(+), 17 deletions(-) diff --git a/benchmarks/CMakeLists.txt b/benchmarks/CMakeLists.txt index f45a82a55..37c4edb53 100644 --- a/benchmarks/CMakeLists.txt +++ b/benchmarks/CMakeLists.txt @@ -40,10 +40,7 @@ if(NOT WIN32) add_subdirectory(develop) add_subdirectory(union_find) endif() - -if(Kokkos_VERSION VERSION_GREATER_EQUAL 4.3.99) - add_subdirectory(triangulated_surface_distance) -endif() +add_subdirectory(triangulated_surface_distance) if (ARBORX_ENABLE_MPI) add_subdirectory(distributed_tree_driver) diff --git a/benchmarks/triangulated_surface_distance/generator.hpp b/benchmarks/triangulated_surface_distance/generator.hpp index 62b4e2cf1..74bbdc13b 100644 --- a/benchmarks/triangulated_surface_distance/generator.hpp +++ b/benchmarks/triangulated_surface_distance/generator.hpp @@ -11,9 +11,26 @@ #include -#include #include +#if KOKKOS_VERSION >= 40400 +#include +using Array = Kokkos::Array; +#else +template +struct Array +{ + static constexpr auto size() { return N; } + KOKKOS_FUNCTION constexpr T &operator[](int i) { return _data[i]; } + KOKKOS_FUNCTION constexpr T const &operator[](int i) const + { + return _data[i]; + } + + T _data[N]; +}; +#endif + static auto icosahedron() { auto a = Kokkos::numbers::phi_v; @@ -36,7 +53,7 @@ static auto icosahedron() vertices.push_back(Point{b, -a, 0}); vertices.push_back(Point{-b, -a, 0}); - std::vector> triangles; + std::vector> triangles; triangles.push_back({2, 1, 0}); triangles.push_back({1, 2, 3}); @@ -62,8 +79,8 @@ static auto icosahedron() return std::make_tuple(vertices, triangles); } -void convertTriangles2EdgeForm(std::vector> &triangles, - std::vector> &edges) +void convertTriangles2EdgeForm(std::vector> &triangles, + std::vector> &edges) { std::map, int> hash; @@ -94,8 +111,8 @@ void convertTriangles2EdgeForm(std::vector> &triangles, template void convertTriangles2VertexForm( ExecutionSpace const &space, - Kokkos::View *, MemorySpace> &triangles, - Kokkos::View *, MemorySpace> const &edges) + Kokkos::View *, MemorySpace> &triangles, + Kokkos::View *, MemorySpace> const &edges) { int const num_triangles = triangles.size(); Kokkos::parallel_for( @@ -126,8 +143,8 @@ template void subdivide(ExecutionSpace const &space, Kokkos::View *, MemorySpace> &vertices, - Kokkos::View *, MemorySpace> &edges, - Kokkos::View *, MemorySpace> &triangles) + Kokkos::View *, MemorySpace> &edges, + Kokkos::View *, MemorySpace> &triangles) { using Point = ArborX::ExperimentalHyperGeometry::Point<3>; @@ -138,7 +155,7 @@ void subdivide(ExecutionSpace const &space, Kokkos::resize(space, vertices, vertices.size() + edges.size()); // Each edge is split in two, and each triangle adds three internal edges - Kokkos::View *, MemorySpace> new_edges( + Kokkos::View *, MemorySpace> new_edges( Kokkos::view_alloc(space, Kokkos::WithoutInitializing, "Benchmark::edges"), 2 * num_edges + 3 * num_triangles); @@ -157,7 +174,7 @@ void subdivide(ExecutionSpace const &space, new_edges(2 * i + 1) = {w, new_vindex}; }); - Kokkos::View *, MemorySpace> new_triangles( + Kokkos::View *, MemorySpace> new_triangles( Kokkos::view_alloc(space, Kokkos::WithoutInitializing, "Benchmark::triangles"), 4 * num_triangles); @@ -234,7 +251,7 @@ auto buildTriangles(ExecutionSpace const &space, float radius, auto [vertices_v, triangles_v] = icosahedron(); // Convert to edge form - std::vector> edges_v; + std::vector> edges_v; convertTriangles2EdgeForm(triangles_v, edges_v); auto vertices = vec2view(vertices_v); diff --git a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp index c8aecd9dc..66610c514 100644 --- a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp +++ b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp @@ -14,7 +14,6 @@ #include #include -#include #include #include @@ -32,7 +31,7 @@ template struct Triangles { Kokkos::View _points; - Kokkos::View *, MemorySpace> _triangle_vertices; + Kokkos::View *, MemorySpace> _triangle_vertices; }; template From b4cb7eede4c0b95512ce3c8c70abe9d89e7a24cc Mon Sep 17 00:00:00 2001 From: Andrey Prokopenko Date: Thu, 2 May 2024 19:16:38 -0400 Subject: [PATCH 5/5] Few more fixes --- .../generator.hpp | 32 +++++++++---------- .../triangulated_surface_distance.cpp | 2 +- 2 files changed, 17 insertions(+), 17 deletions(-) diff --git a/benchmarks/triangulated_surface_distance/generator.hpp b/benchmarks/triangulated_surface_distance/generator.hpp index 74bbdc13b..0e9fb8c72 100644 --- a/benchmarks/triangulated_surface_distance/generator.hpp +++ b/benchmarks/triangulated_surface_distance/generator.hpp @@ -14,11 +14,11 @@ #include #if KOKKOS_VERSION >= 40400 -#include -using Array = Kokkos::Array; +#include +using KokkosArray = Kokkos::Kokkos; #else template -struct Array +struct KokkosArray { static constexpr auto size() { return N; } KOKKOS_FUNCTION constexpr T &operator[](int i) { return _data[i]; } @@ -31,7 +31,7 @@ struct Array }; #endif -static auto icosahedron() +auto icosahedron() { auto a = Kokkos::numbers::phi_v; auto b = 1.f; @@ -53,7 +53,7 @@ static auto icosahedron() vertices.push_back(Point{b, -a, 0}); vertices.push_back(Point{-b, -a, 0}); - std::vector> triangles; + std::vector> triangles; triangles.push_back({2, 1, 0}); triangles.push_back({1, 2, 3}); @@ -79,8 +79,8 @@ static auto icosahedron() return std::make_tuple(vertices, triangles); } -void convertTriangles2EdgeForm(std::vector> &triangles, - std::vector> &edges) +void convertTriangles2EdgeForm(std::vector> &edges, + std::vector> &triangles) { std::map, int> hash; @@ -111,8 +111,8 @@ void convertTriangles2EdgeForm(std::vector> &triangles, template void convertTriangles2VertexForm( ExecutionSpace const &space, - Kokkos::View *, MemorySpace> &triangles, - Kokkos::View *, MemorySpace> const &edges) + Kokkos::View *, MemorySpace> const &edges, + Kokkos::View *, MemorySpace> &triangles) { int const num_triangles = triangles.size(); Kokkos::parallel_for( @@ -143,8 +143,8 @@ template void subdivide(ExecutionSpace const &space, Kokkos::View *, MemorySpace> &vertices, - Kokkos::View *, MemorySpace> &edges, - Kokkos::View *, MemorySpace> &triangles) + Kokkos::View *, MemorySpace> &edges, + Kokkos::View *, MemorySpace> &triangles) { using Point = ArborX::ExperimentalHyperGeometry::Point<3>; @@ -155,7 +155,7 @@ void subdivide(ExecutionSpace const &space, Kokkos::resize(space, vertices, vertices.size() + edges.size()); // Each edge is split in two, and each triangle adds three internal edges - Kokkos::View *, MemorySpace> new_edges( + Kokkos::View *, MemorySpace> new_edges( Kokkos::view_alloc(space, Kokkos::WithoutInitializing, "Benchmark::edges"), 2 * num_edges + 3 * num_triangles); @@ -174,7 +174,7 @@ void subdivide(ExecutionSpace const &space, new_edges(2 * i + 1) = {w, new_vindex}; }); - Kokkos::View *, MemorySpace> new_triangles( + Kokkos::View *, MemorySpace> new_triangles( Kokkos::view_alloc(space, Kokkos::WithoutInitializing, "Benchmark::triangles"), 4 * num_triangles); @@ -251,8 +251,8 @@ auto buildTriangles(ExecutionSpace const &space, float radius, auto [vertices_v, triangles_v] = icosahedron(); // Convert to edge form - std::vector> edges_v; - convertTriangles2EdgeForm(triangles_v, edges_v); + std::vector> edges_v; + convertTriangles2EdgeForm(edges_v, triangles_v); auto vertices = vec2view(vertices_v); auto edges = vec2view(edges_v); @@ -263,7 +263,7 @@ auto buildTriangles(ExecutionSpace const &space, float radius, projectVerticesToSphere(space, vertices, radius); - convertTriangles2VertexForm(space, triangles, edges); + convertTriangles2VertexForm(space, edges, triangles); return std::make_pair(vertices, triangles); } diff --git a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp index 66610c514..f30db26d3 100644 --- a/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp +++ b/benchmarks/triangulated_surface_distance/triangulated_surface_distance.cpp @@ -31,7 +31,7 @@ template struct Triangles { Kokkos::View _points; - Kokkos::View *, MemorySpace> _triangle_vertices; + Kokkos::View *, MemorySpace> _triangle_vertices; }; template