commit bec3cbe6c9a67e4691370a8a07e3c85b9d7ca194
parent 36cd445d4fc396de04505c2da1f060044839cf99
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Wed, 1 Jul 2020 16:48:44 +0200
Add and test the sdis_estimator_accum function
Diffstat:
5 files changed, 108 insertions(+), 0 deletions(-)
diff --git a/src/sdis.h b/src/sdis.h
@@ -936,6 +936,11 @@ sdis_estimator_get_rng_state
* from an estimator buffer */
struct ssp_rng** rng_state);
+SDIS_API res_T
+sdis_estimator_accum
+ (struct sdis_estimator* dst,
+ const struct sdis_estimator* src);
+
/*******************************************************************************
* The green function saves the estimation of the propagator
******************************************************************************/
diff --git a/src/sdis_estimator.c b/src/sdis_estimator.c
@@ -207,6 +207,74 @@ sdis_estimator_get_rng_state
return RES_OK;
}
+res_T
+sdis_estimator_accum
+ (struct sdis_estimator* dst,
+ const struct sdis_estimator* src)
+{
+ res_T res = RES_OK;
+ struct sdis_heat_path* dst_paths = NULL;
+ const struct sdis_heat_path* src_paths = NULL;
+ size_t ndst_paths = 0;
+ size_t nsrc_paths = 0;
+ size_t i;
+
+ if(!dst || !src) { res = RES_BAD_ARG; goto error; }
+
+ if(dst->type != src->type) {
+ log_err(dst->dev, "%s: estimators are mutually incompatibles.\n",
+ FUNC_NAME);
+ res = RES_OK;
+ goto error;
+ }
+
+ /* Begin the function by the copy of the registered paths since if an error
+ * occurs on dynamic array allocation, nothing has to be rewind */
+ ndst_paths = darray_heat_path_size_get(&dst->paths);
+ nsrc_paths = darray_heat_path_size_get(&src->paths);
+ res = darray_heat_path_resize(&dst->paths, ndst_paths + nsrc_paths);
+ if(res != RES_OK) {
+ log_err(dst->dev,
+ "%s: could not allocate the list of registered paths -- %s.\n",
+ FUNC_NAME, res_to_cstr(res));
+ goto error;
+ }
+ dst_paths = darray_heat_path_data_get(&dst->paths);
+ src_paths = darray_heat_path_cdata_get(&src->paths);
+ FOR_EACH(i, 0, nsrc_paths) {
+ res = heat_path_copy(&dst_paths[ndst_paths + i], &src_paths[i]);
+ ASSERT(res == RES_OK);
+ }
+
+ switch(dst->type) {
+ case SDIS_ESTIMATOR_TEMPERATURE:
+ accum_add(&dst->temperature, &dst->temperature, &src->temperature);
+ break;
+ case SDIS_ESTIMATOR_FLUX:
+ FOR_EACH(i, 0, FLUX_NAMES_COUNT__) {
+ accum_add(&dst->fluxes[i], &dst->fluxes[i], &src->fluxes[i]);
+ }
+ break;
+ default: FATAL("Unreachable code.\n"); break;
+ }
+
+ accum_add(&dst->realisation_time, &dst->realisation_time, &src->realisation_time);
+ dst->nrealisations += src->nrealisations;
+ dst->nfailures += src->nfailures;
+
+ /* Relese the rng state on the `dst' estimator since nothing can be said onto
+ * it after the accumulation */
+ if(dst->rng) {
+ SSP(rng_ref_put(dst->rng));
+ dst->rng = NULL;
+ }
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
/*******************************************************************************
* Local functions
******************************************************************************/
diff --git a/src/sdis_misc.h b/src/sdis_misc.h
@@ -39,6 +39,16 @@ static const struct accum ACCUM_NULL = ACCUM_NULL__;
#define BOLTZMANN_CONSTANT 5.6696e-8 /* W/m^2/K^4 */
+static INLINE struct accum*
+accum_add(struct accum* dst, const struct accum* a, const struct accum* b)
+{
+ ASSERT(dst && a && b);
+ dst->sum = a->sum + b->sum;
+ dst->sum2 = a->sum2 + b->sum2;
+ dst->count = a->count + b->count;
+ return dst;
+}
+
static INLINE void
sum_accums
(const struct accum accums[],
diff --git a/src/test_sdis_solve_probe.c b/src/test_sdis_solve_probe.c
@@ -249,6 +249,7 @@ main(int argc, char** argv)
{
struct mem_allocator allocator;
struct sdis_mc T = SDIS_MC_NULL;
+ struct sdis_mc T2 = SDIS_MC_NULL;
struct sdis_mc F = SDIS_MC_NULL;
struct sdis_mc time = SDIS_MC_NULL;
struct sdis_device* dev = NULL;
@@ -450,6 +451,18 @@ main(int argc, char** argv)
OK(sdis_green_function_ref_put(green));
OK(sdis_green_function_ref_put(green));
+ OK(sdis_estimator_ref_put(estimator2));
+ OK(sdis_estimator_get_rng_state(estimator, &rng_state));
+ solve_args.rng_state = rng_state;
+ OK(sdis_solve_probe(scn, &solve_args, &estimator2));
+
+ BA(sdis_estimator_accum(NULL, estimator2));
+ BA(sdis_estimator_accum(estimator, NULL));
+ OK(sdis_estimator_accum(estimator, estimator2));
+
+ OK(sdis_estimator_get_temperature(estimator, &T2));
+ CHK(eq_eps(T2.E, ref, T2.SE));
+
OK(sdis_estimator_ref_put(estimator));
OK(sdis_estimator_ref_put(estimator2));
@@ -460,6 +473,7 @@ main(int argc, char** argv)
CHK(n == 0);
OK(sdis_estimator_ref_put(estimator));
+ solve_args.rng_state = NULL;
solve_args.nrealisations = N_dump;
solve_args.register_paths = SDIS_HEAT_PATH_ALL;
diff --git a/src/test_sdis_solve_probe2.c b/src/test_sdis_solve_probe2.c
@@ -146,6 +146,7 @@ main(int argc, char** argv)
{
struct mem_allocator allocator;
struct sdis_mc T = SDIS_MC_NULL;
+ struct sdis_mc T2 = SDIS_MC_NULL;
struct sdis_mc time = SDIS_MC_NULL;
struct sdis_device* dev = NULL;
struct sdis_data* data = NULL;
@@ -269,6 +270,16 @@ main(int argc, char** argv)
check_green_function(green);
check_estimator_eq(estimator, estimator2);
+ /* Pursue the estimation */
+ OK(sdis_estimator_ref_put(estimator2));
+ OK(sdis_estimator_get_rng_state(estimator, &solve_args.rng_state));
+ solve_args.nrealisations *= 3;
+ OK(sdis_solve_probe(scn, &solve_args, &estimator2));
+ OK(sdis_estimator_accum(estimator, estimator2));
+ OK(sdis_estimator_get_temperature(estimator, &T2));
+ CHK(eq_eps(T2.E, ref, 3*T2.SE));
+ CHK(eq_eps(T2.SE, T.SE*0.5, 1e-4));
+
/* Release data */
OK(sdis_estimator_ref_put(estimator));
OK(sdis_estimator_ref_put(estimator2));
