star-mc

commit cb47445d8a3ff2cdcad63f5b082ba59c032b315e
parent 6729159d8ca838ed01a8c385c3ecd02e89eb62d5
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Mon, 27 Sep 2021 15:24:23 +0200

Add progress log to the solve functions

Diffstat:
M
src/smc_estimator.c
 | 
135
++++++++++++++++++++++++++++++++++++++++++++++++++-----------------------------
M
src/test_smc_light_path.c
 | 
4
++--

2 files changed, 88 insertions(+), 51 deletions(-)
diff --git a/src/smc_estimator.c b/src/smc_estimator.c
@@ -33,6 +33,7 @@
 #include "smc_device_c.h"
 #include "smc_type_c.h"
 
+#include <rsys/logger.h>
 #include <rsys/mem_allocator.h>
 
 #include <limits.h>
@@ -157,9 +158,11 @@ smc_solve
   void** accums = NULL;
   void** accums_sqr = NULL;
   size_t* nsamples = NULL;
-  res_T res = RES_OK;
   size_t nfailed = 0;
-  int cancel = 0;
+  int progress = 0;
+  ATOMIC cancel = 0;
+  ATOMIC nsolved_realisations = 0;
+  res_T res = RES_OK;
 
   if(!dev || !integrator || !out_estimator || !check_integrator(integrator)) {
     res = RES_BAD_ARG;
@@ -198,30 +201,46 @@ smc_solve
   if(res != RES_OK) goto error;
 
   /* Parallel evaluation of the simulation */
-  #pragma omp parallel shared(cancel, nfailed)
-  {
+  logger_print(dev->logger, LOG_OUTPUT, "Solving: %3d%%\r", progress);
+  #pragma omp parallel for schedule(static)
+  for(i = 0; i < (int64_t)integrator->max_realisations; ++i) {
     const int ithread = omp_get_thread_num();
-    #pragma omp for schedule(static)
-    for(i = 0; i < (int64_t)integrator->max_realisations; ++i) {
-      if(!cancel) {
-        const res_T res_local = integrator->integrand
-          (vals[ithread], dev->rngs[ithread], (unsigned)ithread, ctx);
-        if(res_local == RES_OK) {
-          /* call succeded */
-          integrator->type->add(accums[ithread], accums[ithread], vals[ithread]);
-          integrator->type->mul(vals[ithread], vals[ithread], vals[ithread]);
-          integrator->type->add(accums_sqr[ithread], accums_sqr[ithread], vals[ithread]);
-          ++nsamples[ithread];
-        } else {
-          #pragma omp atomic
-          ++nfailed;
-          if (nfailed > integrator->max_failures) {
-            cancel = 1;
-          }
+    int64_t n = 0;
+    int pcent = 0;
+    res_T res_local = RES_OK;
+
+    if(ATOMIC_GET(&cancel)) continue;
+
+    res_local = integrator->integrand
+        (vals[ithread], dev->rngs[ithread], (unsigned)ithread, ctx);
+
+    if(res_local != RES_OK) {
+      #pragma omp critical
+      {
+        nfailed += 1;
+        if(nfailed >= integrator->max_failures) {
+          ATOMIC_SET(&cancel, 1);
         }
       }
+      continue;
+    }
+
+    /* call succeded */
+    integrator->type->add(accums[ithread], accums[ithread], vals[ithread]);
+    integrator->type->mul(vals[ithread], vals[ithread], vals[ithread]);
+    integrator->type->add(accums_sqr[ithread], accums_sqr[ithread], vals[ithread]);
+    ++nsamples[ithread];
+
+    n = ATOMIC_INCR(&nsolved_realisations);
+    pcent = (int)
+      ((double)n * 100.0 / (double)integrator->max_realisations + 0.5/*round*/);
+    #pragma omp critical
+    if(pcent > progress) {
+      progress = pcent;
+      logger_print(dev->logger, LOG_OUTPUT, "Solving: %3d%%\r", progress);
     }
   }
+  logger_print(dev->logger, LOG_OUTPUT, "Solving: %3d%%\n", progress);
 
   /* Merge the parallel estimation into the final estimator */
   FOR_EACH(i, 0, (int64_t)nthreads) {
@@ -260,11 +279,13 @@ smc_solve_N
    const size_t sizeof_ctx,
    struct smc_estimator* estimators[])
 {
-  res_T res = RES_OK;
+  void** vals = NULL;
   int64_t i;
   unsigned nthreads = 0;
-  void** vals = NULL;
-  int cancel = 0;
+  int progress = 0;
+  ATOMIC cancel = 0;
+  ATOMIC nsolved = 0;
+  res_T res = RES_OK;
 
   if(!estimators) {
     res = RES_BAD_ARG;
@@ -297,36 +318,52 @@ smc_solve_N
   }
 
   /* Parallel estimation of N simulations */
-  #pragma omp parallel shared(cancel)
-  {
+  logger_print(dev->logger, LOG_OUTPUT, "Solving: %3d%%\r", progress);
+  #pragma omp parallel for schedule(static, 1)
+  for(i = 0; i < (int64_t)count; ++i) {
+    size_t istep;
+    int64_t n = 0;
+    int pcent = 0;
     const int ithread = omp_get_thread_num();
-    #pragma omp for schedule(static, 1)
-    for(i = 0; i < (int64_t)count; ++i) {
-      size_t istep;
-      FOR_EACH(istep, 0, integrator->max_realisations) {
-        if(!cancel) {
-          const res_T res_local = integrator->integrand
-            (vals[ithread], dev->rngs[ithread], (unsigned)ithread,
-             (char*)ctx + (size_t)i*sizeof_ctx);
-          if(res_local == RES_OK) {
-            /* call succeded */
-            integrator->type->add
-              (estimators[i]->value, estimators[i]->value, vals[ithread]);
-            integrator->type->mul(vals[ithread], vals[ithread], vals[ithread]);
-            integrator->type->add
-              (estimators[i]->square_value, estimators[i]->square_value, vals[ithread]);
-            ++estimators[i]->nsamples;
-          } else {
-            #pragma omp atomic
-            ++estimators[i]->nfailed;
-            if (estimators[i]->nfailed > integrator->max_failures) {
-              cancel = 1;
-            }
-          }
+    res_T res_local = RES_OK;
+
+    if(ATOMIC_GET(&cancel)) continue;
+
+    FOR_EACH(istep, 0, integrator->max_realisations) {
+      if(ATOMIC_GET(&cancel)) break;
+
+      res_local = integrator->integrand
+        (vals[ithread], dev->rngs[ithread], (unsigned)ithread,
+         (char*)ctx + (size_t)i*sizeof_ctx);
+
+      if(res_local != RES_OK) {
+        ++estimators[i]->nfailed;
+        if(estimators[i]->nfailed > integrator->max_failures) {
+          ATOMIC_SET(&cancel, 1);
         }
+        break;
       }
+
+      /* call succeded */
+      integrator->type->add
+        (estimators[i]->value, estimators[i]->value, vals[ithread]);
+      integrator->type->mul
+        (vals[ithread], vals[ithread], vals[ithread]);
+      integrator->type->add
+        (estimators[i]->square_value, estimators[i]->square_value, vals[ithread]);
+      ++estimators[i]->nsamples;
+    }
+
+    n = ATOMIC_INCR(&nsolved);
+    pcent = (int)((double)n * 100.0 / (double)count + 0.5/*round*/);
+    #pragma omp critical
+    if(pcent > progress) {
+      progress = pcent;
+      logger_print(dev->logger, LOG_OUTPUT, "Solving: %3d%%\r", progress);
     }
+    
   }
+  logger_print(dev->logger, LOG_OUTPUT, "Solving: %3d%%\n", progress);
 
 exit:
   if(vals) {
diff --git a/src/test_smc_light_path.c b/src/test_smc_light_path.c
@@ -378,7 +378,7 @@ main(int argc, char** argv)
 
   integrator.integrand = light_path_integrator;
   integrator.type = &smc_float;
-  integrator.max_realisations = 4;
+  integrator.max_realisations = 64;
 
   contexts = MEM_CALLOC
     (&allocator, IMG_WIDTH*IMG_HEIGHT, sizeof(struct integrator_context));
@@ -418,7 +418,7 @@ main(int argc, char** argv)
     CHK(smc_estimator_ref_put(estimators[iestimator]) == RES_OK);
   }}
 
-  image_write_ppm_stream(&img, 0, stdout);
+  image_write_ppm_stream(&img, 0, stderr);
 
   CHK(image_release(&img) == RES_OK);
   MEM_RM(&allocator, contexts);