solstice-solver

commit 6e99b59312815340851f3a4a7188ff6369ae6b73
parent 592d9279520ac5efb2a3da2ba656571da9324ceb
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Tue, 12 Sep 2017 11:22:49 +0200

BugFix: global cos was not including shadowed realisations.

Diffstat:
M
src/ssol_solver.c
 | 
15
++++++++-------
M
src/test_ssol_solver9.c
 | 
5
++---

2 files changed, 10 insertions(+), 10 deletions(-)
diff --git a/src/ssol_solver.c b/src/ssol_solver.c
@@ -930,20 +930,21 @@ trace_radiative_path
       ssol_medium_copy(&in_medium, &out_medium);
     }
 
-    /* Now that the sample ends successfully,
-     * record MC weights and register the remaining power as missing */
-    ACCUM_WEIGHT(pt.mc_samp->cos_factor, pt.cos_factor);
-    ACCUM_WEIGHT(thread_ctx->cos_factor, pt.cos_factor);
+    /* Register the remaining flux as missing */
     ACCUM_WEIGHT(thread_ctx->missing, pt.outgoing_flux);
     pt.energy_loss -= pt.outgoing_flux;
-    #undef ACCUM_WEIGHT
 
-    /* Check conservation of energy */
-    ASSERT((double)depth*DBL_EPSILON*pt.initial_flux >= fabs(pt.energy_loss));
     if(tracker) {
       path.type = hit_a_receiver ? SSOL_PATH_SUCCESS : SSOL_PATH_MISSING;
     }
   }
+  /* Now that the sample ends successfully, record MC weights */
+  ACCUM_WEIGHT(pt.mc_samp->cos_factor, pt.cos_factor);
+  ACCUM_WEIGHT(thread_ctx->cos_factor, pt.cos_factor);
+  #undef ACCUM_WEIGHT
+
+  /* Check conservation of energy at the realisation level */
+  ASSERT((double)depth*DBL_EPSILON*pt.initial_flux >= fabs(pt.energy_loss));
 
   if(tracker) {
     res = path_register_and_clear(&thread_ctx->paths, &path);
diff --git a/src/test_ssol_solver9.c b/src/test_ssol_solver9.c
@@ -141,9 +141,8 @@ main(int argc, char** argv)
 #define DNI_TGT_S (DNI * TGT_X * TGT_Y)
 #define DNI_S (DNI * SZ * SZ)
   CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK);
-  printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE);
-  printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE);
-  printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE);
+  PRINT_GLOBAL(mc_global);
+  CHECK(eq_eps(mc_global.cos_factor.E, 1./3., 3 * mc_global.cos_factor.SE), 1);
   CHECK(eq_eps(mc_global.shadowed.E, DNI_S, 3 * mc_global.shadowed.SE), 1);
   CHECK(eq_eps(mc_global.missing.E, MMAX(DNI_S, DNI_TGT_S), 
     3 * mc_global.missing.SE), 1);