commit 05fd610446aef20b2b554a377d95ae3b37c7b047
parent f04d37f62f0ab35a2f0092bc94263155171db313
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Wed, 23 Sep 2020 16:35:20 +0200
Implement the suvm_tetrahedral_mesh_create function
Diffstat:
3 files changed, 486 insertions(+), 3 deletions(-)
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
@@ -46,7 +46,8 @@ set(VERSION_PATCH 0)
set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH})
set(SUVM_FILES_SRC
- suvm_device.c)
+ suvm_device.c
+ suvm_volume.c)
set(SUVM_FILES_INC
suvm_device.h)
set(SUVM_FILES_INC_API
diff --git a/src/suvm.h b/src/suvm.h
@@ -87,8 +87,8 @@ suvm_device_ref_put
SUVM_API res_T
suvm_tetrahedral_mesh_create
(struct suvm_device* dev,
- const size_t ntethras, /* #tetrahedrals */
- void (*get_indices)(const size_t itetha, size_t ids[4], void* ctx),
+ const size_t ntetras, /* #tetrahedrals */
+ void (*get_indices)(const size_t itetra, size_t ids[4], void* ctx),
const struct suvm_data* tetra_data, /* NULL <=> no tetra data */
const size_t nverts, /* #vertices */
void (*get_position)(const size_t ivert, const double pos[3], void* ctx),
diff --git a/src/suvm_volume.c b/src/suvm_volume.c
@@ -0,0 +1,482 @@
+/* Copyright (C) 2020 |Meso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "suvm.h"
+#include "suvm_c.h"
+#include "suvm_device.h"
+
+#include <rsys/cstr.h>
+#include <rsys/dynamic_array_double.h>
+#include <rsys/dynamic_array_size_t.h>
+#include <rsys/ref_count.h>
+
+#include <embree3/rtcore_device.h>
+#include <embree3/rtcore_builder.h>
+
+/* Generate the dynamic array of RTCBuildPrimitive */
+#define DARRAY_NAME rtc_prim
+#define DARRAY_DATA struct RTCBuildPrimitive
+#define DARRAY_ALIGNMENT ALIGNOF(struct RTCBuildPrimitive)
+#include <rsys/dynamic_array.h>
+
+struct buffer {
+ void* mem; /* Raw memory block storing the data */
+ size_t elmt_size; /* Size in bytes of one data */
+ size_t elmt_stride; /* Size in bytes between two consecutive data */
+ size_t elmt_alignment; /* Data alignment */
+ size_t size; /* #data */
+ struct mem_allocator* allocator;
+};
+static const struct buffer BUFFER_NULL;
+
+enum node_type { NODE_INNER, NODE_LEAF };
+
+/* Generic node */
+struct node {
+ enum node_type type;
+};
+
+/* Inner node */
+struct ALIGN(16) node_inner {
+ float low[2][3];
+ float upp[2][3];
+ struct node* children[2];
+ struct node node;
+};
+
+/* Leaf node */
+struct ALIGN(16) node_leaf {
+ float low[3];
+ unsigned int geom_id;
+ float upp[3];
+ unsigned int prim_id;
+ struct node node;
+};
+
+struct suvm_volume {
+ struct darray_size_t indices; /* List of size_t[4] */
+ struct darray_double positions; /* List of double[3] */
+ struct buffer prim_data; /* Per primitive data */
+ struct buffer vert_data; /* Per vertex data */
+
+ RTCBVH bvh;
+ struct node* bvh_root;
+
+ struct suvm_device* dev;
+ ref_T ref;
+};
+
+/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+static INLINE void
+buffer_release(struct buffer* buf)
+{
+ ASSERT(buf);
+ if(buf->mem) MEM_RM(buf->allocator, buf->mem);
+ *buf = BUFFER_NULL;
+}
+
+static res_T
+buffer_init_from_data
+ (struct mem_allocator* mem_allocator, /* May be NULL */
+ struct buffer* buf,
+ const struct suvm_data* data,
+ const size_t ndata)
+{
+ struct mem_allocator* allocator = NULL;
+ size_t idata;
+ res_T res = RES_OK;
+ ASSERT(buf && data && data->get && ndata);
+
+ allocator = mem_allocator ? mem_allocator : &mem_default_allocator;
+
+ *buf = BUFFER_NULL;
+ if(!data->size || !IS_POW2(data->alignment)) {
+ res = RES_BAD_ARG;
+ goto error;
+ }
+ buf->elmt_size = data->size;
+ buf->elmt_alignment = data->alignment;
+ buf->elmt_stride = ALIGN_SIZE(data->size, data->alignment);
+ buf->size = ndata;
+ buf->allocator = allocator;
+
+ /* Allocate the required memory */
+ buf->mem = MEM_ALLOC_ALIGNED
+ (buf->allocator, buf->elmt_stride*ndata, buf->elmt_alignment);
+ if(!buf->mem) {
+ res = RES_MEM_ERR;
+ goto error;
+ }
+
+ /* Fill the buffer with the submitted data */
+ FOR_EACH(idata, 0, ndata) {
+ void* elmt = (char*)buf->mem + idata*buf->elmt_stride;
+ data->get(idata, 0, elmt);
+ }
+
+exit:
+ return res;
+error:
+ buffer_release(buf);
+ goto exit;
+}
+
+static INLINE size_t
+volume_get_primitives_count(const struct suvm_volume* vol)
+{
+ ASSERT(vol);
+ return darray_size_t_size_get(&vol->indices) / 4/*#indices per primitive*/;
+}
+
+static res_T
+setup_tetrahedral_mesh
+ (struct suvm_volume* vol,
+ const size_t ntetras,
+ void (*get_indices)(const size_t itetra, size_t ids[4], void* ctx),
+ const struct suvm_data* tetra_data,
+ const size_t nverts,
+ void (*get_position)(const size_t ivert, const double pos[3], void* ctx),
+ const struct suvm_data* vert_data,
+ void* context)
+{
+ size_t itetra;
+ size_t ivert;
+ res_T res = RES_OK;
+ ASSERT(vol && ntetras && get_indices && nverts && get_position);
+
+ /* Store the tetrahedral indices */
+ res = darray_size_t_resize(&vol->indices, ntetras*4/*#vertices per tetra*/);
+ if(res != RES_OK) {
+ log_err(vol->dev, "Could not allocate list of tetrahedra indices -- %s.\n",
+ res_to_cstr(res));
+ goto error;
+ }
+ FOR_EACH(itetra, 0, ntetras) {
+ size_t* tetra = darray_size_t_data_get(&vol->indices)+itetra*4;
+ get_indices(itetra, tetra, context);
+ }
+
+ /* Store the vertex positions */
+ res = darray_double_resize(&vol->positions, nverts*3/*#coords per vertex*/);
+ if(res != RES_OK) {
+ log_err(vol->dev, "Could not allocate the list tetrahedra vertices -- %s.\n",
+ res_to_cstr(res));
+ goto error;
+ }
+ FOR_EACH(ivert, 0, nverts) {
+ double* vert = darray_double_data_get(&vol->positions)+ivert*3;
+ get_position(ivert, vert, context);
+ }
+
+ /* Store the per tetrahedral data */
+ if(tetra_data && tetra_data->get) {
+ res = buffer_init_from_data
+ (vol->dev->allocator, &vol->prim_data, tetra_data, ntetras);
+ if(res != RES_OK) {
+ log_err(vol->dev,
+ "Could not setup the per volumetric primitive data -- %s.\n",
+ res_to_cstr(res));
+ goto error;
+ }
+ }
+
+ /* Store the per vertex data */
+ if(vert_data && vert_data->get) {
+ res = buffer_init_from_data
+ (vol->dev->allocator, &vol->vert_data, vert_data, nverts);
+ if(res != RES_OK) {
+ log_err(vol->dev, "Could not setup the per vertex data -- %s.\n",
+ res_to_cstr(res));
+ goto error;
+ }
+ }
+exit:
+ return res;
+error:
+ darray_size_t_purge(&vol->indices);
+ darray_double_purge(&vol->positions);
+ buffer_release(&vol->prim_data);
+ buffer_release(&vol->vert_data);
+ goto exit;
+}
+
+static INLINE void*
+rtc_node_inner_create
+ (RTCThreadLocalAllocator allocator, unsigned int nchildren, void* ctx)
+{
+ struct node_inner* inner = NULL;
+ ASSERT(nchildren == 2);
+ (void)ctx, (void)nchildren;
+ inner = rtcThreadLocalAlloc(allocator, sizeof(*inner), 16);
+ if(!inner) return NULL;
+ inner->node.type = NODE_INNER;
+ return &inner->node;
+}
+
+static INLINE void
+rtc_node_inner_set_children
+ (void* ptr, void* children[2], unsigned nchildren, void* ctx)
+{
+ struct node_inner* inner = CONTAINER_OF(ptr, struct node_inner, node);
+ struct node* node = ptr;
+ ASSERT(node && node->type == NODE_INNER && children && nchildren == 2);
+ (void)ctx, (void)nchildren, (void)node;
+ inner->children[0] = children[0];
+ inner->children[1] = children[1];
+}
+
+static INLINE void
+rtc_node_inner_set_bounds
+ (void* ptr,
+ const struct RTCBounds* bounds[2],
+ unsigned nchildren,
+ void* ctx)
+{
+ struct node_inner* inner = CONTAINER_OF(ptr, struct node_inner, node);
+ struct node* node = ptr;
+ ASSERT(node && node->type == NODE_INNER && bounds && nchildren == 2);
+ (void)ctx, (void)nchildren, (void)node;
+
+ /* Setup the AABB of the 1st child */
+ inner->low[0][0] = bounds[0]->lower_x;
+ inner->low[0][1] = bounds[0]->lower_y;
+ inner->low[0][2] = bounds[0]->lower_z;
+ inner->upp[0][0] = bounds[0]->upper_x;
+ inner->upp[0][1] = bounds[0]->upper_y;
+ inner->upp[0][2] = bounds[0]->upper_z;
+
+ /* Setup the AABB of the 2nd child */
+ inner->low[1][0] = bounds[1]->lower_x;
+ inner->low[1][1] = bounds[1]->lower_y;
+ inner->low[1][2] = bounds[1]->lower_z;
+ inner->upp[1][0] = bounds[1]->upper_x;
+ inner->upp[1][1] = bounds[1]->upper_y;
+ inner->upp[1][2] = bounds[1]->upper_z;
+}
+
+static INLINE void*
+rtc_node_leaf_create
+ (RTCThreadLocalAllocator allocator,
+ const struct RTCBuildPrimitive* prim,
+ size_t nprims,
+ void* ctx)
+{
+ struct node_leaf* leaf = NULL;
+ ASSERT(prim && nprims == 1);
+ (void)ctx, (void)nprims;
+
+ leaf = rtcThreadLocalAlloc(allocator, sizeof(*leaf), 16);
+ if(!leaf) return NULL;
+
+ leaf->low[0] = prim->lower_x;
+ leaf->low[1] = prim->lower_y;
+ leaf->low[2] = prim->lower_z;
+ leaf->upp[0] = prim->upper_x;
+ leaf->upp[1] = prim->upper_y;
+ leaf->upp[2] = prim->upper_z;
+ leaf->geom_id = prim->geomID;
+ leaf->prim_id = prim->primID;
+ leaf->node.type = NODE_LEAF;
+ return leaf;
+}
+
+static res_T
+build_bvh(struct suvm_volume* vol)
+{
+ struct darray_rtc_prim rtc_prims;
+ struct RTCBuildArguments args;
+ size_t iprim, nprims;
+ res_T res = RES_OK;
+ ASSERT(vol);
+
+ nprims = volume_get_primitives_count(vol);
+
+ /* Create the BVH */
+ vol->bvh = rtcNewBVH(vol->dev->rtc);
+ if(!vol->bvh) {
+ const enum RTCError rtc_err = rtcGetDeviceError(vol->dev->rtc);
+ log_err(vol->dev, "Could not create the BVH -- %s.\n",
+ rtc_error_string(rtc_err));
+ res = rtc_error_to_res_T(rtc_err);
+ goto error;
+ }
+
+ /* Allocate the array of geometric primitives */
+ darray_rtc_prim_init(vol->dev->allocator, &rtc_prims);
+ res = darray_rtc_prim_resize(&rtc_prims, nprims);
+ if(res != RES_OK) goto error;
+
+ /* Setup the primitive array */
+ FOR_EACH(iprim, 0, nprims) {
+ const size_t* ids = darray_size_t_cdata_get(&vol->indices) + iprim*4;
+ const double* verts[4];
+ struct RTCBuildPrimitive* prim = darray_rtc_prim_data_get(&rtc_prims)+iprim;
+ double low[3], upp[3];
+
+ /* Fetch the tetrahedral vertices */
+ verts[0] = darray_double_cdata_get(&vol->positions) + ids[0]*3/*#coords*/;
+ verts[1] = darray_double_cdata_get(&vol->positions) + ids[1]*3/*#coords*/;
+ verts[2] = darray_double_cdata_get(&vol->positions) + ids[2]*3/*#coords*/;
+ verts[3] = darray_double_cdata_get(&vol->positions) + ids[3]*3/*#coords*/;
+
+ /* Compute the tetrahedral AABB */
+ low[0] = MMIN(MMIN(verts[0][0],verts[1][0]), MMIN(verts[2][0],verts[3][0]));
+ low[1] = MMIN(MMIN(verts[0][1],verts[1][1]), MMIN(verts[2][1],verts[3][1]));
+ low[2] = MMIN(MMIN(verts[0][2],verts[1][2]), MMIN(verts[2][2],verts[3][2]));
+ upp[0] = MMAX(MMAX(verts[0][0],verts[1][0]), MMAX(verts[2][0],verts[3][0]));
+ upp[1] = MMAX(MMAX(verts[0][1],verts[1][1]), MMAX(verts[2][1],verts[3][1]));
+ upp[2] = MMAX(MMAX(verts[0][2],verts[1][2]), MMAX(verts[2][2],verts[3][2]));
+
+ /* Setup the build primitive */
+ prim->lower_x = (float)low[0];
+ prim->lower_y = (float)low[1];
+ prim->lower_z = (float)low[2];
+ prim->upper_x = (float)upp[0];
+ prim->upper_x = (float)upp[1];
+ prim->upper_x = (float)upp[2];
+ prim->geomID = 0;
+ prim->primID = (unsigned int)iprim;
+ }
+
+ /* Setup the build arguments */
+ args = rtcDefaultBuildArguments();
+ args.byteSize = sizeof(args);
+ args.buildQuality = RTC_BUILD_QUALITY_MEDIUM;
+ args.buildFlags = RTC_BUILD_FLAG_NONE;
+ args.maxBranchingFactor = 2;
+ args.maxDepth = 1024;
+ args.sahBlockSize = 1;
+ args.minLeafSize = 1;
+ args.maxLeafSize = 1;
+ args.traversalCost = 1.f;
+ args.intersectionCost = 10.f;
+ args.bvh = vol->bvh;
+ args.primitives = darray_rtc_prim_data_get(&rtc_prims);
+ args.primitiveCount = nprims;
+ args.primitiveArrayCapacity = darray_rtc_prim_capacity(&rtc_prims);
+ args.createNode = rtc_node_inner_create;
+ args.setNodeChildren = rtc_node_inner_set_children;
+ args.setNodeBounds = rtc_node_inner_set_bounds;
+ args.createLeaf = rtc_node_leaf_create;
+ args.splitPrimitive = NULL;
+ args.buildProgress = NULL;
+ args.userPtr = vol;
+
+ /* Build the BVH */
+ vol->bvh_root = rtcBuildBVH(&args);
+ if(!vol->bvh_root) {
+ const enum RTCError rtc_err = rtcGetDeviceError(vol->dev->rtc);
+ log_err(vol->dev, "Error building the BVH -- %s.\n",
+ rtc_error_string(rtc_err));
+ res = rtc_error_to_res_T(rtc_err);
+ goto error;
+ }
+
+exit:
+ darray_rtc_prim_release(&rtc_prims);
+ return res;
+error:
+ if(vol->bvh) {
+ rtcReleaseBVH(vol->bvh);
+ vol->bvh = NULL;
+ }
+ goto exit;
+}
+
+static void
+volume_release(ref_T* ref)
+{
+ struct suvm_volume* vol = NULL;
+ struct suvm_device* dev = NULL;
+ ASSERT(ref);
+ vol = CONTAINER_OF(ref, struct suvm_volume, ref);
+ dev = vol->dev;
+ darray_size_t_release(&vol->indices);
+ darray_double_release(&vol->positions);
+ buffer_release(&vol->prim_data);
+ buffer_release(&vol->vert_data);
+ if(vol->bvh) rtcReleaseBVH(vol->bvh);
+ MEM_RM(dev->allocator, vol);
+ SUVM(device_ref_put(dev));
+}
+
+/*******************************************************************************
+ * Exported functions
+ ******************************************************************************/
+res_T
+suvm_tetrahedral_mesh_create
+ (struct suvm_device* dev,
+ const size_t ntetras,
+ void (*get_indices)(const size_t itetha, size_t ids[4], void* ctx),
+ const struct suvm_data* tetra_data,
+ const size_t nverts,
+ void (*get_position)(const size_t ivert, const double pos[3], void* ctx),
+ const struct suvm_data* vert_data,
+ void* context,
+ struct suvm_volume** out_vol)
+{
+ struct suvm_volume* vol = NULL;
+ res_T res = RES_OK;
+
+ if(!dev || !ntetras || !get_indices || !nverts || !get_position || !out_vol) {
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ vol = MEM_CALLOC(dev->allocator, 1, sizeof(*vol));
+ if(!vol) {
+ res = RES_MEM_ERR;
+ goto error;
+ }
+ ref_init(&vol->ref);
+ SUVM(device_ref_get(dev));
+ vol->dev = dev;
+ darray_size_t_init(dev->allocator, &vol->indices);
+ darray_double_init(dev->allocator, &vol->positions);
+
+ /* Locally copy the volumetric mesh data */
+ res = setup_tetrahedral_mesh(vol, ntetras, get_indices, tetra_data, nverts,
+ get_position, vert_data, context);
+ if(res != RES_OK) goto error;
+
+ res = build_bvh(vol);
+ if(res != RES_OK) goto error;
+
+exit:
+ if(out_vol) *out_vol = vol;
+ return res;
+error:
+ if(vol) { SUVM(volume_ref_put(vol)); vol = NULL; }
+ goto exit;
+}
+
+res_T
+suvm_volume_ref_get(struct suvm_volume* vol)
+{
+ if(!vol) return RES_BAD_ARG;
+ ref_get(&vol->ref);
+ return RES_OK;
+}
+
+res_T
+suvm_volume_ref_put(struct suvm_volume* vol)
+{
+ if(!vol) return RES_BAD_ARG;
+ ref_put(&vol->ref, volume_release);
+ return RES_OK;
+}
+
