Add missing free buffers

CMakeLists.txt

@@ -2,7 +2,7 @@ cmake_minimum_required(VERSION 3.20 FATAL_ERROR)
 project(
   ViennaRay
   LANGUAGES CXX
-  VERSION 3.8.0)
+  VERSION 3.8.1)
 
 # --------------------------------------------------------------------------------------------------------
 # Library switches
@@ -92,7 +92,7 @@ include("cmake/cpm.cmake")
 
 CPMAddPackage(
   NAME ViennaCore
-  VERSION 1.6.3
+  VERSION 1.7.0
   GIT_REPOSITORY "https://github.com/ViennaTools/ViennaCore"
   OPTIONS "VIENNACORE_USE_GPU ${VIENNARAY_USE_GPU}")
 

README.md

@@ -63,7 +63,7 @@ We recommend using [CPM.cmake](https://github.com/cpm-cmake/CPM.cmake) to consum
 * Installation with CPM
 
   ```cmake
-  CPMAddPackage("gh:viennatools/viennaray@3.8.0") # Use the latest release version
+  CPMAddPackage("gh:viennatools/viennaray@3.8.1") # Use the latest release version
   ```
 
 * With a local installation

cmake/generate_ptx.cmake

@@ -50,6 +50,7 @@ function(generate_kernel generated_files)
   cuda_include_directories(${ViennaCore_SOURCE_DIR}/include/viennacore)
   cuda_include_directories(${VIENNARAY_GPU_INCLUDE})
   cuda_include_directories(${OptiX_INCLUDE_DIR})
+  add_compile_definitions(VIENNACORE_COMPILE_GPU)
 
   cuda_compile_ptx(generated_ptx_files ${cu_source_files} ${cmake_options} ${options})
 

gpu/include/raygDiskGeometry.hpp

@@ -164,7 +164,7 @@ template <int D> struct DiskGeometry {
     CudaBuffer compactedSizeBuffer;
     compactedSizeBuffer.alloc(sizeof(uint64_t));
 
-    OptixAccelEmitDesc emitDesc;
+    OptixAccelEmitDesc emitDesc = {};
     emitDesc.type = OPTIX_PROPERTY_TYPE_COMPACTED_SIZE;
     emitDesc.result = compactedSizeBuffer.dPointer();
 
@@ -194,6 +194,8 @@ template <int D> struct DiskGeometry {
     outputBuffer.free(); // << the UNcompacted, temporary output buffer
     tempBuffer.free();
     compactedSizeBuffer.free();
+    d_aabbBuffer.free();
+    d_aabbBoundaryBuffer.free();
 
     launchParams.traversable = asHandle;
   }

gpu/include/raygLineGeometry.hpp

@@ -64,9 +64,9 @@ struct LineGeometry {
     }
 
     // Send AABB boxes to GPU
-    CudaBuffer d_aabbBuffer;
-    d_aabbBuffer.allocUpload(aabb);
-    CUdeviceptr d_aabb = d_aabbBuffer.dPointer();
+    CudaBuffer aabbBuffer;
+    aabbBuffer.allocUpload(aabb);
+    CUdeviceptr d_aabb = aabbBuffer.dPointer();
 
     // line inputs
     lineInput[0] = {};
@@ -141,7 +141,7 @@ struct LineGeometry {
     CudaBuffer compactedSizeBuffer;
     compactedSizeBuffer.alloc(sizeof(uint64_t));
 
-    OptixAccelEmitDesc emitDesc;
+    OptixAccelEmitDesc emitDesc = {};
     emitDesc.type = OPTIX_PROPERTY_TYPE_COMPACTED_SIZE;
     emitDesc.result = compactedSizeBuffer.dPointer();
 
@@ -171,6 +171,8 @@ struct LineGeometry {
     outputBuffer.free(); // << the UNcompacted, temporary output buffer
     tempBuffer.free();
     compactedSizeBuffer.free();
+    aabbBuffer.free();
+    d_aabbBoundaryBuffer.free();
 
     launchParams.traversable = asHandle;
   }
@@ -206,6 +208,7 @@ struct LineGeometry {
   void freeBuffers() {
     geometryNodesBuffer.free();
     geometryLinesBuffer.free();
+    geometryNormalsBuffer.free();
     boundaryNodesBuffer.free();
     boundaryLinesBuffer.free();
     asBuffer.free();

gpu/include/raygPerRayData.hpp

@@ -2,8 +2,7 @@
 
 #include <optix.h>
 
-#include "raygRNG.hpp"
-
+#include <vcRNG.hpp>
 #include <vcVectorType.hpp>
 
 #include <stdint.h>
@@ -27,7 +26,7 @@ struct PerRayData {
   float load = 0.f;
 
   // RNG
-  RNGState RNGstate;
+  CudaRNG RNGstate;
 
   // Hit data
   unsigned int numBoundaryHits = 0;

gpu/include/raygReflection.hpp

@@ -6,9 +6,9 @@
 #include <curand.h>
 
 #include "raygPerRayData.hpp"
-#include "raygRNG.hpp"
 #include "raygSBTRecords.hpp"
 
+#include <vcRNG.hpp>
 #include <vcVectorType.hpp>
 
 namespace viennaray::gpu {
@@ -64,7 +64,7 @@ specularReflection(PerRayData *prd, const Vec3Df &geoNormal) {
   prd->dir = prd->dir - (2 * DotProduct(prd->dir, geoNormal)) * geoNormal;
 }
 
-static __device__ Vec3Df PickRandomPointOnUnitSphere(RNGState *state) {
+static __device__ Vec3Df PickRandomPointOnUnitSphere(CudaRNG *state) {
   const float4 u = curand_uniform4(state); // (0,1]
   const float z = 1.0f - 2.0f * u.x;       // uniform in [-1,1]
   const float r2 = fmaxf(0.0f, 1.0f - z * z);

gpu/include/raygTrace.hpp

@@ -241,6 +241,10 @@ template <class T, int D> class Trace {
           .print();
     }
     cellDataBuffer_ = passedCellDataBuffer;
+#ifndef NDEBUG
+    // In debug mode, we set the buffer as reference to avoid accidental frees
+    cellDataBuffer_.isRef = true;
+#endif
     numCellData = numData;
   }
 

gpu/include/raygTriangleGeometry.hpp

@@ -145,7 +145,7 @@ struct TriangleGeometry {
     CudaBuffer compactedSizeBuffer;
     compactedSizeBuffer.alloc(sizeof(uint64_t));
 
-    OptixAccelEmitDesc emitDesc;
+    OptixAccelEmitDesc emitDesc = {};
     emitDesc.type = OPTIX_PROPERTY_TYPE_COMPACTED_SIZE;
     emitDesc.result = compactedSizeBuffer.dPointer();
 
@@ -269,6 +269,7 @@ struct TriangleGeometry {
   void freeBuffers() {
     geometryIndexBuffer.free();
     geometryVertexBuffer.free();
+    geometryNormalBuffer.free();
     boundaryIndexBuffer.free();
     boundaryVertexBuffer.free();
     asBuffer.free();

gpu/pipelines/GeneralPipelineDisk.cu

@@ -8,7 +8,6 @@
 #include <raygCallableConfig.hpp>
 #include <raygLaunchParams.hpp>
 #include <raygPerRayData.hpp>
-#include <raygRNG.hpp>
 #include <raygReflection.hpp>
 #include <raygSBTRecords.hpp>
 #include <raygSource.hpp>

gpu/pipelines/GeneralPipelineLine.cu

@@ -8,7 +8,6 @@
 #include <raygCallableConfig.hpp>
 #include <raygLaunchParams.hpp>
 #include <raygPerRayData.hpp>
-#include <raygRNG.hpp>
 #include <raygReflection.hpp>
 #include <raygSBTRecords.hpp>
 #include <raygSource.hpp>

gpu/pipelines/GeneralPipelineTriangle.cu

@@ -8,7 +8,6 @@
 #include <raygCallableConfig.hpp>
 #include <raygLaunchParams.hpp>
 #include <raygPerRayData.hpp>
-#include <raygRNG.hpp>
 #include <raygReflection.hpp>
 #include <raygSBTRecords.hpp>
 #include <raygSource.hpp>

include/viennaray/rayTraceKernel.hpp

@@ -222,10 +222,9 @@ template <typename NumericType, int D, GeometryType geoType> class TraceKernel {
 
           // Calculate point of impact
           const auto &ray = rayHit.ray;
-          const auto hitPoint =
-              Vec3D<rtcNumericType>{ray.org_x + ray.dir_x * ray.tfar,
-                                    ray.org_y + ray.dir_y * ray.tfar,
-                                    ray.org_z + ray.dir_z * ray.tfar};
+          const auto hitPoint = Vec3Df{ray.org_x + ray.dir_x * ray.tfar,
+                                       ray.org_y + ray.dir_y * ray.tfar,
+                                       ray.org_z + ray.dir_z * ray.tfar};
 
           const auto rayDir =
               Vec3D<NumericType>{ray.dir_x, ray.dir_y, ray.dir_z};
@@ -239,6 +238,7 @@ template <typename NumericType, int D, GeometryType geoType> class TraceKernel {
               if (hitFromBack) {
                 // if hitFromBack == true, then the ray hits the back of a disk
                 // the second time. In this case we discard the ray.
+                ++raysTerminated;
                 break;
               }
               hitFromBack = true;
@@ -263,8 +263,7 @@ template <typename NumericType, int D, GeometryType geoType> class TraceKernel {
                                  // origins of hit disks
             {                    // distance on first disk hit
               const auto &disk = geometry_.getPrimRef(rayHit.hit.primID);
-              const auto &diskOrigin =
-                  *reinterpret_cast<Vec3D<rtcNumericType> const *>(&disk);
+              const auto &diskOrigin = *reinterpret_cast<Vec3Df const *>(&disk);
               impactDistances.push_back(
                   Distance(hitPoint, diskOrigin) +
                   1e-6f); // add eps to avoid division by 0
@@ -323,6 +322,7 @@ template <typename NumericType, int D, GeometryType geoType> class TraceKernel {
           }
           if (++numReflections > config_.maxReflections) {
             // terminate ray if too many reflections
+            ++raysTerminated;
             break;
           }
           reflect = rejectionControl(rayWeight, initialRayWeight, rngState);

include/viennaray/rayUtil.hpp

@@ -267,14 +267,15 @@ template <typename NumericType>
 [[nodiscard]] static Vec3D<NumericType>
 pickRandomPointOnUnitSphere(RNG &rngState) {
   static thread_local std::uniform_real_distribution<NumericType> uniDist(
-      NumericType(0), NumericType(1));
-  NumericType x, y, z, x2py2;
+      NumericType(-1), NumericType(1));
+  NumericType x, y, z;
+  double x2py2;
   do {
-    x = 2 * uniDist(rngState) - 1.;
-    y = 2 * uniDist(rngState) - 1.;
+    x = uniDist(rngState);
+    y = uniDist(rngState);
     x2py2 = x * x + y * y;
   } while (x2py2 >= 1.);
-  NumericType tmp = 2 * std::sqrt(1. - x2py2);
+  double tmp = 2. * std::sqrt(1. - x2py2);
   x *= tmp;
   y *= tmp;
   z = 1. - 2 * x2py2;
@@ -297,9 +298,7 @@ template <typename T>
     return B;
   }
   const T invLen = T(1) / std::sqrt(len2);
-  u[0] *= invLen;
-  u[1] *= invLen;
-  u[2] *= invLen;
+  u = u * invLen;
   B[0] = u;
 
   // 2) choose a helper vector not collinear with u
@@ -313,14 +312,12 @@ template <typename T>
     h = Vec3D<T>{T(0), -u[2], u[1]};
   }
 
-  // 3) v = normalized(h)
-  auto v = h;
-  Normalize(v);
-  B[1] = v;
+  // 3) normalize h
+  Normalize(h);
+  B[1] = h;
 
   // 4) w = u × v  (already unit-length up to tiny FP error)
-  auto w = CrossProduct(u, v);
-  B[2] = w;
+  B[2] = CrossProduct(u, h);
 
   return B;
 }