star-gf

commit f73c4fa7b729b19d304c973f22fb3bf9bfccce19
parent 2816d41307f8ef530049c5a123dae2dc401809f1
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Thu, 24 Sep 2015 14:21:33 +0200

Conditionally use a sparse data structure or not

Choose at compile time whether or not the Gebhart Factor matrix is
sparsely stored.

Diffstat:
M
src/sgf.c
 | 
75
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++-----
M
src/test_sgf_cube.c
 | 
1
-

2 files changed, 70 insertions(+), 6 deletions(-)
diff --git a/src/sgf.c b/src/sgf.c
@@ -48,10 +48,21 @@ struct accum {
 #define HTABLE_KEY uint64_t
 #include <rsys/hash_table.h>
 
+#define DARRAY_NAME accum
+#define DARRAY_DATA struct accum
+#include <rsys/dynamic_array.h>
+
+#define SPARSE
+
 #define CELL_KEY(From, To) (((uint64_t)From<<32) | (uint64_t)(To))
 
 struct sgf_device {
+#ifdef SPARSE
   struct htable_accum matrix;
+#else
+  struct darray_accum matrix;
+  size_t nprims;
+#endif
   struct darray_size_t nsteps;
   struct mutex* mutex;
 
@@ -88,7 +99,7 @@ accum_weight
   uint64_t key;
   res_T res = RES_OK;
   ASSERT(dev);
-
+#ifdef SPARSE
   key = CELL_KEY(ifrom, ito);
   mutex_lock(dev->mutex);
   acc = htable_accum_find(&dev->matrix, &key);
@@ -102,6 +113,14 @@ accum_weight
     res = htable_accum_set(&dev->matrix, &key, &tmp);
   }
   mutex_unlock(dev->mutex);
+#else
+  key = ifrom * dev->nprims + ito;
+  acc = darray_accum_data_get(&dev->matrix) + key;
+  mutex_lock(dev->mutex);
+  acc->radiative_flux += weight;
+  acc->sqr_radiative_flux += weight * weight;
+  mutex_unlock(dev->mutex);
+#endif
   return res;
 }
 
@@ -244,7 +263,11 @@ device_release(ref_T* ref)
   struct sgf_device* dev;
   ASSERT(ref);
   dev = CONTAINER_OF(ref, struct sgf_device, ref);
+#ifdef SPARSE
   htable_accum_release(&dev->matrix);
+#else
+  darray_accum_release(&dev->matrix);
+#endif
   darray_size_t_release(&dev->nsteps);
   if(dev->mutex) mutex_destroy(dev->mutex);
   MEM_RM(dev->allocator, dev);
@@ -286,7 +309,11 @@ sgf_device_create
   dev->allocator = mem_allocator;
   dev->logger = logger;
   dev->verbose = verbose;
+#ifdef SPARSE
   htable_accum_init(dev->allocator, &dev->matrix);
+#else
+  darray_accum_init(dev->allocator, &dev->matrix);
+#endif
   darray_size_t_init(dev->allocator, &dev->nsteps);
 
   dev->mutex = mutex_create();
@@ -363,7 +390,21 @@ sgf_begin_integrate(struct sgf_device* dev, struct sgf_scene_desc* desc)
   }
 
   memset(darray_size_t_data_get(&dev->nsteps), 0, nprims * sizeof(size_t));
+#ifdef SPARSE
   htable_accum_clear(&dev->matrix);
+#else
+  res = darray_accum_resize(&dev->matrix, nprims*nprims);
+  if(res != RES_OK) {
+    if(dev->verbose) {
+      logger_print(dev->logger, LOG_ERROR,
+"Couldn't allocate the Gebhart Factor integration data structure.\n");
+    }
+    S3D(scene_end_session(desc->scene));
+    return res;
+  }
+  memset(darray_accum_data_get(&dev->matrix), 0, nprims*nprims*sizeof(struct accum));
+  dev->nprims = nprims;
+#endif
   S3D(scene_ref_get(dev->scn_desc.scene));
   dev->is_integrating = 1;
 
@@ -394,7 +435,13 @@ sgf_end_integrate(struct sgf_device* dev)
 
   S3D(scene_ref_put(dev->scn_desc.scene));
   dev->is_integrating = 0;
+
+#ifdef SPARSE
   htable_accum_clear(&dev->matrix);
+#else
+  darray_accum_clear(&dev->matrix);
+#endif
+  darray_size_t_clear(&dev->nsteps);
   return RES_OK;
 }
 
@@ -433,7 +480,10 @@ sgf_integrate
     goto error;
   }
 
+  mutex_lock(dev->mutex);
   darray_size_t_data_get(&dev->nsteps)[primitive_id] += steps_count;
+  mutex_unlock(dev->mutex);
+
   FOR_EACH(i, 0, steps_count) {
     res = gebhart_radiative_path(dev, primitive_id, rng);
     ASSERT(res == RES_OK);
@@ -453,7 +503,10 @@ sgf_get
    struct sgf_estimator* estimator)
 {
   uint64_t key;
+#ifdef SPARSE
   struct accum* acc;
+#endif
+  struct accum acc_cp;
   size_t nprims;
 
   if(!dev || !estimator)
@@ -478,19 +531,31 @@ sgf_get
     return RES_BAD_ARG;
   }
 
+#ifdef SPARSE
   key = CELL_KEY(ifrom, ito);
+  mutex_lock(dev->mutex);
   acc = htable_accum_find(&dev->matrix, &key);
+  if(acc) acc_cp = *acc;
+  mutex_unlock(dev->mutex);
   if(!acc) {
     estimator->E = estimator->V = estimator->SE = 0.0;
     estimator->nsteps = 0;
-  } else {
+  } else
+#else
+  key = ifrom * dev->nprims + ito;
+  mutex_lock(dev->mutex);
+  acc_cp = darray_accum_data_get(&dev->matrix)[key];
+  mutex_unlock(dev->mutex);
+#endif
+  {
     size_t nsteps = darray_size_t_data_get(&dev->nsteps)[ifrom];
-    estimator->E = acc->radiative_flux / (double)nsteps;
+    estimator->E = acc_cp.radiative_flux / (double)nsteps;
     estimator->V =
-        acc->sqr_radiative_flux / (double)nsteps
-      - (acc->radiative_flux * acc->radiative_flux) / (double)(nsteps * nsteps);
+        acc_cp.sqr_radiative_flux / (double)nsteps
+      - (acc_cp.radiative_flux * acc_cp.radiative_flux) / (double)(nsteps * nsteps);
     estimator->SE = sqrt(estimator->V / (double)nsteps);
     estimator->nsteps = nsteps;
   }
   return RES_OK;
 }
+
diff --git a/src/test_sgf_cube.c b/src/test_sgf_cube.c
@@ -121,7 +121,6 @@ main(int argc, char** argv)
   CHECK(ssp_rng_proxy_create(&allocator, &ssp_rng_threefry, nbuckets, &proxy), RES_OK);
   CHECK(sgf_device_create(NULL, NULL, 1, &sgf), RES_OK);
 
-
   attribs[0].type = S3D_FLOAT3;
   attribs[0].usage = S3D_POSITION;
   attribs[0].get = get_pos;