From cb67959498f83a7c52ae1a9c11d315e1dc6afe79 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Mon, 27 Apr 2026 13:09:07 +0100 Subject: [PATCH 01/22] Added cuBQL as a submodule --- .gitmodules | 3 +++ vendor/cuBQL | 1 + 2 files changed, 4 insertions(+) create mode 160000 vendor/cuBQL diff --git a/.gitmodules b/.gitmodules index 77e2a79a..dd0e33bd 100644 --- a/.gitmodules +++ b/.gitmodules @@ -20,3 +20,6 @@ [submodule "vendor/GPRT"] path = vendor/GPRT url = https://github.com/gprt-org/GPRT.git +[submodule "vendor/cuBQL"] + path = vendor/cuBQL + url = https://github.com/NVIDIA/cuBQL diff --git a/vendor/cuBQL b/vendor/cuBQL new file mode 160000 index 00000000..5a651b37 --- /dev/null +++ b/vendor/cuBQL @@ -0,0 +1 @@ +Subproject commit 5a651b3787e0cfae123eb427088fa73921e5a5fe From d024cef4017de47517012417ae554b95155b5083 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Tue, 28 Apr 2026 15:05:36 +0100 Subject: [PATCH 02/22] Added required CMake wiring + a CMakePresets.json for offload flags --- CMakeLists.txt | 34 +++++++++++++++++++++++++++++++--- CMakePresets.json | 45 +++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 76 insertions(+), 3 deletions(-) create mode 100644 CMakePresets.json diff --git a/CMakeLists.txt b/CMakeLists.txt index fb762521..3b95d647 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -7,6 +7,7 @@ option(XDG_ENABLE_LIBMESH "Enable support for the libMesh mesh library" OFF option(XDG_LINK_MPI "Link with MPI (for dependency compatibility)" OFF) option(XDG_ENABLE_EMBREE "Enable support for the Embree ray tracing library" ON) option(XDG_ENABLE_GPRT "Enable support for the GPRT ray tracing library" OFF) +option(XDG_ENABLE_CUBQL "Enable support for the cuBQL ray tracing library" OFF) option(XDG_BUILD_TESTS "Enable C++ unit testing" ON) option(XDG_BUILD_TOOLS "Enable tools and miniapps" ON) @@ -20,6 +21,10 @@ if(NOT CMAKE_BUILD_TYPE) set(CMAKE_BUILD_TYPE Release CACHE STRING "Choose build type" FORCE) endif() +if(DEFINED XDG_CMAKE_PRESET) + message(STATUS "XDG CMake preset: ${XDG_CMAKE_PRESET}") +endif() + # Compiler options (things in this section may not be platform-portable) set(CMAKE_POSITION_INDEPENDENT_CODE ON) @@ -111,6 +116,12 @@ if(XDG_ENABLE_GPRT) ) endif() +if(XDG_ENABLE_CUBQL) + list(APPEND VENDOR_PATHS + vendor/cuBQL + ) +endif() + if(GIT_FOUND AND EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/.git") option(XDG_GIT_SUBMODULE "Check submodules during build" ON) if(XDG_GIT_SUBMODULE) @@ -175,11 +186,12 @@ if (NOT XDG_ENABLE_MOAB AND NOT XDG_ENABLE_LIBMESH) endif() # Ensure at least one ray tracing backend is enabled -if (NOT XDG_ENABLE_EMBREE AND NOT XDG_ENABLE_GPRT) +if (NOT XDG_ENABLE_EMBREE AND NOT XDG_ENABLE_GPRT AND NOT XDG_ENABLE_CUBQL) message(FATAL_ERROR "No ray tracing backend enabled. Enable at least one of:\n" " -DXDG_ENABLE_EMBREE=ON\n" - " -DXDG_ENABLE_GPRT=ON") + " -DXDG_ENABLE_GPRT=ON\n" + " -DXDG_ENABLE_CUBQL=ON") endif() # GPRT @@ -188,6 +200,11 @@ if (XDG_ENABLE_GPRT) add_subdirectory(vendor/GPRT) endif() +if (XDG_ENABLE_CUBQL) + set(CUBQL_OMP ON CACHE BOOL "Build cuBQL with OpenMP target offload" FORCE) + add_subdirectory(vendor/cuBQL) +endif() + list(APPEND xdg_sources src/geometry/measure.cpp src/geometry/plucker.cpp @@ -221,6 +238,12 @@ dbl_deviceCode endif() +if (XDG_ENABLE_CUBQL) +list(APPEND xdg_sources +src/cuBQL/ray_tracer.cpp +) +endif() + if (XDG_ENABLE_LIBMESH) list(APPEND xdg_sources src/libmesh/mesh_manager.cpp @@ -281,7 +304,12 @@ if (${CMAKE_BUILD_TYPE} MATCHES "Debug") endif() # attempt to find OpenMP and include it if found -find_package(OpenMP) +if (XDG_ENABLE_CUBQL) + find_package(OpenMP REQUIRED) +else() + find_package(OpenMP) +endif() + if (OpenMP_CXX_FOUND) target_link_libraries(xdg PRIVATE OpenMP::OpenMP_CXX) target_compile_definitions(xdg PRIVATE XDG_HAVE_OPENMP) diff --git a/CMakePresets.json b/CMakePresets.json new file mode 100644 index 00000000..d1e26cb8 --- /dev/null +++ b/CMakePresets.json @@ -0,0 +1,45 @@ +{ + "version": 3, + "configurePresets": [ + { + "name": "base", + "generator": "Unix Makefiles", + "binaryDir": "${sourceDir}/build/${presetName}", + "cacheVariables": { + "CMAKE_CXX_FLAGS": "$env{PRESET_CXX_FLAGS} $env{LLVM_CXX_FLAGS} $env{COMMON_CXX_FLAGS}", + "XDG_CMAKE_PRESET": "${presetName}" + }, + "environment": { + "COMMON_CXX_FLAGS": "" + } + }, + { + "name": "llvm", + "inherits": ["base"], + "cacheVariables": { + "CMAKE_C_COMPILER": "clang", + "CMAKE_CXX_COMPILER": "clang++" + }, + "environment": { + "LLVM_CXX_FLAGS": "-fopenmp -fopenmp-cuda-mode" + } + }, + { + "name": "cubql_llvm_ada", + "inherits": ["llvm"], + "displayName": "cuBQL LLVM OpenMP offload RTX 2000 Ada", + "cacheVariables": { + "XDG_ENABLE_CUBQL": "ON" + }, + "environment": { + "PRESET_CXX_FLAGS": "-fopenmp-targets=nvptx64 -Xopenmp-target -march=sm_80" + } + } + ], + "buildPresets": [ + { + "name": "cubql_llvm_ada", + "configurePreset": "cubql_llvm_ada" + } + ] +} From 5bbaf150b9bc77d2d204c5a064a6838285d4dd6a Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Tue, 28 Apr 2026 15:06:37 +0100 Subject: [PATCH 03/22] Added constructor setup for a cubql backend --- include/xdg/constants.h | 8 +++++--- include/xdg/ray_tracers.h | 6 +++++- src/xdg.cpp | 16 +++++++++++++++- 3 files changed, 25 insertions(+), 5 deletions(-) diff --git a/include/xdg/constants.h b/include/xdg/constants.h index 38df9db1..e143d537 100644 --- a/include/xdg/constants.h +++ b/include/xdg/constants.h @@ -54,7 +54,8 @@ enum class MeshLibrary { // Ray Tracing library identifier enum class RTLibrary { EMBREE, - GPRT + GPRT, + CUBQL }; static const std::map MESH_LIB_TO_STR = @@ -67,7 +68,8 @@ static const std::map MESH_LIB_TO_STR = static const std::map RT_LIB_TO_STR = { {RTLibrary::EMBREE, "EMBREE"}, - {RTLibrary::GPRT, "GPRT"} + {RTLibrary::GPRT, "GPRT"}, + {RTLibrary::CUBQL, "CUBQL"} }; // Mesh identifer type @@ -148,4 +150,4 @@ struct formatter : fmt::formatter { } -#endif // include guard \ No newline at end of file +#endif // include guard diff --git a/include/xdg/ray_tracers.h b/include/xdg/ray_tracers.h index b4816571..2031baf0 100644 --- a/include/xdg/ray_tracers.h +++ b/include/xdg/ray_tracers.h @@ -5,4 +5,8 @@ #ifdef XDG_ENABLE_GPRT #include "xdg/gprt/ray_tracer.h" -#endif \ No newline at end of file +#endif + +#ifdef XDG_ENABLE_CUBQL +#include "xdg/cuBQL/ray_tracer.h" +#endif diff --git a/src/xdg.cpp b/src/xdg.cpp index 371a9a1c..83816e6c 100644 --- a/src/xdg.cpp +++ b/src/xdg.cpp @@ -31,6 +31,14 @@ XDG::XDG(std::shared_ptr mesh_manager, RTLibrary ray_tracing_lib) #else fatal_error("This build was not compiled with GPRT support (XDG_ENABLE_GPRT=OFF)."); #endif + + case RTLibrary::CUBQL: + #ifdef XDG_ENABLE_CUBQL + set_ray_tracing_interface(std::make_shared()); + break; + #else + fatal_error("This build was not compiled with cuBQL support (XDG_ENABLE_CUBQL=OFF)."); + #endif } } @@ -84,6 +92,9 @@ std::shared_ptr XDG::create(MeshLibrary mesh_lib, RTLibrary ray_tracing_lib #ifdef XDG_ENABLE_GPRT if (ray_tracing_lib == RTLibrary::GPRT) return std::make_shared(); #endif + #ifdef XDG_ENABLE_CUBQL + if (ray_tracing_lib == RTLibrary::CUBQL) return std::make_shared(); + #endif // If no supported ray tracing library throw an error std::string msg = fmt::format("Invalid ray tracing library '{}'. Supported:", RT_LIB_TO_STR.at(ray_tracing_lib)); @@ -93,6 +104,9 @@ std::shared_ptr XDG::create(MeshLibrary mesh_lib, RTLibrary ray_tracing_lib #ifdef XDG_ENABLE_GPRT msg += " GPRT"; #endif + #ifdef XDG_ENABLE_CUBQL + msg += " CUBQL"; + #endif fatal_error(msg); }; @@ -325,4 +339,4 @@ double XDG::measure_volume_area(MeshID volume) const return area; } -} // namespace xdg \ No newline at end of file +} // namespace xdg From 42d8b947ffc7a81a0c45a514f0fadde30e9d2bbf Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Tue, 28 Apr 2026 15:07:12 +0100 Subject: [PATCH 04/22] Added cubql ray tracer header + implmentation stubs --- include/xdg/cuBQL/ray_tracer.h | 68 +++++++++++++++++++++ src/cuBQL/ray_tracer.cpp | 104 +++++++++++++++++++++++++++++++++ 2 files changed, 172 insertions(+) create mode 100644 include/xdg/cuBQL/ray_tracer.h create mode 100644 src/cuBQL/ray_tracer.cpp diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h new file mode 100644 index 00000000..eae80009 --- /dev/null +++ b/include/xdg/cuBQL/ray_tracer.h @@ -0,0 +1,68 @@ +#ifndef _XDG_CUBQL_RAY_TRACING_INTERFACE_H +#define _XDG_CUBQL_RAY_TRACING_INTERFACE_H + +#include +#include +#include + +#include "xdg/constants.h" +#include "xdg/geometry_data.h" +#include "xdg/mesh_manager_interface.h" +#include "xdg/ray.h" +#include "xdg/ray_tracing_interface.h" + + + +namespace xdg { + +class CuBQLRayTracer : public RayTracer { +public: + CuBQLRayTracer(); + ~CuBQLRayTracer() override; + + RTLibrary library() const override { return RTLibrary::CUBQL; } + + void init() override; + + std::pair + register_volume(const std::shared_ptr& mesh_manager, + MeshID volume) override; + + TreeID create_surface_tree(const std::shared_ptr& mesh_manager, + MeshID volume) override; + + TreeID create_element_tree(const std::shared_ptr& mesh_manager, + MeshID volume) override; + + void create_global_surface_tree() override; + + void create_global_element_tree() override; + + MeshID find_element(const Position& point) const override; + + MeshID find_element(TreeID tree, const Position& point) const override; + + bool point_in_volume(TreeID tree, + const Position& point, + const Direction* direction = nullptr, + const std::vector* exclude_primitives = nullptr) const override; + + std::pair ray_fire(TreeID tree, + const Position& origin, + const Direction& direction, + const double dist_limit = INFTY, + HitOrientation orientation = HitOrientation::EXITING, + std::vector* const exclude_primitives = nullptr) override; + + std::pair closest(TreeID tree, + const Position& origin) override; + + bool occluded(TreeID tree, + const Position& origin, + const Direction& direction, + double& dist) const override; +}; + +} // namespace xdg + +#endif // include guard diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp new file mode 100644 index 00000000..05de9ea7 --- /dev/null +++ b/src/cuBQL/ray_tracer.cpp @@ -0,0 +1,104 @@ +#include "xdg/cuBQL/ray_tracer.h" +#include "xdg/error.h" + +#include "cuBQL/bvh.h" +#include "cuBQL/builder/omp.h" + +namespace xdg { + +CuBQLRayTracer::CuBQLRayTracer() = default; + +CuBQLRayTracer::~CuBQLRayTracer() = default; + +void CuBQLRayTracer::init() +{ +} + +std::pair +CuBQLRayTracer::register_volume(const std::shared_ptr& mesh_manager, + MeshID volume) +{ + TreeID surface_tree = create_surface_tree(mesh_manager, volume); + TreeID element_tree = create_element_tree(mesh_manager, volume); + return {surface_tree, element_tree}; +} + +TreeID +CuBQLRayTracer::create_surface_tree(const std::shared_ptr&, + MeshID) +{ + // Testing cubql compilation + cuBQL::bvh_t bvh; + + + return TREE_NONE; +} + +TreeID +CuBQLRayTracer::create_element_tree(const std::shared_ptr&, + MeshID) +{ + fatal_error("Element trees not currently supported with cuBQL ray tracer"); + return TREE_NONE; +} + +void CuBQLRayTracer::create_global_surface_tree() +{ + fatal_error("Global surface trees not currently supported with cuBQL ray tracer"); +} + +void CuBQLRayTracer::create_global_element_tree() +{ + fatal_error("Global element trees not currently supported with cuBQL ray tracer"); +} + +MeshID CuBQLRayTracer::find_element(const Position&) const +{ + fatal_error("Element queries not currently supported with cuBQL ray tracer"); + return ID_NONE; +} + +MeshID CuBQLRayTracer::find_element(TreeID, const Position&) const +{ + fatal_error("Element queries not currently supported with cuBQL ray tracer"); + return ID_NONE; +} + +bool CuBQLRayTracer::point_in_volume(TreeID, + const Position&, + const Direction*, + const std::vector*) const +{ + fatal_error("Point-in-volume queries not currently supported with cuBQL ray tracer"); + return false; +} + +std::pair +CuBQLRayTracer::ray_fire(TreeID, + const Position&, + const Direction&, + const double, + HitOrientation, + std::vector* const) +{ + fatal_error("Ray-fire queries not currently supported with cuBQL ray tracer"); + return {INFTY, ID_NONE}; +} + +std::pair +CuBQLRayTracer::closest(TreeID, const Position&) +{ + fatal_error("Closest queries not currently supported with cuBQL ray tracer"); + return {INFTY, ID_NONE}; +} + +bool CuBQLRayTracer::occluded(TreeID, + const Position&, + const Direction&, + double&) const +{ + fatal_error("Occlusion queries not currently supported with cuBQL ray tracer"); + return false; +} + +} // namespace xdg From e3ef067917bede5f2e7772578fbc6a6bb7e36409 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Tue, 28 Apr 2026 15:36:30 +0100 Subject: [PATCH 05/22] Added a nvhpc preset --- CMakePresets.json | 23 +++++++++++++++++++++++ 1 file changed, 23 insertions(+) diff --git a/CMakePresets.json b/CMakePresets.json index d1e26cb8..09b86db1 100644 --- a/CMakePresets.json +++ b/CMakePresets.json @@ -34,12 +34,35 @@ "environment": { "PRESET_CXX_FLAGS": "-fopenmp-targets=nvptx64 -Xopenmp-target -march=sm_80" } + }, + { + "name": "nvhpc", + "inherits": ["base"], + "cacheVariables": { + "CMAKE_C_COMPILER": "nvc", + "CMAKE_CXX_COMPILER": "nvc++" + } + }, + { + "name": "cubql_nvhpc_ada", + "inherits": ["nvhpc"], + "displayName": "cuBQL NVHPC OpenMP offload RTX 2000 Ada", + "cacheVariables": { + "XDG_ENABLE_CUBQL": "ON" + }, + "environment": { + "PRESET_CXX_FLAGS": "-mp=gpu -Minfo=mp -gpu=cc89" + } } ], "buildPresets": [ { "name": "cubql_llvm_ada", "configurePreset": "cubql_llvm_ada" + }, + { + "name": "cubql_nvhpc_ada", + "configurePreset": "cubql_nvhpc_ada" } ] } From f83cefc829b078295627a887b0e32ab9ff2b6f33 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Tue, 28 Apr 2026 17:27:44 +0100 Subject: [PATCH 06/22] Added some guards to ensure that cubql path stays on OMP and doesnt switch to cuda path --- src/cuBQL/ray_tracer.cpp | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index 05de9ea7..fa9c3a02 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -1,6 +1,11 @@ #include "xdg/cuBQL/ray_tracer.h" #include "xdg/error.h" +// Guards to prevent CUDA headers from being included in host code, which causes failed compilation with LLVM-clang +#if defined(__CUDA_ARCH__) && !defined(__CUDACC__) +#undef __CUDA_ARCH__ +#endif + #include "cuBQL/bvh.h" #include "cuBQL/builder/omp.h" From 3b977852c59fd892793562aa3381c8e9351784cd Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Wed, 29 Apr 2026 16:13:17 +0100 Subject: [PATCH 07/22] Added a TODO comment to fix the loading of libomptarget.so --- CMakeLists.txt | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/CMakeLists.txt b/CMakeLists.txt index 3b95d647..00ab66e9 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -348,6 +348,22 @@ if (XDG_ENABLE_GPRT) endif() target_link_libraries(xdg PRIVATE fmt::fmt) +if (XDG_ENABLE_CUBQL) + target_compile_definitions(xdg PUBLIC XDG_ENABLE_CUBQL) + target_link_libraries(xdg $) + # TODO: Stop relying on LD_LIBRARY_PATH for LLVM OpenMP offload runtimes. + # As a temporary measure whilst figuring out our way around cuBQL this is + # okay but in the long run we should aim for a more robust solutions here. + # Clang injects libomptarget when offload flags are supplied by presets, but + # CMake does not currently add that compiler runtime directory to RPATH. + # Add targeted BUILD_RPATH handling for libomptarget/libomp, and decide + # whether install RPATH should remain environment-module based or be opt-in. +endif() + +target_link_libraries(xdg fmt::fmt) + +#TODO: refactor to only be required when tools are enabled +target_link_libraries(xdg indicators::indicators) # ========================== # Link ray tracing libraries From 250fe625fce87a32e05d420a55748a2627ce8ba5 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Thu, 30 Apr 2026 12:53:01 +0100 Subject: [PATCH 08/22] Working towards creation of a surface BLAS from mesh manager data --- include/xdg/cuBQL/ray_tracer.h | 4 ++ src/cuBQL/ray_tracer.cpp | 75 +++++++++++++++++++++++++++++++--- 2 files changed, 74 insertions(+), 5 deletions(-) diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h index eae80009..b2950f51 100644 --- a/include/xdg/cuBQL/ray_tracer.h +++ b/include/xdg/cuBQL/ray_tracer.h @@ -2,6 +2,7 @@ #define _XDG_CUBQL_RAY_TRACING_INTERFACE_H #include +#include #include #include @@ -61,6 +62,9 @@ class CuBQLRayTracer : public RayTracer { const Position& origin, const Direction& direction, double& dist) const override; + +private: + }; } // namespace xdg diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index fa9c3a02..91f7e281 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -6,6 +6,7 @@ #undef __CUDA_ARCH__ #endif +#include #include "cuBQL/bvh.h" #include "cuBQL/builder/omp.h" @@ -29,13 +30,77 @@ CuBQLRayTracer::register_volume(const std::shared_ptr& mesh_manager } TreeID -CuBQLRayTracer::create_surface_tree(const std::shared_ptr&, - MeshID) +CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_manager, + MeshID volume_id) { - // Testing cubql compilation - cuBQL::bvh_t bvh; - + int gpu_id = 0; // TODO - how to manage GPU IDs in a multi-GPU system? + + cuBQL::box3d *d_boxes = nullptr; // box3d is an alias for box_t + int num_boxes = 0; + + cuBQL::bvh3d bvh; // bvh3d is an alias for BinaryBVH + // bvh_t is an alias for BinaryBVH, so bvh_t is also BinaryBVH + + cuBQL::BuildConfig buildParams; // It looks like spatial median is the only supported method for omp builder + // Looks like the cuBQL::gpuBuilder is only pulled in when cubql is built with CUDA support + + cuBQL::build_omp_target(bvh, d_boxes, num_boxes, buildParams, gpu_id); // build bvh on GPU 0 using OpenMP target offloading + + SurfaceTreeID tree = next_surface_tree_id(); + surface_trees_.push_back(tree); + auto volume_surfaces = mesh_manager->get_volume_surfaces(volume_id); + + for (const auto &surf : volume_surfaces) { + auto num_faces = mesh_manager->num_surface_faces(surf); + auto vertices = mesh_manager->get_surface_vertices(surf); + auto indices = mesh_manager->get_surface_connectivity(surf); + + // Get host side storage for vertices using cubql friendly types + std::vector h_vertices; + h_vertices.reserve(vertices.size()); + for (const auto& vertex : vertices) { + h_vertices.emplace_back(vertex.x, vertex.y, vertex.z); + } + + // Get host side storage for indices using cubql friendly types + std::vector h_indices; + h_indices.reserve(indices.size() / 3); + for (size_t i = 0; i < indices.size(); i += 3) { + h_indices.emplace_back(indices[i], indices[i + 1], indices[i + 2]); + } + + // copy vertices and indices to device + auto* d_vertices = static_cast + (omp_target_alloc(h_vertices.size() * sizeof(cuBQL::vec3d), gpu_id)); + + auto* d_indices = static_cast + (omp_target_alloc(h_indices.size() * sizeof(cuBQL::vec3i), gpu_id)); + + // Create device storage for triangle AABBs to be computed in parallel on GPU + auto* d_aabbs = static_cast + (omp_target_alloc(h_indices.size() * sizeof(cuBQL::box3d), gpu_id)); + + // Create AABBs for each triangle in parallel on GPU + // TODO - Should this be its own function? + // TODO - How can this be extended for tets and other element types? + #pragma omp target device(gpu_id) is_device_ptr(d_vertices, d_indices, d_aabbs) + #pragma omp teams distribute parallel for + for (uint32_t primID = 0; primID < num_faces; ++primID) { + cuBQL::vec3i indices = d_indices[primID]; + + cuBQL::vec3d A = d_vertices[indices.x]; + cuBQL::vec3d B = d_vertices[indices.y]; + cuBQL::vec3d C = d_vertices[indices.z]; + + cuBQL::box3d aabb; + aabb.extend(A); + aabb.extend(B); + aabb.extend(C); + + d_aabbs[primID] = aabb; + } + } return TREE_NONE; } From 41cd31ced99cc97b3ee6f4b4bbf4cfaeeffc4ceb Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Thu, 30 Apr 2026 15:08:23 +0100 Subject: [PATCH 09/22] Copy vertices+indices to device and build bvh --- src/cuBQL/ray_tracer.cpp | 28 ++++++++++++++++++++-------- 1 file changed, 20 insertions(+), 8 deletions(-) diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index 91f7e281..aa8e724d 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -35,23 +35,27 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man { int gpu_id = 0; // TODO - how to manage GPU IDs in a multi-GPU system? + int host_id = omp_get_initial_device(); - cuBQL::box3d *d_boxes = nullptr; // box3d is an alias for box_t - int num_boxes = 0; - cuBQL::bvh3d bvh; // bvh3d is an alias for BinaryBVH - // bvh_t is an alias for BinaryBVH, so bvh_t is also BinaryBVH + // cuBQL::box3d *d_boxes = nullptr; // box3d is an alias for box_t + // int num_boxes = 0; + + // cuBQL::bvh3d bvh; // bvh3d is an alias for BinaryBVH + // // bvh_t is an alias for BinaryBVH, so bvh_t is also BinaryBVH cuBQL::BuildConfig buildParams; // It looks like spatial median is the only supported method for omp builder - // Looks like the cuBQL::gpuBuilder is only pulled in when cubql is built with CUDA support + // // Looks like the cuBQL::gpuBuilder is only pulled in when cubql is built with CUDA support - cuBQL::build_omp_target(bvh, d_boxes, num_boxes, buildParams, gpu_id); // build bvh on GPU 0 using OpenMP target offloading + // cuBQL::build_omp_target(bvh, d_boxes, num_boxes, buildParams, gpu_id); // build bvh on GPU 0 using OpenMP target offloading SurfaceTreeID tree = next_surface_tree_id(); surface_trees_.push_back(tree); auto volume_surfaces = mesh_manager->get_volume_surfaces(volume_id); for (const auto &surf : volume_surfaces) { + cuBQL::bvh3d bvh; // bvh3d is an alias for BinaryBVH + auto num_faces = mesh_manager->num_surface_faces(surf); auto vertices = mesh_manager->get_surface_vertices(surf); auto indices = mesh_manager->get_surface_connectivity(surf); @@ -73,9 +77,13 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man // copy vertices and indices to device auto* d_vertices = static_cast (omp_target_alloc(h_vertices.size() * sizeof(cuBQL::vec3d), gpu_id)); + omp_target_memcpy(d_vertices, h_vertices.data(), + h_vertices.size() * sizeof(cuBQL::vec3d), 0, 0, gpu_id, host_id); auto* d_indices = static_cast (omp_target_alloc(h_indices.size() * sizeof(cuBQL::vec3i), gpu_id)); + omp_target_memcpy(d_indices, h_indices.data(), + h_indices.size() * sizeof(cuBQL::vec3i), 0, 0, gpu_id, host_id); // Create device storage for triangle AABBs to be computed in parallel on GPU auto* d_aabbs = static_cast @@ -86,7 +94,7 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man // TODO - How can this be extended for tets and other element types? #pragma omp target device(gpu_id) is_device_ptr(d_vertices, d_indices, d_aabbs) #pragma omp teams distribute parallel for - for (uint32_t primID = 0; primID < num_faces; ++primID) { + for (uint32_t primID = 0; primID < h_indices.size(); ++primID) { cuBQL::vec3i indices = d_indices[primID]; cuBQL::vec3d A = d_vertices[indices.x]; @@ -99,7 +107,11 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man aabb.extend(C); d_aabbs[primID] = aabb; - } + } // This is our AABB population kernel written with OpenMP target offloading pragmas to run in parallel on the GPU + + // Construct the bvh on the gpu using the AABBs with openmp pathway + cuBQL::build_omp_target(bvh, d_aabbs, num_faces, buildParams, gpu_id); + } return TREE_NONE; } From 4031f75ae4036ff016f8339eec6754133898631a Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 1 May 2026 15:32:08 +0100 Subject: [PATCH 10/22] Added objects for storage of bvhs --- include/xdg/cuBQL/ray_tracer.h | 10 ++++++++++ src/cuBQL/ray_tracer.cpp | 22 ++++++++++------------ 2 files changed, 20 insertions(+), 12 deletions(-) diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h index b2950f51..a44c57b6 100644 --- a/include/xdg/cuBQL/ray_tracer.h +++ b/include/xdg/cuBQL/ray_tracer.h @@ -12,6 +12,13 @@ #include "xdg/ray.h" #include "xdg/ray_tracing_interface.h" +// Guards to prevent CUDA headers from being included in host code, which causes failed compilation with LLVM-clang +#if defined(__CUDA_ARCH__) && !defined(__CUDACC__) +#undef __CUDA_ARCH__ +#endif + +#include "cuBQL/bvh.h" + namespace xdg { @@ -64,6 +71,9 @@ class CuBQLRayTracer : public RayTracer { double& dist) const override; private: + std::vector surface_bvhs_; // BVH for each surface tree +// std::unordered_map surface_volume_tree_to_bvh_map; // Map from surface tree IDs to their corresponding cuBQL BVH structures + }; diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index aa8e724d..dc71cde8 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -1,15 +1,9 @@ #include "xdg/cuBQL/ray_tracer.h" #include "xdg/error.h" -// Guards to prevent CUDA headers from being included in host code, which causes failed compilation with LLVM-clang -#if defined(__CUDA_ARCH__) && !defined(__CUDACC__) -#undef __CUDA_ARCH__ -#endif - #include -#include "cuBQL/bvh.h" #include "cuBQL/builder/omp.h" - +#include "cuBQL/traversal/shrinkingRadiusQuery.h" namespace xdg { CuBQLRayTracer::CuBQLRayTracer() = default; @@ -111,7 +105,7 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man // Construct the bvh on the gpu using the AABBs with openmp pathway cuBQL::build_omp_target(bvh, d_aabbs, num_faces, buildParams, gpu_id); - + surface_bvhs_.push_back(bvh); // Store the BVH for this surface tree } return TREE_NONE; } @@ -157,12 +151,16 @@ bool CuBQLRayTracer::point_in_volume(TreeID, std::pair CuBQLRayTracer::ray_fire(TreeID, - const Position&, - const Direction&, - const double, - HitOrientation, + const Position& origin, + const Direction& direction, + const double tmax, + HitOrientation hitOrientation, std::vector* const) { + + // For a closest hit query this is the traversal template we want to make use of + cuBQL::shrinkingRadiusQuery::forEachPrim(prim_lambda, node_lambda, surface_bvhs_[0], tmax) + fatal_error("Ray-fire queries not currently supported with cuBQL ray tracer"); return {INFTY, ID_NONE}; } From a70df01cdcd589a18b40897fc0da8a89d03bae67 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 1 May 2026 16:51:37 +0100 Subject: [PATCH 11/22] Prototyping with codex. Working ray-fire. MAKE SURE TO RESET TO PREVIOUS COMMIT --- include/xdg/cuBQL/ray_tracer.h | 26 ++- src/cuBQL/ray_tracer.cpp | 345 +++++++++++++++++++++++++++++++-- tests/test_ray_fire.cpp | 3 +- tests/util.h | 20 +- tools/ray_fire.cpp | 6 +- 5 files changed, 378 insertions(+), 22 deletions(-) diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h index a44c57b6..b12681f6 100644 --- a/include/xdg/cuBQL/ray_tracer.h +++ b/include/xdg/cuBQL/ray_tracer.h @@ -71,10 +71,28 @@ class CuBQLRayTracer : public RayTracer { double& dist) const override; private: - std::vector surface_bvhs_; // BVH for each surface tree -// std::unordered_map surface_volume_tree_to_bvh_map; // Map from surface tree IDs to their corresponding cuBQL BVH structures - - + struct CuBQLSurfaceBVH { + MeshID surface {ID_NONE}; + cuBQL::bvh3d bvh; + + cuBQL::vec3d* d_vertices {nullptr}; + cuBQL::vec3i* d_indices {nullptr}; + cuBQL::vec3d* d_normals {nullptr}; + MeshID* d_primitive_refs {nullptr}; + + uint32_t num_vertices {0}; + uint32_t num_faces {0}; + int gpu_id {0}; + }; + + struct CuBQLRayHit { + double distance {INFTY}; + MeshID surface {ID_NONE}; + MeshID primitive {ID_NONE}; + }; + + std::vector surface_bvhs_; + std::unordered_map> tree_to_surface_bvh_indices_; }; } // namespace xdg diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index dc71cde8..9da83d6a 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -3,12 +3,154 @@ #include #include "cuBQL/builder/omp.h" -#include "cuBQL/traversal/shrinkingRadiusQuery.h" +#include "cuBQL/math/Ray.h" +#include "cuBQL/queries/triangleData/Triangle.h" +#include "cuBQL/queries/triangleData/math/rayTriangleIntersections.h" +#include "cuBQL/traversal/rayQueries.h" + namespace xdg { +struct CuBQLPluckerIntersectionResult { + bool hit {false}; + double t {0.0}; +}; + +static inline __cubql_both double cubql_ray_hit_tolerance(double t) +{ + constexpr double tolerance = 64.0 * 2.2204460492503131e-16; + return tolerance * (1.0 + cuBQL::abst(t)); +} + +static inline __cubql_both bool cubql_plucker_first(cuBQL::vec3d a, + cuBQL::vec3d b) +{ + if (a.x < b.x) return true; + if (a.x > b.x) return false; + + if (a.y < b.y) return true; + if (a.y > b.y) return false; + + return a.z < b.z; +} + +static inline __cubql_both double cubql_plucker_edge_test(cuBQL::vec3d vertex_a, + cuBQL::vec3d vertex_b, + cuBQL::vec3d ray, + cuBQL::vec3d ray_normal) +{ + double pip; + if (cubql_plucker_first(vertex_a, vertex_b)) { + const cuBQL::vec3d edge = vertex_b - vertex_a; + const cuBQL::vec3d edge_normal = cross(edge, vertex_a); + pip = dot(ray, edge_normal) + dot(ray_normal, edge); + } else { + const cuBQL::vec3d edge = vertex_a - vertex_b; + const cuBQL::vec3d edge_normal = cross(edge, vertex_b); + pip = dot(ray, edge_normal) + dot(ray_normal, edge); + pip = -pip; + } + + constexpr double dbl_zero_tol = 20.0 * 2.2204460492503131e-16; + if (cuBQL::abst(pip) < dbl_zero_tol) { + pip = 0.0; + } + return pip; +} + +static inline __cubql_both CuBQLPluckerIntersectionResult +cubql_plucker_ray_tri_intersect(cuBQL::vec3d vertices[3], + cuBQL::vec3d origin, + cuBQL::vec3d direction, + double t_max, + double t_min) +{ + const cuBQL::vec3d ray_a = direction; + const cuBQL::vec3d ray_b = cross(direction, origin); + + const double plucker_coord_0 = + cubql_plucker_edge_test(vertices[0], vertices[1], ray_a, ray_b); + + const double plucker_coord_1 = + cubql_plucker_edge_test(vertices[1], vertices[2], ray_a, ray_b); + + if ((0.0 < plucker_coord_0 && 0.0 > plucker_coord_1) || + (0.0 > plucker_coord_0 && 0.0 < plucker_coord_1)) { + return {}; + } + + const double plucker_coord_2 = + cubql_plucker_edge_test(vertices[2], vertices[0], ray_a, ray_b); + + if ((0.0 < plucker_coord_1 && 0.0 > plucker_coord_2) || + (0.0 > plucker_coord_1 && 0.0 < plucker_coord_2) || + (0.0 < plucker_coord_0 && 0.0 > plucker_coord_2) || + (0.0 > plucker_coord_0 && 0.0 < plucker_coord_2)) { + return {}; + } + + if (plucker_coord_0 == 0.0 && + plucker_coord_1 == 0.0 && + plucker_coord_2 == 0.0) { + return {}; + } + + const double inverse_sum = + 1.0 / (plucker_coord_0 + plucker_coord_1 + plucker_coord_2); + + const cuBQL::vec3d intersection = + plucker_coord_0 * inverse_sum * vertices[2] + + plucker_coord_1 * inverse_sum * vertices[0] + + plucker_coord_2 * inverse_sum * vertices[1]; + + int idx = 0; + double max_abs_dir = 0.0; + for (int i = 0; i < 3; ++i) { + if (cuBQL::abst(direction[i]) > max_abs_dir) { + idx = i; + max_abs_dir = cuBQL::abst(direction[i]); + } + } + + double dist_out = (intersection[idx] - origin[idx]) / direction[idx]; + + const double u = plucker_coord_2 * inverse_sum; + const double v = plucker_coord_0 * inverse_sum; + if (u < 0.0 || v < 0.0 || (u + v) > 1.0) { + return {}; + } + + if (dist_out < t_min || dist_out > t_max) { + return {}; + } + + return {true, dist_out}; +} + CuBQLRayTracer::CuBQLRayTracer() = default; -CuBQLRayTracer::~CuBQLRayTracer() = default; +CuBQLRayTracer::~CuBQLRayTracer() +{ + for (auto& surface_bvh : surface_bvhs_) { + cuBQL::omp::Context context(surface_bvh.gpu_id); + + if (surface_bvh.bvh.nodes || surface_bvh.bvh.primIDs) { + cuBQL::omp::freeBVH(surface_bvh.bvh, &context); + } + + if (surface_bvh.d_vertices) { + omp_target_free(surface_bvh.d_vertices, surface_bvh.gpu_id); + } + if (surface_bvh.d_indices) { + omp_target_free(surface_bvh.d_indices, surface_bvh.gpu_id); + } + if (surface_bvh.d_normals) { + omp_target_free(surface_bvh.d_normals, surface_bvh.gpu_id); + } + if (surface_bvh.d_primitive_refs) { + omp_target_free(surface_bvh.d_primitive_refs, surface_bvh.gpu_id); + } + } +} void CuBQLRayTracer::init() { @@ -45,11 +187,10 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man SurfaceTreeID tree = next_surface_tree_id(); surface_trees_.push_back(tree); + auto& surface_bvh_indices = tree_to_surface_bvh_indices_[tree]; auto volume_surfaces = mesh_manager->get_volume_surfaces(volume_id); for (const auto &surf : volume_surfaces) { - cuBQL::bvh3d bvh; // bvh3d is an alias for BinaryBVH - auto num_faces = mesh_manager->num_surface_faces(surf); auto vertices = mesh_manager->get_surface_vertices(surf); auto indices = mesh_manager->get_surface_connectivity(surf); @@ -68,6 +209,26 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man h_indices.emplace_back(indices[i], indices[i + 1], indices[i + 2]); } + // Get storage for primitive refs so a hit can be mapped back to a mesh face. + std::vector h_primitive_refs; + h_primitive_refs.reserve(num_faces); + + std::vector h_normals; + h_normals.reserve(num_faces); + + auto [forward_parent, reverse_parent] = mesh_manager->get_parent_volumes(surf); + for (const auto& face : mesh_manager->get_surface_faces(surf)) { + h_primitive_refs.push_back(face); + + auto normal = mesh_manager->face_normal(face); + if (volume_id == reverse_parent) { + normal = -normal; + } else if (volume_id != forward_parent) { + fatal_error("Volume {} is not a parent of surface {}", volume_id, surf); + } + h_normals.emplace_back(normal.x, normal.y, normal.z); + } + // copy vertices and indices to device auto* d_vertices = static_cast (omp_target_alloc(h_vertices.size() * sizeof(cuBQL::vec3d), gpu_id)); @@ -79,16 +240,27 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man omp_target_memcpy(d_indices, h_indices.data(), h_indices.size() * sizeof(cuBQL::vec3i), 0, 0, gpu_id, host_id); + auto* d_normals = static_cast + (omp_target_alloc(h_normals.size() * sizeof(cuBQL::vec3d), gpu_id)); + omp_target_memcpy(d_normals, h_normals.data(), + h_normals.size() * sizeof(cuBQL::vec3d), 0, 0, gpu_id, host_id); + + auto* d_primitive_refs = static_cast + (omp_target_alloc(h_primitive_refs.size() * sizeof(MeshID), gpu_id)); + omp_target_memcpy(d_primitive_refs, h_primitive_refs.data(), + h_primitive_refs.size() * sizeof(MeshID), 0, 0, gpu_id, host_id); + // Create device storage for triangle AABBs to be computed in parallel on GPU auto* d_aabbs = static_cast (omp_target_alloc(h_indices.size() * sizeof(cuBQL::box3d), gpu_id)); + const auto num_primitives = static_cast(h_indices.size()); // Create AABBs for each triangle in parallel on GPU // TODO - Should this be its own function? // TODO - How can this be extended for tets and other element types? #pragma omp target device(gpu_id) is_device_ptr(d_vertices, d_indices, d_aabbs) #pragma omp teams distribute parallel for - for (uint32_t primID = 0; primID < h_indices.size(); ++primID) { + for (uint32_t primID = 0; primID < num_primitives; ++primID) { cuBQL::vec3i indices = d_indices[primID]; cuBQL::vec3d A = d_vertices[indices.x]; @@ -104,17 +276,33 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man } // This is our AABB population kernel written with OpenMP target offloading pragmas to run in parallel on the GPU // Construct the bvh on the gpu using the AABBs with openmp pathway + cuBQL::bvh3d bvh; // bvh3d is an alias for BinaryBVH cuBQL::build_omp_target(bvh, d_aabbs, num_faces, buildParams, gpu_id); - surface_bvhs_.push_back(bvh); // Store the BVH for this surface tree + + omp_target_free(d_aabbs, gpu_id); + + CuBQLSurfaceBVH surface_bvh; + surface_bvh.surface = surf; + surface_bvh.bvh = bvh; + surface_bvh.d_vertices = d_vertices; + surface_bvh.d_indices = d_indices; + surface_bvh.d_normals = d_normals; + surface_bvh.d_primitive_refs = d_primitive_refs; + surface_bvh.num_vertices = h_vertices.size(); + surface_bvh.num_faces = num_faces; + surface_bvh.gpu_id = gpu_id; + + surface_bvh_indices.push_back(surface_bvhs_.size()); + surface_bvhs_.push_back(surface_bvh); } - return TREE_NONE; + return tree; } TreeID CuBQLRayTracer::create_element_tree(const std::shared_ptr&, MeshID) { - fatal_error("Element trees not currently supported with cuBQL ray tracer"); + warning("Element trees not currently supported with cuBQL ray tracer"); return TREE_NONE; } @@ -150,19 +338,148 @@ bool CuBQLRayTracer::point_in_volume(TreeID, } std::pair -CuBQLRayTracer::ray_fire(TreeID, +CuBQLRayTracer::ray_fire(TreeID tree, const Position& origin, const Direction& direction, const double tmax, HitOrientation hitOrientation, - std::vector* const) + std::vector* const exclude_primitives) { + int gpu_id = 0; // TODO - how to manage GPU IDs in a multi-GPU system? + int host_id = omp_get_initial_device(); - // For a closest hit query this is the traversal template we want to make use of - cuBQL::shrinkingRadiusQuery::forEachPrim(prim_lambda, node_lambda, surface_bvhs_[0], tmax) + MeshID* d_exclude_primitives = nullptr; + int exclude_count = 0; + if (exclude_primitives && !exclude_primitives->empty()) { + exclude_count = static_cast(exclude_primitives->size()); + d_exclude_primitives = static_cast + (omp_target_alloc(exclude_count * sizeof(MeshID), gpu_id)); + omp_target_memcpy(d_exclude_primitives, + exclude_primitives->data(), + exclude_count * sizeof(MeshID), + 0, + 0, + gpu_id, + host_id); + } - fatal_error("Ray-fire queries not currently supported with cuBQL ray tracer"); - return {INFTY, ID_NONE}; + double closest_distance = tmax; + MeshID closest_surface = ID_NONE; + MeshID closest_primitive = ID_NONE; + auto* d_surface_hit = static_cast + (omp_target_alloc(sizeof(CuBQLRayHit), gpu_id)); + + const auto& surface_bvh_indices = tree_to_surface_bvh_indices_.at(tree); + for (const auto surface_bvh_index : surface_bvh_indices) { + const auto& surface_bvh = surface_bvhs_.at(surface_bvh_index); + + cuBQL::bvh3d bvh = surface_bvh.bvh; + const cuBQL::vec3d* d_vertices = surface_bvh.d_vertices; + const cuBQL::vec3i* d_indices = surface_bvh.d_indices; + const cuBQL::vec3d* d_normals = surface_bvh.d_normals; + const MeshID* d_primitive_refs = surface_bvh.d_primitive_refs; + + CuBQLRayHit surface_hit; + surface_hit.distance = closest_distance + cubql_ray_hit_tolerance(closest_distance); + omp_target_memcpy(d_surface_hit, + &surface_hit, + sizeof(CuBQLRayHit), + 0, + 0, + gpu_id, + host_id); + + const int orientation = static_cast(hitOrientation); + const cuBQL::vec3d ray_origin(origin.x, origin.y, origin.z); + const cuBQL::vec3d ray_direction(direction.x, direction.y, direction.z); + + #pragma omp target device(gpu_id) \ + is_device_ptr(d_vertices, d_indices, d_normals, d_primitive_refs, d_exclude_primitives, d_surface_hit) + { + cuBQL::ray3d ray(ray_origin, ray_direction, 0.0, d_surface_hit->distance); + + auto intersect_prim = [=, &ray] + (uint32_t prim_id) -> double + { + const MeshID primitive_ref = d_primitive_refs[prim_id]; + + for (int i = 0; i < exclude_count; ++i) { + if (d_exclude_primitives[i] == primitive_ref) { + return ray.tMax; + } + } + + const cuBQL::vec3i index = d_indices[prim_id]; + cuBQL::vec3d vertices[3] = { + d_vertices[index.x], + d_vertices[index.y], + d_vertices[index.z] + }; + + const double normal_dot_direction = dot(d_normals[prim_id], ray.direction); + bool culled = false; + if (orientation == static_cast(HitOrientation::EXITING)) { + culled = normal_dot_direction < 0.0; + } else if (orientation == static_cast(HitOrientation::ENTERING)) { + culled = normal_dot_direction >= 0.0; + } + if (culled) { + return ray.tMax; + } + + auto intersection = cubql_plucker_ray_tri_intersect(vertices, + ray.origin, + ray.direction, + ray.tMax, + ray.tMin); + if (intersection.hit) { + d_surface_hit->distance = intersection.t; + d_surface_hit->primitive = primitive_ref; + ray.tMax = intersection.t; + } + + return ray.tMax; + }; + + cuBQL::shrinkingRayQuery::forEachPrim(intersect_prim, bvh, ray); + } + + omp_target_memcpy(&surface_hit, + d_surface_hit, + sizeof(CuBQLRayHit), + 0, + 0, + host_id, + gpu_id); + + const double hit_tolerance = cubql_ray_hit_tolerance(closest_distance); + const bool closer_hit = surface_hit.distance < closest_distance - hit_tolerance; + const bool tied_hit = cuBQL::abst(surface_hit.distance - closest_distance) <= hit_tolerance; + const bool lower_surface_tie = closest_surface == ID_NONE || + surface_bvh.surface < closest_surface; + + if (surface_hit.primitive != ID_NONE && (closer_hit || (tied_hit && lower_surface_tie))) { + closest_distance = surface_hit.distance; + closest_surface = surface_bvh.surface; + closest_primitive = surface_hit.primitive; + } + } + + omp_target_free(d_surface_hit, gpu_id); + + if (d_exclude_primitives) { + omp_target_free(d_exclude_primitives, gpu_id); + } + + if (closest_surface == ID_NONE) { + return {INFTY, ID_NONE}; + } + + if (exclude_primitives) { + exclude_primitives->push_back(closest_primitive); + } + + return {closest_distance, closest_surface}; } std::pair diff --git a/tests/test_ray_fire.cpp b/tests/test_ray_fire.cpp index be816c36..70a5353f 100644 --- a/tests/test_ray_fire.cpp +++ b/tests/test_ray_fire.cpp @@ -18,7 +18,8 @@ using namespace xdg::test; TEMPLATE_TEST_CASE("Ray Fire on MeshMock (per-backend sections)", "[rayfire][mock]", Embree_Raytracer, - GPRT_Raytracer) + GPRT_Raytracer, + CuBQL_Raytracer) { // Generate one test run per enabled backend constexpr auto rt_backend = TestType::value; diff --git a/tests/util.h b/tests/util.h index a9841086..86657d9d 100644 --- a/tests/util.h +++ b/tests/util.h @@ -17,7 +17,7 @@ using LibMesh_Interface = std::integral_constant; using GPRT_Raytracer = std::integral_constant; - +using CuBQL_Raytracer = std::integral_constant; } // namespace xdg::test namespace Catch { @@ -49,6 +49,13 @@ inline bool ray_tracer_available(xdg::RTLibrary rt) { #else return false; #endif + + case xdg::RTLibrary::CUBQL: + #ifdef XDG_ENABLE_CUBQL + return true; + #else + return false; + #endif } return false; @@ -81,6 +88,12 @@ inline void check_ray_tracer_supported(xdg::RTLibrary rt) { if (!ray_tracer_available(rt)) { SKIP(fmt::format("{} backend unavailable; skipping.", xdg::RT_LIB_TO_STR.at(rt))); } + #endif + #ifndef XDG_ENABLE_CUBQL + if(!ray_tracer_available(xdg::RTLibrary::CUBQL)) { + SKIP(fmt::format("{} backend unavailable; skipping.", xdg::RT_LIB_TO_STR.at(rt))); + } + #endif } inline void check_mesh_library_supported(xdg::MeshLibrary mesh) { @@ -118,5 +131,10 @@ create_raytracer(xdg::RTLibrary rt) { return std::make_shared(); #endif + #ifdef XDG_ENABLE_CUBQL + if (rt == xdg::RTLibrary::CUBQL) + return std::make_shared(); + #endif + return nullptr; } diff --git a/tools/ray_fire.cpp b/tools/ray_fire.cpp index c75bbe82..d827b335 100644 --- a/tools/ray_fire.cpp +++ b/tools/ray_fire.cpp @@ -40,7 +40,7 @@ int main(int argc, char** argv) { .default_value("MOAB"); args.add_argument("-r", "--rt-library") - .help("Ray tracing library to use. One of (EMBREE, GPRT)") + .help("Ray tracing library to use. One of (EMBREE, GPRT, CUBQL)") .default_value("EMBREE"); try { @@ -60,6 +60,8 @@ if (rt_str == "EMBREE") rt_lib = RTLibrary::EMBREE; else if (rt_str == "GPRT") rt_lib = RTLibrary::GPRT; +else if (rt_str == "CUBQL") + rt_lib = RTLibrary::CUBQL; else fatal_error("Invalid ray tracing library '{}' specified", rt_str); @@ -68,7 +70,7 @@ if (mesh_str == "MOAB") mesh_lib = MeshLibrary::MOAB; else if (mesh_str == "LIBMESH") { mesh_lib = MeshLibrary::LIBMESH; - if (rt_lib == RTLibrary::GPRT) + if (rt_lib == RTLibrary::GPRT || rt_lib == RTLibrary::CUBQL) fatal_error("LibMesh is not currently supported with GPRT"); } else From 1766d71ca62216d4e84796aed99d906438f04e76 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Wed, 6 May 2026 14:42:19 +0100 Subject: [PATCH 12/22] Added some additional comments and cleanups --- include/xdg/cuBQL/ray_tracer.h | 2 +- src/cuBQL/ray_tracer.cpp | 37 ++++++++++++++++++++++++---------- 2 files changed, 27 insertions(+), 12 deletions(-) diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h index b12681f6..f0a59b04 100644 --- a/include/xdg/cuBQL/ray_tracer.h +++ b/include/xdg/cuBQL/ray_tracer.h @@ -85,7 +85,7 @@ class CuBQLRayTracer : public RayTracer { int gpu_id {0}; }; - struct CuBQLRayHit { + struct CubqlHit { double distance {INFTY}; MeshID surface {ID_NONE}; MeshID primitive {ID_NONE}; diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index 9da83d6a..bc9680ac 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -220,6 +220,9 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man for (const auto& face : mesh_manager->get_surface_faces(surf)) { h_primitive_refs.push_back(face); + // TODO: This pre-orients normals for the current volume tree. If cuBQL + // starts reusing a single surface BVH/BLAS across both parent volumes, + // normals must be flipped at intersection time based on the queried volume. auto normal = mesh_manager->face_normal(face); if (volume_id == reverse_parent) { normal = -normal; @@ -366,8 +369,8 @@ CuBQLRayTracer::ray_fire(TreeID tree, double closest_distance = tmax; MeshID closest_surface = ID_NONE; MeshID closest_primitive = ID_NONE; - auto* d_surface_hit = static_cast - (omp_target_alloc(sizeof(CuBQLRayHit), gpu_id)); + auto* d_surface_hit = static_cast + (omp_target_alloc(sizeof(CubqlHit), gpu_id)); const auto& surface_bvh_indices = tree_to_surface_bvh_indices_.at(tree); for (const auto surface_bvh_index : surface_bvh_indices) { @@ -379,11 +382,11 @@ CuBQLRayTracer::ray_fire(TreeID tree, const cuBQL::vec3d* d_normals = surface_bvh.d_normals; const MeshID* d_primitive_refs = surface_bvh.d_primitive_refs; - CuBQLRayHit surface_hit; + CubqlHit surface_hit; surface_hit.distance = closest_distance + cubql_ray_hit_tolerance(closest_distance); omp_target_memcpy(d_surface_hit, &surface_hit, - sizeof(CuBQLRayHit), + sizeof(CubqlHit), 0, 0, gpu_id, @@ -396,9 +399,17 @@ CuBQLRayTracer::ray_fire(TreeID tree, #pragma omp target device(gpu_id) \ is_device_ptr(d_vertices, d_indices, d_normals, d_primitive_refs, d_exclude_primitives, d_surface_hit) { - cuBQL::ray3d ray(ray_origin, ray_direction, 0.0, d_surface_hit->distance); - - auto intersect_prim = [=, &ray] + cuBQL::ray3d ray; + ray.origin = ray_origin; + ray.direction = ray_direction; + ray.tMin = 0.0; + ray.tMax = d_surface_hit->distance; + + // lambda for primitive intersection + // - culls excluded primitives + // - culls backfacing or frontfacing hits based on hitOrientation + // - calls plucker intersection + auto xdg_plucker_intersect_prim = [=, &ray] (uint32_t prim_id) -> double { const MeshID primitive_ref = d_primitive_refs[prim_id]; @@ -441,23 +452,27 @@ CuBQLRayTracer::ray_fire(TreeID tree, return ray.tMax; }; - cuBQL::shrinkingRayQuery::forEachPrim(intersect_prim, bvh, ray); + // Traversal template from cuBQL + cuBQL::shrinkingRayQuery::forEachPrim(xdg_plucker_intersect_prim, bvh, ray); } omp_target_memcpy(&surface_hit, d_surface_hit, - sizeof(CuBQLRayHit), + sizeof(CubqlHit), 0, 0, host_id, gpu_id); - + + // Temporary cuBQL-only closest-hit reduction while each surface BVH is queried + // independently. Embree/GPRT delegate this to a single scene/TLAS traversal; + // once cuBQL does the same, this tolerance/tie-break logic should go away. const double hit_tolerance = cubql_ray_hit_tolerance(closest_distance); const bool closer_hit = surface_hit.distance < closest_distance - hit_tolerance; const bool tied_hit = cuBQL::abst(surface_hit.distance - closest_distance) <= hit_tolerance; const bool lower_surface_tie = closest_surface == ID_NONE || surface_bvh.surface < closest_surface; - + if (surface_hit.primitive != ID_NONE && (closer_hit || (tied_hit && lower_surface_tie))) { closest_distance = surface_hit.distance; closest_surface = surface_bvh.surface; From 69895a095869f24e9bfb768c0e11d73ddd2f39a3 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Thu, 7 May 2026 13:50:01 +0100 Subject: [PATCH 13/22] Cleaned up the plucker intersection function for cuBQL --- CMakeLists.txt | 8 +++ include/xdg/geometry/dp_math.h | 22 +++++- src/cuBQL/ray_tracer.cpp | 124 +++------------------------------ 3 files changed, 38 insertions(+), 116 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 00ab66e9..ab5b282b 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -242,6 +242,14 @@ if (XDG_ENABLE_CUBQL) list(APPEND xdg_sources src/cuBQL/ray_tracer.cpp ) + +# We need a precompile definition to switch to using the cuBQL math types in the shared +# plucker intersection code. The compile definition is used in dp__math.h +set_source_files_properties( + src/cuBQL/ray_tracer.cpp + PROPERTIES COMPILE_DEFINITIONS XDG_DP_MATH_CUBQL +) + endif() if (XDG_ENABLE_LIBMESH) diff --git a/include/xdg/geometry/dp_math.h b/include/xdg/geometry/dp_math.h index 5b162d3f..5dfedffd 100644 --- a/include/xdg/geometry/dp_math.h +++ b/include/xdg/geometry/dp_math.h @@ -22,7 +22,27 @@ namespace dp { static const double INFTY = 1.7976931348623157e+308; // std::numeric_limits::max() is not available in slang } -#else +// TODO - Is this the right way to handle this for cubql openmp target offload compilation? +// In theory, we can compile the C++ pathway that Embree uses but will the vec3da types be +// omptarget friendly? +// For now I have defined this separate compilation pathway which is enabled with a new +// precompile definition only set when when compiling the CuBQLRayTracer that maps to cuBQL's math types +#elif defined(XDG_DP_MATH_CUBQL) +#include "cuBQL/math/common.h" + +// C++ compilation with cuBQL math types, map dp::vec3 -> cuBQL::vec3d +namespace dp { + typedef cuBQL::vec3d vec3; + + inline double dot(vec3 a, vec3 b) { return cuBQL::dot(a, b); } + inline vec3 cross(vec3 a, vec3 b) { return cuBQL::cross(a, b); } + inline double abs(double a) { return cuBQL::abst(a); } + + static constexpr double DBL_ZERO_TOL = 20.0 * 2.2204460492503131e-16; // same as 20 * std::numeric_limits::epsilon() + constexpr double INFTY {std::numeric_limits::max()}; +} + +#else #include "xdg/vec3da.h" // C++ compilation map dp::vec3 -> xdg::Vec3da diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index bc9680ac..ea0e2216 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -1,5 +1,6 @@ #include "xdg/cuBQL/ray_tracer.h" #include "xdg/error.h" +#include "xdg/geometry/plucker.h" #include #include "cuBQL/builder/omp.h" @@ -10,122 +11,12 @@ namespace xdg { -struct CuBQLPluckerIntersectionResult { - bool hit {false}; - double t {0.0}; -}; - static inline __cubql_both double cubql_ray_hit_tolerance(double t) { constexpr double tolerance = 64.0 * 2.2204460492503131e-16; return tolerance * (1.0 + cuBQL::abst(t)); } -static inline __cubql_both bool cubql_plucker_first(cuBQL::vec3d a, - cuBQL::vec3d b) -{ - if (a.x < b.x) return true; - if (a.x > b.x) return false; - - if (a.y < b.y) return true; - if (a.y > b.y) return false; - - return a.z < b.z; -} - -static inline __cubql_both double cubql_plucker_edge_test(cuBQL::vec3d vertex_a, - cuBQL::vec3d vertex_b, - cuBQL::vec3d ray, - cuBQL::vec3d ray_normal) -{ - double pip; - if (cubql_plucker_first(vertex_a, vertex_b)) { - const cuBQL::vec3d edge = vertex_b - vertex_a; - const cuBQL::vec3d edge_normal = cross(edge, vertex_a); - pip = dot(ray, edge_normal) + dot(ray_normal, edge); - } else { - const cuBQL::vec3d edge = vertex_a - vertex_b; - const cuBQL::vec3d edge_normal = cross(edge, vertex_b); - pip = dot(ray, edge_normal) + dot(ray_normal, edge); - pip = -pip; - } - - constexpr double dbl_zero_tol = 20.0 * 2.2204460492503131e-16; - if (cuBQL::abst(pip) < dbl_zero_tol) { - pip = 0.0; - } - return pip; -} - -static inline __cubql_both CuBQLPluckerIntersectionResult -cubql_plucker_ray_tri_intersect(cuBQL::vec3d vertices[3], - cuBQL::vec3d origin, - cuBQL::vec3d direction, - double t_max, - double t_min) -{ - const cuBQL::vec3d ray_a = direction; - const cuBQL::vec3d ray_b = cross(direction, origin); - - const double plucker_coord_0 = - cubql_plucker_edge_test(vertices[0], vertices[1], ray_a, ray_b); - - const double plucker_coord_1 = - cubql_plucker_edge_test(vertices[1], vertices[2], ray_a, ray_b); - - if ((0.0 < plucker_coord_0 && 0.0 > plucker_coord_1) || - (0.0 > plucker_coord_0 && 0.0 < plucker_coord_1)) { - return {}; - } - - const double plucker_coord_2 = - cubql_plucker_edge_test(vertices[2], vertices[0], ray_a, ray_b); - - if ((0.0 < plucker_coord_1 && 0.0 > plucker_coord_2) || - (0.0 > plucker_coord_1 && 0.0 < plucker_coord_2) || - (0.0 < plucker_coord_0 && 0.0 > plucker_coord_2) || - (0.0 > plucker_coord_0 && 0.0 < plucker_coord_2)) { - return {}; - } - - if (plucker_coord_0 == 0.0 && - plucker_coord_1 == 0.0 && - plucker_coord_2 == 0.0) { - return {}; - } - - const double inverse_sum = - 1.0 / (plucker_coord_0 + plucker_coord_1 + plucker_coord_2); - - const cuBQL::vec3d intersection = - plucker_coord_0 * inverse_sum * vertices[2] + - plucker_coord_1 * inverse_sum * vertices[0] + - plucker_coord_2 * inverse_sum * vertices[1]; - - int idx = 0; - double max_abs_dir = 0.0; - for (int i = 0; i < 3; ++i) { - if (cuBQL::abst(direction[i]) > max_abs_dir) { - idx = i; - max_abs_dir = cuBQL::abst(direction[i]); - } - } - - double dist_out = (intersection[idx] - origin[idx]) / direction[idx]; - - const double u = plucker_coord_2 * inverse_sum; - const double v = plucker_coord_0 * inverse_sum; - if (u < 0.0 || v < 0.0 || (u + v) > 1.0) { - return {}; - } - - if (dist_out < t_min || dist_out > t_max) { - return {}; - } - - return {true, dist_out}; -} - CuBQLRayTracer::CuBQLRayTracer() = default; CuBQLRayTracer::~CuBQLRayTracer() @@ -438,11 +329,14 @@ CuBQLRayTracer::ray_fire(TreeID tree, return ray.tMax; } - auto intersection = cubql_plucker_ray_tri_intersect(vertices, - ray.origin, - ray.direction, - ray.tMax, - ray.tMin); + // Reuse same plucker intersection function applied to other ray tracers + auto intersection = plucker_ray_tri_intersect(vertices, + ray.origin, + ray.direction, + ray.tMax, + ray.tMin, + false, + 0); if (intersection.hit) { d_surface_hit->distance = intersection.t; d_surface_hit->primitive = primitive_ref; From 764086d85c23e28697c3c6239e6c1e1b64fa8f15 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Thu, 7 May 2026 13:54:27 +0100 Subject: [PATCH 14/22] Removed constexpr EARLY_EXIT as nvlink couldn't resolve it --- include/xdg/geometry/plucker.h | 16 +++++++--------- 1 file changed, 7 insertions(+), 9 deletions(-) diff --git a/include/xdg/geometry/plucker.h b/include/xdg/geometry/plucker.h index acf0ebfa..8b7581c9 100644 --- a/include/xdg/geometry/plucker.h +++ b/include/xdg/geometry/plucker.h @@ -25,8 +25,6 @@ struct PluckerIntersectionResult { double t = 0.0; // Distance along the ray to the intersection point }; -static constexpr PluckerIntersectionResult EXIT_EARLY = {false, 0.0}; - /* Function to return the vertex with the lowest coordinates. To force the same ray-edge computation, the Plücker test needs to use consistent edge representation. This would be more simple with MOAB handles instead of @@ -81,7 +79,7 @@ inline PluckerIntersectionResult plucker_ray_tri_intersect(dp::vec3 vertices[3], // If orientation is set, confirm that sign of plucker_coordinate indicate // correct orientation of intersection if (useOrientation && orientation * plucker_coord0 > 0) { - return EXIT_EARLY; + return {false, 0.0}; } // Determine the value of the second Plucker coordinate from edge 1 @@ -92,13 +90,13 @@ inline PluckerIntersectionResult plucker_ray_tri_intersect(dp::vec3 vertices[3], // correct orientation of intersection if (useOrientation) { if (orientation * plucker_coord1 > 0) { - return EXIT_EARLY; + return {false, 0.0}; } // If the orientation is not specified, all plucker_coords must be the same // sign or zero. } else if ((0.0 < plucker_coord0 && 0.0 > plucker_coord1) || (0.0 > plucker_coord0 && 0.0 < plucker_coord1)) { - return EXIT_EARLY; + return {false, 0.0}; } // Determine the value of the third Plucker coordinate from edge 2 @@ -109,7 +107,7 @@ inline PluckerIntersectionResult plucker_ray_tri_intersect(dp::vec3 vertices[3], // correct orientation of intersection if (useOrientation) { if (orientation * plucker_coord2 > 0) { - return EXIT_EARLY; + return {false, 0.0}; } // If the orientation is not specified, all plucker_coords must be the same // sign or zero. @@ -117,12 +115,12 @@ inline PluckerIntersectionResult plucker_ray_tri_intersect(dp::vec3 vertices[3], (0.0 > plucker_coord1 && 0.0 < plucker_coord2) || (0.0 < plucker_coord0 && 0.0 > plucker_coord2) || (0.0 > plucker_coord0 && 0.0 < plucker_coord2)) { - return EXIT_EARLY; + return {false, 0.0}; } // check for coplanar case to avoid dividing by zero if (0.0 == plucker_coord0 && 0.0 == plucker_coord1 && 0.0 == plucker_coord2) { - return EXIT_EARLY; + return {false, 0.0}; } // get the distance to intersection @@ -155,7 +153,7 @@ inline PluckerIntersectionResult plucker_ray_tri_intersect(dp::vec3 vertices[3], } // is the intersection within distance limits? - if (dist_out < tMin || dist_out > tMax) return EXIT_EARLY; + if (dist_out < tMin || dist_out > tMax) return {false, 0.0}; return {true, dist_out}; } From 21c9fd2471243f4a7b01b579e7a6b8bb975c2e39 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Thu, 14 May 2026 14:01:19 +0100 Subject: [PATCH 15/22] Fixes for rebase from main --- CMakeLists.txt | 6 +----- vendor/GPRT | 2 +- 2 files changed, 2 insertions(+), 6 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index ab5b282b..13a61fcf 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -355,7 +355,6 @@ if (XDG_ENABLE_GPRT) target_link_options(xdg PRIVATE -Wl,--unresolved-symbols=ignore-in-shared-libs) endif() -target_link_libraries(xdg PRIVATE fmt::fmt) if (XDG_ENABLE_CUBQL) target_compile_definitions(xdg PUBLIC XDG_ENABLE_CUBQL) target_link_libraries(xdg $) @@ -368,10 +367,7 @@ if (XDG_ENABLE_CUBQL) # whether install RPATH should remain environment-module based or be opt-in. endif() -target_link_libraries(xdg fmt::fmt) - -#TODO: refactor to only be required when tools are enabled -target_link_libraries(xdg indicators::indicators) +target_link_libraries(xdg PRIVATE fmt::fmt) # ========================== # Link ray tracing libraries diff --git a/vendor/GPRT b/vendor/GPRT index f1e95e41..7b55a1ae 160000 --- a/vendor/GPRT +++ b/vendor/GPRT @@ -1 +1 @@ -Subproject commit f1e95e4188cde591547d6b4a33a70bf2afaeec59 +Subproject commit 7b55a1ae432b226be18592f6f2e7f8980c8df9ce From 68c51ae5097e6ad4940463a194f0eb9ce5eb68e6 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Thu, 14 May 2026 16:35:19 +0100 Subject: [PATCH 16/22] Updated CMakeLists to make proper use of PRIVATE after rebasing recent changes from main --- CMakeLists.txt | 2 +- tests/CMakeLists.txt | 3 +++ 2 files changed, 4 insertions(+), 1 deletion(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 13a61fcf..d05306a1 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -357,7 +357,7 @@ endif() if (XDG_ENABLE_CUBQL) target_compile_definitions(xdg PUBLIC XDG_ENABLE_CUBQL) - target_link_libraries(xdg $) + target_link_libraries(xdg PRIVATE $) # TODO: Stop relying on LD_LIBRARY_PATH for LLVM OpenMP offload runtimes. # As a temporary measure whilst figuring out our way around cuBQL this is # okay but in the long run we should aim for a more robust solutions here. diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 222626f2..0009dbfb 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -53,6 +53,9 @@ foreach(test ${TEST_NAMES}) if (XDG_ENABLE_MOAB) target_link_libraries(${test} PRIVATE MOAB) endif() + if (XDG_ENABLE_CUBQL) + target_link_libraries(${test} PRIVATE $) + endif() set_target_properties(${test} PROPERTIES BUILD_RPATH "$") catch_discover_tests(${test} From b742aebe08e390fa2a6661a2ffcae37c5abca219 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 15 May 2026 13:42:54 +0100 Subject: [PATCH 17/22] Fix runtime check for available rt libraries --- tests/util.h | 6 ------ 1 file changed, 6 deletions(-) diff --git a/tests/util.h b/tests/util.h index 86657d9d..074ff899 100644 --- a/tests/util.h +++ b/tests/util.h @@ -88,12 +88,6 @@ inline void check_ray_tracer_supported(xdg::RTLibrary rt) { if (!ray_tracer_available(rt)) { SKIP(fmt::format("{} backend unavailable; skipping.", xdg::RT_LIB_TO_STR.at(rt))); } - #endif - #ifndef XDG_ENABLE_CUBQL - if(!ray_tracer_available(xdg::RTLibrary::CUBQL)) { - SKIP(fmt::format("{} backend unavailable; skipping.", xdg::RT_LIB_TO_STR.at(rt))); - } - #endif } inline void check_mesh_library_supported(xdg::MeshLibrary mesh) { From 6cf0adb911b51932d2d488d79ab701bddc2d229b Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 15 May 2026 14:21:24 +0100 Subject: [PATCH 18/22] Reducing the number of queries in cross check for now since cubql pathway a little slow --- tests/test_ray_tracer_cross_check.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_ray_tracer_cross_check.cpp b/tests/test_ray_tracer_cross_check.cpp index c75b0273..3cde9f73 100644 --- a/tests/test_ray_tracer_cross_check.cpp +++ b/tests/test_ray_tracer_cross_check.cpp @@ -57,7 +57,7 @@ TEST_CASE("Test Pincell RT libraries Cross-Check ray_fire queries", "[moab][rayf } srand48(12345); // set fixed seed for rng - std::vector directions(1000); + std::vector directions(100); for (auto &dir : directions) { dir = rand_dir(); } From f486c3b351338124c832ec20c85d4e32a7ae052e Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 15 May 2026 14:24:03 +0100 Subject: [PATCH 19/22] Reworked normal/sense handling so normals are calculated and sense handling is performed at intersection time Also had to change the create_global_*_trees stubs for cuBQL to return a warning rather than fatal_error Cross check verified to produce consistent results across all three ray tracers --- include/xdg/cuBQL/ray_tracer.h | 4 ++- src/cuBQL/ray_tracer.cpp | 49 ++++++++++++++-------------------- 2 files changed, 23 insertions(+), 30 deletions(-) diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h index f0a59b04..528e16ce 100644 --- a/include/xdg/cuBQL/ray_tracer.h +++ b/include/xdg/cuBQL/ray_tracer.h @@ -73,11 +73,12 @@ class CuBQLRayTracer : public RayTracer { private: struct CuBQLSurfaceBVH { MeshID surface {ID_NONE}; + MeshID forward_parent {ID_NONE}; + MeshID reverse_parent {ID_NONE}; cuBQL::bvh3d bvh; cuBQL::vec3d* d_vertices {nullptr}; cuBQL::vec3i* d_indices {nullptr}; - cuBQL::vec3d* d_normals {nullptr}; MeshID* d_primitive_refs {nullptr}; uint32_t num_vertices {0}; @@ -93,6 +94,7 @@ class CuBQLRayTracer : public RayTracer { std::vector surface_bvhs_; std::unordered_map> tree_to_surface_bvh_indices_; + std::unordered_map surface_tree_to_volume_; }; } // namespace xdg diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index ea0e2216..22406943 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -34,9 +34,6 @@ CuBQLRayTracer::~CuBQLRayTracer() if (surface_bvh.d_indices) { omp_target_free(surface_bvh.d_indices, surface_bvh.gpu_id); } - if (surface_bvh.d_normals) { - omp_target_free(surface_bvh.d_normals, surface_bvh.gpu_id); - } if (surface_bvh.d_primitive_refs) { omp_target_free(surface_bvh.d_primitive_refs, surface_bvh.gpu_id); } @@ -78,6 +75,7 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man SurfaceTreeID tree = next_surface_tree_id(); surface_trees_.push_back(tree); + surface_tree_to_volume_[tree] = volume_id; auto& surface_bvh_indices = tree_to_surface_bvh_indices_[tree]; auto volume_surfaces = mesh_manager->get_volume_surfaces(volume_id); @@ -104,23 +102,8 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man std::vector h_primitive_refs; h_primitive_refs.reserve(num_faces); - std::vector h_normals; - h_normals.reserve(num_faces); - - auto [forward_parent, reverse_parent] = mesh_manager->get_parent_volumes(surf); for (const auto& face : mesh_manager->get_surface_faces(surf)) { h_primitive_refs.push_back(face); - - // TODO: This pre-orients normals for the current volume tree. If cuBQL - // starts reusing a single surface BVH/BLAS across both parent volumes, - // normals must be flipped at intersection time based on the queried volume. - auto normal = mesh_manager->face_normal(face); - if (volume_id == reverse_parent) { - normal = -normal; - } else if (volume_id != forward_parent) { - fatal_error("Volume {} is not a parent of surface {}", volume_id, surf); - } - h_normals.emplace_back(normal.x, normal.y, normal.z); } // copy vertices and indices to device @@ -134,11 +117,6 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man omp_target_memcpy(d_indices, h_indices.data(), h_indices.size() * sizeof(cuBQL::vec3i), 0, 0, gpu_id, host_id); - auto* d_normals = static_cast - (omp_target_alloc(h_normals.size() * sizeof(cuBQL::vec3d), gpu_id)); - omp_target_memcpy(d_normals, h_normals.data(), - h_normals.size() * sizeof(cuBQL::vec3d), 0, 0, gpu_id, host_id); - auto* d_primitive_refs = static_cast (omp_target_alloc(h_primitive_refs.size() * sizeof(MeshID), gpu_id)); omp_target_memcpy(d_primitive_refs, h_primitive_refs.data(), @@ -175,12 +153,19 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man omp_target_free(d_aabbs, gpu_id); + auto [forward_parent, reverse_parent] = mesh_manager->get_parent_volumes(surf); + CuBQLSurfaceBVH surface_bvh; + surface_bvh.surface = surf; + surface_bvh.forward_parent = forward_parent; + surface_bvh.reverse_parent = reverse_parent; + if (volume_id != forward_parent && volume_id != reverse_parent) { + fatal_error("Volume {} is not a parent of surface {}", volume_id, surf); + } surface_bvh.bvh = bvh; surface_bvh.d_vertices = d_vertices; surface_bvh.d_indices = d_indices; - surface_bvh.d_normals = d_normals; surface_bvh.d_primitive_refs = d_primitive_refs; surface_bvh.num_vertices = h_vertices.size(); surface_bvh.num_faces = num_faces; @@ -202,12 +187,12 @@ CuBQLRayTracer::create_element_tree(const std::shared_ptr&, void CuBQLRayTracer::create_global_surface_tree() { - fatal_error("Global surface trees not currently supported with cuBQL ray tracer"); + warning("Global surface trees not currently supported with cuBQL ray tracer"); } void CuBQLRayTracer::create_global_element_tree() { - fatal_error("Global element trees not currently supported with cuBQL ray tracer"); + warning("Global element trees not currently supported with cuBQL ray tracer"); } MeshID CuBQLRayTracer::find_element(const Position&) const @@ -263,6 +248,7 @@ CuBQLRayTracer::ray_fire(TreeID tree, auto* d_surface_hit = static_cast (omp_target_alloc(sizeof(CubqlHit), gpu_id)); + const MeshID query_volume = surface_tree_to_volume_.at(tree); const auto& surface_bvh_indices = tree_to_surface_bvh_indices_.at(tree); for (const auto surface_bvh_index : surface_bvh_indices) { const auto& surface_bvh = surface_bvhs_.at(surface_bvh_index); @@ -270,8 +256,8 @@ CuBQLRayTracer::ray_fire(TreeID tree, cuBQL::bvh3d bvh = surface_bvh.bvh; const cuBQL::vec3d* d_vertices = surface_bvh.d_vertices; const cuBQL::vec3i* d_indices = surface_bvh.d_indices; - const cuBQL::vec3d* d_normals = surface_bvh.d_normals; const MeshID* d_primitive_refs = surface_bvh.d_primitive_refs; + const bool reverse_sense = query_volume == surface_bvh.reverse_parent; CubqlHit surface_hit; surface_hit.distance = closest_distance + cubql_ray_hit_tolerance(closest_distance); @@ -288,7 +274,7 @@ CuBQLRayTracer::ray_fire(TreeID tree, const cuBQL::vec3d ray_direction(direction.x, direction.y, direction.z); #pragma omp target device(gpu_id) \ - is_device_ptr(d_vertices, d_indices, d_normals, d_primitive_refs, d_exclude_primitives, d_surface_hit) + is_device_ptr(d_vertices, d_indices, d_primitive_refs, d_exclude_primitives, d_surface_hit) { cuBQL::ray3d ray; ray.origin = ray_origin; @@ -318,7 +304,12 @@ CuBQLRayTracer::ray_fire(TreeID tree, d_vertices[index.z] }; - const double normal_dot_direction = dot(d_normals[prim_id], ray.direction); + cuBQL::vec3d normal = cuBQL::cross(vertices[1] - vertices[0], vertices[2] - vertices[0]); + if (reverse_sense) { + normal = -normal; + } + + const double normal_dot_direction = dot(normal, ray.direction); bool culled = false; if (orientation == static_cast(HitOrientation::EXITING)) { culled = normal_dot_direction < 0.0; From 1895ee0ae70b3a92303d8ec8b6474990a646e727 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 15 May 2026 14:39:48 +0100 Subject: [PATCH 20/22] Added an omp_target_device_probe and put it into the same header as the vk_device probe --- ...ulkan_probe.h => available_device_probe.h} | 48 ++++++++++++++++++- src/cuBQL/ray_tracer.cpp | 8 +++- src/gprt/ray_tracer.cpp | 2 +- tests/util.h | 4 +- 4 files changed, 57 insertions(+), 5 deletions(-) rename include/xdg/{gprt/vulkan_probe.h => available_device_probe.h} (64%) diff --git a/include/xdg/gprt/vulkan_probe.h b/include/xdg/available_device_probe.h similarity index 64% rename from include/xdg/gprt/vulkan_probe.h rename to include/xdg/available_device_probe.h index 8e0ff9cf..8aa6c5a1 100644 --- a/include/xdg/gprt/vulkan_probe.h +++ b/include/xdg/available_device_probe.h @@ -1,8 +1,12 @@ #pragma once -#ifdef XDG_ENABLE_GPRT #include "xdg/error.h" +// -------------------------------------------------------------------------------------- +// Vulkan probe functions to check for ray tracing capable devices at runtime +// -------------------------------------------------------------------------------------- +#ifdef XDG_ENABLE_GPRT + #include #include #include @@ -90,5 +94,47 @@ inline bool system_has_vk_device() missing); return false; } +#endif + +// -------------------------------------------------------------------------------------- +// OpenMP target probe functions to check for devices capable of running cuBQL at runtime +// -------------------------------------------------------------------------------------- + +#ifdef XDG_ENABLE_CUBQL + +#include + +inline bool system_has_omp_target_device() +{ + const int device_count = omp_get_num_devices(); + if (device_count <= 0) { + warning("No OpenMP target devices found; cuBQL ray tracer unavailable."); + return false; + } + + const int host_id = omp_get_initial_device(); + for (int device_id = 0; device_id < device_count; ++device_id) { + int value = 1; + void* d_value = omp_target_alloc(sizeof(int), device_id); + if (!d_value) continue; + + const int copy_result = omp_target_memcpy(d_value, + &value, + sizeof(int), + 0, + 0, + device_id, + host_id); + omp_target_free(d_value, device_id); + + if (copy_result == 0) { + write_message("Found OpenMP target device {}.", device_id); + return true; + } + } + + warning("OpenMP target devices were found, but none accepted target allocation; cuBQL ray tracer unavailable."); + return false; +} #endif diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index 22406943..4e3e7add 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -1,6 +1,7 @@ #include "xdg/cuBQL/ray_tracer.h" #include "xdg/error.h" #include "xdg/geometry/plucker.h" +#include "xdg/available_device_probe.h" #include #include "cuBQL/builder/omp.h" @@ -17,7 +18,12 @@ static inline __cubql_both double cubql_ray_hit_tolerance(double t) return tolerance * (1.0 + cuBQL::abst(t)); } -CuBQLRayTracer::CuBQLRayTracer() = default; +CuBQLRayTracer::CuBQLRayTracer() +{ + if (!system_has_omp_target_device()) { + fatal_error("No OpenMP target capable device found; cannot initialize cuBQL ray tracer."); + } +} CuBQLRayTracer::~CuBQLRayTracer() { diff --git a/src/gprt/ray_tracer.cpp b/src/gprt/ray_tracer.cpp index 15d21a7b..9db035ec 100644 --- a/src/gprt/ray_tracer.cpp +++ b/src/gprt/ray_tracer.cpp @@ -1,6 +1,6 @@ #include "xdg/gprt/ray_tracer.h" #include "gprt/gprt.h" -#include "xdg/gprt/vulkan_probe.h" +#include "xdg/available_device_probe.h" namespace xdg { diff --git a/tests/util.h b/tests/util.h index 074ff899..96946aa2 100644 --- a/tests/util.h +++ b/tests/util.h @@ -8,7 +8,7 @@ #include "xdg/constants.h" #include "xdg/ray_tracers.h" #include "xdg/mesh_managers.h" -#include "xdg/gprt/vulkan_probe.h" +#include "xdg/available_device_probe.h" namespace xdg::test { @@ -52,7 +52,7 @@ inline bool ray_tracer_available(xdg::RTLibrary rt) { case xdg::RTLibrary::CUBQL: #ifdef XDG_ENABLE_CUBQL - return true; + return system_has_omp_target_device(); #else return false; #endif From 3473b763372d2ef1129ae6878e8ae24944bbf6e0 Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 15 May 2026 16:38:16 +0100 Subject: [PATCH 21/22] Introduced a cubql_backend Context object to store host_id and device_id --- include/xdg/cuBQL/cuBQL_backend.h | 30 +++++++++++++++++++++++++++ include/xdg/cuBQL/ray_tracer.h | 17 ++++++---------- src/cuBQL/ray_tracer.cpp | 34 +++++++++++++++++-------------- 3 files changed, 55 insertions(+), 26 deletions(-) create mode 100644 include/xdg/cuBQL/cuBQL_backend.h diff --git a/include/xdg/cuBQL/cuBQL_backend.h b/include/xdg/cuBQL/cuBQL_backend.h new file mode 100644 index 00000000..b348802f --- /dev/null +++ b/include/xdg/cuBQL/cuBQL_backend.h @@ -0,0 +1,30 @@ +#ifndef _XDG_CUBQL_BACKEND_H +#define _XDG_CUBQL_BACKEND_H + +#include +#include +#include +#include + +#include + +// Guards to prevent CUDA headers from being included in host code, which causes +// failed compilation with LLVM-clang. +#if defined(__CUDA_ARCH__) && !defined(__CUDACC__) +#undef __CUDA_ARCH__ +#endif + +#include "xdg/error.h" + +namespace xdg::cubql { + +struct Context { + int gpuID {0}; + int hostID {omp_get_initial_device()}; +}; + + + +} // namespace xdg::cubql + +#endif // include guard diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h index 528e16ce..c737ded6 100644 --- a/include/xdg/cuBQL/ray_tracer.h +++ b/include/xdg/cuBQL/ray_tracer.h @@ -8,19 +8,11 @@ #include "xdg/constants.h" #include "xdg/geometry_data.h" +#include "xdg/cuBQL/triangles.h" #include "xdg/mesh_manager_interface.h" #include "xdg/ray.h" #include "xdg/ray_tracing_interface.h" -// Guards to prevent CUDA headers from being included in host code, which causes failed compilation with LLVM-clang -#if defined(__CUDA_ARCH__) && !defined(__CUDACC__) -#undef __CUDA_ARCH__ -#endif - -#include "cuBQL/bvh.h" - - - namespace xdg { class CuBQLRayTracer : public RayTracer { @@ -71,7 +63,9 @@ class CuBQLRayTracer : public RayTracer { double& dist) const override; private: - struct CuBQLSurfaceBVH { + cubql::Context context_; + + struct CuBQLSurfaceBLAS { MeshID surface {ID_NONE}; MeshID forward_parent {ID_NONE}; MeshID reverse_parent {ID_NONE}; @@ -85,6 +79,7 @@ class CuBQLRayTracer : public RayTracer { uint32_t num_faces {0}; int gpu_id {0}; }; + struct CubqlHit { double distance {INFTY}; @@ -92,7 +87,7 @@ class CuBQLRayTracer : public RayTracer { MeshID primitive {ID_NONE}; }; - std::vector surface_bvhs_; + std::vector surface_bvhs_; std::unordered_map> tree_to_surface_bvh_indices_; std::unordered_map surface_tree_to_volume_; }; diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index 4e3e7add..96705041 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -23,6 +23,9 @@ CuBQLRayTracer::CuBQLRayTracer() if (!system_has_omp_target_device()) { fatal_error("No OpenMP target capable device found; cannot initialize cuBQL ray tracer."); } + + context_.gpuID = 0; // TODO - support selecting among multiple OpenMP target devices. + context_.hostID = omp_get_initial_device(); } CuBQLRayTracer::~CuBQLRayTracer() @@ -64,8 +67,8 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man MeshID volume_id) { - int gpu_id = 0; // TODO - how to manage GPU IDs in a multi-GPU system? - int host_id = omp_get_initial_device(); + const auto& context = context_; + const int gpu_id = context.gpuID; // cuBQL::box3d *d_boxes = nullptr; // box3d is an alias for box_t @@ -116,17 +119,17 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man auto* d_vertices = static_cast (omp_target_alloc(h_vertices.size() * sizeof(cuBQL::vec3d), gpu_id)); omp_target_memcpy(d_vertices, h_vertices.data(), - h_vertices.size() * sizeof(cuBQL::vec3d), 0, 0, gpu_id, host_id); + h_vertices.size() * sizeof(cuBQL::vec3d), 0, 0, gpu_id, context.hostID); auto* d_indices = static_cast (omp_target_alloc(h_indices.size() * sizeof(cuBQL::vec3i), gpu_id)); omp_target_memcpy(d_indices, h_indices.data(), - h_indices.size() * sizeof(cuBQL::vec3i), 0, 0, gpu_id, host_id); + h_indices.size() * sizeof(cuBQL::vec3i), 0, 0, gpu_id, context.hostID); auto* d_primitive_refs = static_cast (omp_target_alloc(h_primitive_refs.size() * sizeof(MeshID), gpu_id)); omp_target_memcpy(d_primitive_refs, h_primitive_refs.data(), - h_primitive_refs.size() * sizeof(MeshID), 0, 0, gpu_id, host_id); + h_primitive_refs.size() * sizeof(MeshID), 0, 0, gpu_id, context.hostID); // Create device storage for triangle AABBs to be computed in parallel on GPU auto* d_aabbs = static_cast @@ -160,15 +163,14 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man omp_target_free(d_aabbs, gpu_id); auto [forward_parent, reverse_parent] = mesh_manager->get_parent_volumes(surf); + if (volume_id != forward_parent && volume_id != reverse_parent) { + fatal_error("Volume {} is not a parent of surface {}", volume_id, surf); + } - CuBQLSurfaceBVH surface_bvh; - + CuBQLSurfaceBLAS surface_bvh; surface_bvh.surface = surf; surface_bvh.forward_parent = forward_parent; surface_bvh.reverse_parent = reverse_parent; - if (volume_id != forward_parent && volume_id != reverse_parent) { - fatal_error("Volume {} is not a parent of surface {}", volume_id, surf); - } surface_bvh.bvh = bvh; surface_bvh.d_vertices = d_vertices; surface_bvh.d_indices = d_indices; @@ -180,6 +182,8 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man surface_bvh_indices.push_back(surface_bvhs_.size()); surface_bvhs_.push_back(surface_bvh); } + + return tree; } @@ -230,8 +234,8 @@ CuBQLRayTracer::ray_fire(TreeID tree, HitOrientation hitOrientation, std::vector* const exclude_primitives) { - int gpu_id = 0; // TODO - how to manage GPU IDs in a multi-GPU system? - int host_id = omp_get_initial_device(); + const auto& context = context_; + const int gpu_id = context.gpuID; MeshID* d_exclude_primitives = nullptr; int exclude_count = 0; @@ -245,7 +249,7 @@ CuBQLRayTracer::ray_fire(TreeID tree, 0, 0, gpu_id, - host_id); + context.hostID); } double closest_distance = tmax; @@ -273,7 +277,7 @@ CuBQLRayTracer::ray_fire(TreeID tree, 0, 0, gpu_id, - host_id); + context.hostID); const int orientation = static_cast(hitOrientation); const cuBQL::vec3d ray_origin(origin.x, origin.y, origin.z); @@ -352,7 +356,7 @@ CuBQLRayTracer::ray_fire(TreeID tree, sizeof(CubqlHit), 0, 0, - host_id, + context.hostID, gpu_id); // Temporary cuBQL-only closest-hit reduction while each surface BVH is queried From 0fd848ecfcc66ecfe7b45d56070079394978e50b Mon Sep 17 00:00:00 2001 From: waqar-ukaea Date: Fri, 15 May 2026 16:59:10 +0100 Subject: [PATCH 22/22] Implemented a version of DPRT's AutoUploadArray for RAII style containers for omp offload --- include/xdg/cuBQL/cuBQL_backend.h | 97 ++++++++++++++++++++++++++++++- include/xdg/cuBQL/ray_tracer.h | 8 +-- src/cuBQL/ray_tracer.cpp | 76 ++++++++---------------- 3 files changed, 123 insertions(+), 58 deletions(-) diff --git a/include/xdg/cuBQL/cuBQL_backend.h b/include/xdg/cuBQL/cuBQL_backend.h index b348802f..6a6c8111 100644 --- a/include/xdg/cuBQL/cuBQL_backend.h +++ b/include/xdg/cuBQL/cuBQL_backend.h @@ -16,15 +16,108 @@ #include "xdg/error.h" -namespace xdg::cubql { +namespace xdg_cubql_backend { struct Context { int gpuID {0}; int hostID {omp_get_initial_device()}; }; +template +struct AutoUploadArray { + AutoUploadArray() = default; + AutoUploadArray(Context context, const InputT* input, std::size_t count) + : context(context), count(count) + { + if (count == 0) return; -} // namespace xdg::cubql + elements = static_cast(omp_target_alloc(count * sizeof(T), context.gpuID)); + if (!elements) { + fatal_error("Failed to allocate cuBQL device array."); + } + needsFree = true; + + int result = 0; + if constexpr (std::is_same_v) { + result = omp_target_memcpy(elements, + input, + count * sizeof(T), + 0, + 0, + context.gpuID, + context.hostID); + } else { + std::vector converted(count); + for (std::size_t i = 0; i < count; ++i) { + converted[i] = T(input[i]); + } + + result = omp_target_memcpy(elements, + converted.data(), + count * sizeof(T), + 0, + 0, + context.gpuID, + context.hostID); + } + + if (result != 0) { + reset(); + fatal_error("Failed to upload cuBQL device array."); + } + } + + ~AutoUploadArray() + { + reset(); + } + + AutoUploadArray(const AutoUploadArray&) = delete; + AutoUploadArray& operator=(const AutoUploadArray&) = delete; + + AutoUploadArray(AutoUploadArray&& other) noexcept + { + *this = std::move(other); + } + + AutoUploadArray& operator=(AutoUploadArray&& other) noexcept + { + if (this != &other) { + reset(); + + context = other.context; + elements = other.elements; + count = other.count; + needsFree = other.needsFree; + + other.elements = nullptr; + other.count = 0; + other.needsFree = false; + } + + return *this; + } + + void reset() + { + if (needsFree && elements) { + omp_target_free(elements, context.gpuID); + } + + elements = nullptr; + count = 0; + needsFree = false; + } + + T* get() const { return elements; } + + T* elements {nullptr}; + std::size_t count {0}; + bool needsFree {false}; + Context context {}; +}; + +} // namespace xdg_cubql_backend #endif // include guard diff --git a/include/xdg/cuBQL/ray_tracer.h b/include/xdg/cuBQL/ray_tracer.h index c737ded6..6c8a6501 100644 --- a/include/xdg/cuBQL/ray_tracer.h +++ b/include/xdg/cuBQL/ray_tracer.h @@ -63,7 +63,7 @@ class CuBQLRayTracer : public RayTracer { double& dist) const override; private: - cubql::Context context_; + xdg_cubql_backend::Context context_; struct CuBQLSurfaceBLAS { MeshID surface {ID_NONE}; @@ -71,9 +71,9 @@ class CuBQLRayTracer : public RayTracer { MeshID reverse_parent {ID_NONE}; cuBQL::bvh3d bvh; - cuBQL::vec3d* d_vertices {nullptr}; - cuBQL::vec3i* d_indices {nullptr}; - MeshID* d_primitive_refs {nullptr}; + xdg_cubql_backend::AutoUploadArray vertices; + xdg_cubql_backend::AutoUploadArray indices; + xdg_cubql_backend::AutoUploadArray primitive_refs; uint32_t num_vertices {0}; uint32_t num_faces {0}; diff --git a/src/cuBQL/ray_tracer.cpp b/src/cuBQL/ray_tracer.cpp index 96705041..3660fb44 100644 --- a/src/cuBQL/ray_tracer.cpp +++ b/src/cuBQL/ray_tracer.cpp @@ -37,15 +37,6 @@ CuBQLRayTracer::~CuBQLRayTracer() cuBQL::omp::freeBVH(surface_bvh.bvh, &context); } - if (surface_bvh.d_vertices) { - omp_target_free(surface_bvh.d_vertices, surface_bvh.gpu_id); - } - if (surface_bvh.d_indices) { - omp_target_free(surface_bvh.d_indices, surface_bvh.gpu_id); - } - if (surface_bvh.d_primitive_refs) { - omp_target_free(surface_bvh.d_primitive_refs, surface_bvh.gpu_id); - } } } @@ -115,21 +106,9 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man h_primitive_refs.push_back(face); } - // copy vertices and indices to device - auto* d_vertices = static_cast - (omp_target_alloc(h_vertices.size() * sizeof(cuBQL::vec3d), gpu_id)); - omp_target_memcpy(d_vertices, h_vertices.data(), - h_vertices.size() * sizeof(cuBQL::vec3d), 0, 0, gpu_id, context.hostID); - - auto* d_indices = static_cast - (omp_target_alloc(h_indices.size() * sizeof(cuBQL::vec3i), gpu_id)); - omp_target_memcpy(d_indices, h_indices.data(), - h_indices.size() * sizeof(cuBQL::vec3i), 0, 0, gpu_id, context.hostID); - - auto* d_primitive_refs = static_cast - (omp_target_alloc(h_primitive_refs.size() * sizeof(MeshID), gpu_id)); - omp_target_memcpy(d_primitive_refs, h_primitive_refs.data(), - h_primitive_refs.size() * sizeof(MeshID), 0, 0, gpu_id, context.hostID); + xdg_cubql_backend::AutoUploadArray d_vertices(context, h_vertices.data(), h_vertices.size()); + xdg_cubql_backend::AutoUploadArray d_indices(context, h_indices.data(), h_indices.size()); + xdg_cubql_backend::AutoUploadArray d_primitive_refs(context, h_primitive_refs.data(), h_primitive_refs.size()); // Create device storage for triangle AABBs to be computed in parallel on GPU auto* d_aabbs = static_cast @@ -139,14 +118,16 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man // Create AABBs for each triangle in parallel on GPU // TODO - Should this be its own function? // TODO - How can this be extended for tets and other element types? - #pragma omp target device(gpu_id) is_device_ptr(d_vertices, d_indices, d_aabbs) + auto* vertices_ptr = d_vertices.get(); + auto* indices_ptr = d_indices.get(); + #pragma omp target device(gpu_id) is_device_ptr(vertices_ptr, indices_ptr, d_aabbs) #pragma omp teams distribute parallel for for (uint32_t primID = 0; primID < num_primitives; ++primID) { - cuBQL::vec3i indices = d_indices[primID]; + cuBQL::vec3i indices = indices_ptr[primID]; - cuBQL::vec3d A = d_vertices[indices.x]; - cuBQL::vec3d B = d_vertices[indices.y]; - cuBQL::vec3d C = d_vertices[indices.z]; + cuBQL::vec3d A = vertices_ptr[indices.x]; + cuBQL::vec3d B = vertices_ptr[indices.y]; + cuBQL::vec3d C = vertices_ptr[indices.z]; cuBQL::box3d aabb; aabb.extend(A); @@ -172,15 +153,15 @@ CuBQLRayTracer::create_surface_tree(const std::shared_ptr& mesh_man surface_bvh.forward_parent = forward_parent; surface_bvh.reverse_parent = reverse_parent; surface_bvh.bvh = bvh; - surface_bvh.d_vertices = d_vertices; - surface_bvh.d_indices = d_indices; - surface_bvh.d_primitive_refs = d_primitive_refs; + surface_bvh.vertices = std::move(d_vertices); + surface_bvh.indices = std::move(d_indices); + surface_bvh.primitive_refs = std::move(d_primitive_refs); surface_bvh.num_vertices = h_vertices.size(); surface_bvh.num_faces = num_faces; surface_bvh.gpu_id = gpu_id; surface_bvh_indices.push_back(surface_bvhs_.size()); - surface_bvhs_.push_back(surface_bvh); + surface_bvhs_.push_back(std::move(surface_bvh)); } @@ -237,19 +218,13 @@ CuBQLRayTracer::ray_fire(TreeID tree, const auto& context = context_; const int gpu_id = context.gpuID; - MeshID* d_exclude_primitives = nullptr; int exclude_count = 0; + xdg_cubql_backend::AutoUploadArray d_exclude_primitives; if (exclude_primitives && !exclude_primitives->empty()) { exclude_count = static_cast(exclude_primitives->size()); - d_exclude_primitives = static_cast - (omp_target_alloc(exclude_count * sizeof(MeshID), gpu_id)); - omp_target_memcpy(d_exclude_primitives, - exclude_primitives->data(), - exclude_count * sizeof(MeshID), - 0, - 0, - gpu_id, - context.hostID); + d_exclude_primitives = xdg_cubql_backend::AutoUploadArray(context, + exclude_primitives->data(), + exclude_count); } double closest_distance = tmax; @@ -264,9 +239,10 @@ CuBQLRayTracer::ray_fire(TreeID tree, const auto& surface_bvh = surface_bvhs_.at(surface_bvh_index); cuBQL::bvh3d bvh = surface_bvh.bvh; - const cuBQL::vec3d* d_vertices = surface_bvh.d_vertices; - const cuBQL::vec3i* d_indices = surface_bvh.d_indices; - const MeshID* d_primitive_refs = surface_bvh.d_primitive_refs; + const cuBQL::vec3d* d_vertices = surface_bvh.vertices.get(); + const cuBQL::vec3i* d_indices = surface_bvh.indices.get(); + const MeshID* d_primitive_refs = surface_bvh.primitive_refs.get(); + const MeshID* d_exclude_primitive_ptr = d_exclude_primitives.get(); const bool reverse_sense = query_volume == surface_bvh.reverse_parent; CubqlHit surface_hit; @@ -284,7 +260,7 @@ CuBQLRayTracer::ray_fire(TreeID tree, const cuBQL::vec3d ray_direction(direction.x, direction.y, direction.z); #pragma omp target device(gpu_id) \ - is_device_ptr(d_vertices, d_indices, d_primitive_refs, d_exclude_primitives, d_surface_hit) + is_device_ptr(d_vertices, d_indices, d_primitive_refs, d_exclude_primitive_ptr, d_surface_hit) { cuBQL::ray3d ray; ray.origin = ray_origin; @@ -302,7 +278,7 @@ CuBQLRayTracer::ray_fire(TreeID tree, const MeshID primitive_ref = d_primitive_refs[prim_id]; for (int i = 0; i < exclude_count; ++i) { - if (d_exclude_primitives[i] == primitive_ref) { + if (d_exclude_primitive_ptr[i] == primitive_ref) { return ray.tMax; } } @@ -377,10 +353,6 @@ CuBQLRayTracer::ray_fire(TreeID tree, omp_target_free(d_surface_hit, gpu_id); - if (d_exclude_primitives) { - omp_target_free(d_exclude_primitives, gpu_id); - } - if (closest_surface == ID_NONE) { return {INFTY, ID_NONE}; }