commit ab90717ccb684c83a240a4617e4ecbd2371f8169 parent b9fed004322bf3939d09eac4231afc6ebefdd8b5 Author: Vincent Forest <vincent.forest@meso-star.com> Date: Tue, 14 May 2019 12:27:48 +0200 Merge branch 'feature_realisation_time' into develop Diffstat:
27 files changed, 413 insertions(+), 147 deletions(-)
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -32,12 +32,12 @@ CMAKE_DEPENDENT_OPTION(ALL_TESTS # Check dependencies ################################################################################ find_package(RCMake 0.4 REQUIRED) -find_package(Star2D 0.3 REQUIRED) +find_package(Star2D 0.3.1 REQUIRED) find_package(Star3D 0.6 REQUIRED) find_package(StarSP 0.8 REQUIRED) find_package(StarEnc 0.2.2 REQUIRED) find_package(StarEnc2D 0.2.2 REQUIRED) -find_package(RSys 0.8 REQUIRED) +find_package(RSys 0.8.1 REQUIRED) find_package(OpenMP 2.0 REQUIRED) set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${RCMAKE_SOURCE_DIR}) @@ -52,6 +52,8 @@ include_directories( ${StarEnc2D_INCLUDE_DIR} ${RSys_INCLUDE_DIR}) +set(CMAKE_C_FLAGS "-fno-omit-frame-pointer") + rcmake_append_runtime_dirs(_runtime_dirs RSys Star3D StarSP StarEnc StarEnc2D) ################################################################################ diff --git a/src/sdis.h b/src/sdis.h @@ -722,6 +722,11 @@ sdis_estimator_get_temperature struct sdis_mc* temperature); SDIS_API res_T +sdis_estimator_get_realisation_time + (const struct sdis_estimator* estimator, + struct sdis_mc* time); + +SDIS_API res_T sdis_estimator_get_convective_flux (const struct sdis_estimator* estimator, struct sdis_mc* flux); diff --git a/src/sdis_Xd_begin.h b/src/sdis_Xd_begin.h @@ -22,14 +22,6 @@ struct green_path_handle; struct sdis_heat_path; -struct accum { - double sum; /* Sum of MC weights */ - double sum2; /* Sum of square MC weights */ - size_t count; /* #accumulated MC weights */ -}; -#define ACCUM_NULL__ {0,0,0} -static const struct accum ACCUM_NULL = ACCUM_NULL__; - struct rwalk_context { struct green_path_handle* green_path; struct sdis_heat_path* heat_path; @@ -39,24 +31,6 @@ struct rwalk_context { #define RWALK_CONTEXT_NULL__ {NULL, NULL, 0, 0} static const struct rwalk_context RWALK_CONTEXT_NULL = RWALK_CONTEXT_NULL__; -static INLINE void -sum_accums - (const struct accum accums[], - const size_t naccums, - struct accum* accum) -{ - struct accum acc = ACCUM_NULL; - size_t i; - ASSERT(accums && naccums && accum); - - FOR_EACH(i, 0, naccums) { - acc.sum += accums[i].sum; - acc.sum2 += accums[i].sum2; - acc.count += accums[i].count; - } - *accum = acc; -} - #endif /* SDIS_XD_BEGIN_H */ #ifdef SDIS_XD_BEGIN_H__ diff --git a/src/sdis_estimator.c b/src/sdis_estimator.c @@ -92,6 +92,15 @@ sdis_estimator_get_temperature } res_T +sdis_estimator_get_realisation_time + (const struct sdis_estimator* estimator, struct sdis_mc* mc) +{ + if(!estimator || !mc) return RES_BAD_ARG; + *mc = estimator->realisation_time; + return RES_OK; +} + +res_T sdis_estimator_get_convective_flux (const struct sdis_estimator* estimator, struct sdis_mc* flux) { diff --git a/src/sdis_estimator_c.h b/src/sdis_estimator_c.h @@ -35,6 +35,7 @@ enum flux_name { struct sdis_estimator { struct sdis_mc temperature; + struct sdis_mc realisation_time; struct sdis_mc fluxes[FLUX_NAMES_COUNT__]; size_t nrealisations; size_t nfailures; @@ -76,19 +77,31 @@ estimator_setup_realisations_count estimator->nfailures = nrealisations - nsuccesses; } +#define SETUP_MC(McName, Sum, Sum2, Count) { \ + (McName).E = (Sum) / (double)(Count); \ + (McName).V = (Sum2) / (double)(Count) - (McName).E*(McName).E; \ + (McName).V = MMAX((McName).V, 0); \ + (McName).SE = sqrt((McName).V / (double)(Count)); \ +} (void)0 + static INLINE void estimator_setup_temperature (struct sdis_estimator* estim, const double sum, const double sum2) { - double N; ASSERT(estim && estim->nrealisations); - N = (double)estim->nrealisations; - estim->temperature.E = sum/N; - estim->temperature.V = sum2/N - estim->temperature.E*estim->temperature.E; - estim->temperature.V = MMAX(estim->temperature.V, 0); - estim->temperature.SE = sqrt(estim->temperature.V/N); + SETUP_MC(estim->temperature, sum, sum2, estim->nrealisations); +} + +static INLINE void +estimator_setup_realisation_time + (struct sdis_estimator* estim, + const double sum, + const double sum2) +{ + ASSERT(estim && estim->nrealisations); + SETUP_MC(estim->realisation_time, sum, sum2, estim->nrealisations); } static INLINE void @@ -98,14 +111,11 @@ estimator_setup_flux const double sum, const double sum2) { - double N; ASSERT(estim && (unsigned)name < FLUX_NAMES_COUNT__ && estim->nrealisations); - N = (double)estim->nrealisations; - estim->fluxes[name].E = sum/N; - estim->fluxes[name].V = sum2/N - estim->fluxes[name].E*estim->fluxes[name].E; - estim->fluxes[name].V = MMAX(estim->fluxes[name].V, 0); - estim->fluxes[name].SE = sqrt(estim->fluxes[name].V/N); + SETUP_MC(estim->fluxes[name], sum, sum2, estim->nrealisations); } +#undef SETUP_MC + #endif /* SDIS_PROBE_ESTIMATOR_C_H */ diff --git a/src/sdis_green.c b/src/sdis_green.c @@ -192,6 +192,8 @@ struct sdis_green_function { size_t npaths_valid; size_t npaths_invalid; + struct accum realisation_time; /* Time per realisation */ + struct ssp_rng_type rng_type; FILE* rng_state; @@ -489,6 +491,8 @@ sdis_green_function_solve /* Setup the estimated temperature */ estimator_setup_realisations_count(estimator, npaths, N); estimator_setup_temperature(estimator, accum, accum2); + estimator_setup_realisation_time + (estimator, green->realisation_time.sum, green->realisation_time.sum2); exit: if(rng) SSP(rng_ref_put(rng)); @@ -797,12 +801,13 @@ error: res_T green_function_finalize (struct sdis_green_function* green, - struct ssp_rng_proxy* proxy) + struct ssp_rng_proxy* proxy, + const struct accum* time) { size_t i, n; res_T res = RES_OK; - if(!green || !proxy) { + if(!green || !proxy || !time) { res = RES_BAD_ARG; goto error; } @@ -821,6 +826,9 @@ green_function_finalize } green->npaths_invalid = n - green->npaths_valid; + ASSERT(green->npaths_valid == time->count); + green->realisation_time = *time; + exit: return res; error: diff --git a/src/sdis_green.h b/src/sdis_green.h @@ -19,6 +19,7 @@ #include <rsys/rsys.h> /* Forward declaration */ +struct accum; struct sdis_green_function; struct ssp_rng_proxy; struct green_path; @@ -53,7 +54,8 @@ green_function_redux_and_clear extern LOCAL_SYM res_T green_function_finalize (struct sdis_green_function* green, - struct ssp_rng_proxy* rng_proxy); /* Proxy RNG used to estimate the function */ + struct ssp_rng_proxy* rng_proxy, /* Proxy RNG used to estimate the function */ + const struct accum* time); /* Accumulator of the realisation time */ extern LOCAL_SYM res_T green_function_create_path diff --git a/src/sdis_misc.h b/src/sdis_misc.h @@ -20,6 +20,14 @@ #include <rsys/float3.h> #include <star/ssp.h> +struct accum { + double sum; /* Sum of MC weights */ + double sum2; /* Sum of square MC weights */ + size_t count; /* #accumulated MC weights */ +}; +#define ACCUM_NULL__ {0,0,0} +static const struct accum ACCUM_NULL = ACCUM_NULL__; + /* Empirical scale factor to apply to the upper bound of the ray range in order * to handle numerical imprecisions */ #define RAY_RANGE_MAX_SCALE 1.001f @@ -31,6 +39,24 @@ #define BOLTZMANN_CONSTANT 5.6696e-8 /* W/m^2/K^4 */ +static INLINE void +sum_accums + (const struct accum accums[], + const size_t naccums, + struct accum* accum) +{ + struct accum acc = ACCUM_NULL; + size_t i; + ASSERT(accums && naccums && accum); + + FOR_EACH(i, 0, naccums) { + acc.sum += accums[i].sum; + acc.sum2 += accums[i].sum2; + acc.count += accums[i].count; + } + *accum = acc; +} + /* Reflect the V wrt the normal N. By convention V points outward the surface */ static FINLINE float* reflect_2d(float res[2], const float V[2], const float N[2]) diff --git a/src/sdis_solve_boundary_Xd.h b/src/sdis_solve_boundary_Xd.h @@ -21,6 +21,8 @@ #include "sdis_realisation.h" #include "sdis_scene_c.h" +#include <rsys/clock_time.h> +#include <rsys/rsys.h> #include <star/ssp.h> #include <omp.h> @@ -97,7 +99,8 @@ XD(solve_boundary) struct sdis_green_function** greens = NULL; struct ssp_rng_proxy* rng_proxy = NULL; struct ssp_rng** rngs = NULL; - struct accum* accums = NULL; + struct accum* acc_temps = NULL; + struct accum* acc_times = NULL; size_t i; size_t view_nprims; int64_t irealisation; @@ -185,9 +188,13 @@ XD(solve_boundary) if(res != RES_OK) goto error; } - /* Create the per thread accumulator */ - accums = MEM_CALLOC(scn->dev->allocator, scn->dev->nthreads, sizeof(*accums)); - if(!accums) { res = RES_MEM_ERR; goto error; } + /* Create the per thread accumulators */ + acc_temps = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_temps)); + if(!acc_temps) { res = RES_MEM_ERR; goto error; } + acc_times = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_times)); + if(!acc_times) { res = RES_MEM_ERR; goto error; } /* Create the per thread green function */ if(out_green) { @@ -208,9 +215,11 @@ XD(solve_boundary) omp_set_num_threads((int)scn->dev->nthreads); #pragma omp parallel for schedule(static) for(irealisation=0; irealisation<(int64_t)nrealisations; ++irealisation) { + struct time t0, t1; const int ithread = omp_get_thread_num(); struct ssp_rng* rng = rngs[ithread]; - struct accum* accum = &accums[ithread]; + struct accum* acc_temp = &acc_temps[ithread]; + struct accum* acc_time = &acc_times[ithread]; struct green_path_handle* pgreen_path = NULL; struct green_path_handle green_path = GREEN_PATH_HANDLE_NULL; struct sdis_heat_path* pheat_path = NULL; @@ -223,9 +232,13 @@ XD(solve_boundary) float st[DIM-1]; double time; res_T res_local = RES_OK; + res_T res_simul = RES_OK; if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */ + /* Begin time registration */ + time_current(&t0); + if(!out_green) { time = sample_time(rng, time_range); if(register_paths) { @@ -266,21 +279,18 @@ XD(solve_boundary) side = sides[prim.prim_id]; /* Invoke the boundary realisation */ - res_local = XD(boundary_realisation)(scn, rng, iprim, uv, time, side, + res_simul = XD(boundary_realisation)(scn, rng, iprim, uv, time, side, fp_to_meter, Tarad, Tref, pgreen_path, pheat_path, &w); - /* Update the MC accumulators */ - if(res_local == RES_OK) { - accum->sum += w; - accum->sum2 += w*w; - ++accum->count; - } else if(res_local != RES_BAD_OP) { - ATOMIC_SET(&res, res_local); + /* Fatal error */ + if(res_simul != RES_OK && res_simul != RES_BAD_OP) { + ATOMIC_SET(&res, res_simul); continue; } + /* Register heat path */ if(pheat_path) { - pheat_path->status = res_local == RES_OK + pheat_path->status = res_simul == RES_OK ? SDIS_HEAT_PATH_SUCCEED : SDIS_HEAT_PATH_FAILED; @@ -292,18 +302,36 @@ XD(solve_boundary) if(res_local != RES_OK) { ATOMIC_SET(&res, res_local); continue; } } } + + /* Stop time registration */ + time_sub(&t0, time_current(&t1), &t0); + + /* Update accumulators */ + if(res_simul == RES_OK) { + const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; + acc_temp->sum += w; acc_temp->sum2 += w*w; ++acc_temp->count; + acc_time->sum += usec; acc_time->sum2 += usec*usec; ++acc_time->count; + } } /* Setup the estimated temperature */ if(out_estimator) { - struct accum acc; - sum_accums(accums, scn->dev->nthreads, &acc); - estimator_setup_realisations_count(estimator, nrealisations, acc.count); - estimator_setup_temperature(estimator, acc.sum, acc.sum2); + struct accum acc_temp; + struct accum acc_time; + + sum_accums(acc_temps, scn->dev->nthreads, &acc_temp); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + ASSERT(acc_temp.count == acc_time.count); + + estimator_setup_realisations_count(estimator, nrealisations, acc_temp.count); + estimator_setup_temperature(estimator, acc_temp.sum, acc_temp.sum2); + estimator_setup_realisation_time(estimator, acc_time.sum, acc_time.sum2); } /* Setup the green function */ if(out_green) { + struct accum acc_time; + /* Redux the per thread green function into the green of the 1st thread */ green = greens[0]; /* Return the green of the 1st thread */ greens[0] = NULL; /* Make invalid the 1st green for 'on exit' clean up*/ @@ -311,7 +339,8 @@ XD(solve_boundary) if(res != RES_OK) goto error; /* Finalize the estimated green */ - res = green_function_finalize(green, rng_proxy); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + res = green_function_finalize(green, rng_proxy, &acc_time); if(res != RES_OK) goto error; } @@ -328,7 +357,8 @@ exit: } MEM_RM(scn->dev->allocator, greens); } - if(accums) MEM_RM(scn->dev->allocator, accums); + if(acc_temps) MEM_RM(scn->dev->allocator, acc_temps); + if(acc_times) MEM_RM(scn->dev->allocator, acc_times); if(scene) SXD(scene_ref_put(scene)); if(shape) SXD(shape_ref_put(shape)); if(view) SXD(scene_view_ref_put(view)); @@ -369,6 +399,7 @@ XD(solve_boundary_flux) struct ssp_rng_proxy* rng_proxy = NULL; struct ssp_rng** rngs = NULL; struct accum* acc_tp = NULL; /* Per thread temperature accumulator */ + struct accum* acc_ti = NULL; /* Per thread realisation time accumulator */ struct accum* acc_fl = NULL; /* Per thread flux accumulator */ struct accum* acc_fc = NULL; /* Per thread convective flux accumulator */ struct accum* acc_fr = NULL; /* Per thread radiative flux accumulator */ @@ -458,6 +489,7 @@ XD(solve_boundary_flux) if(!Dst) { res = RES_MEM_ERR; goto error; } \ } (void)0 ALLOC_ACCUMS(acc_tp); + ALLOC_ACCUMS(acc_ti); ALLOC_ACCUMS(acc_fc); ALLOC_ACCUMS(acc_fl); ALLOC_ACCUMS(acc_fr); @@ -470,9 +502,11 @@ XD(solve_boundary_flux) omp_set_num_threads((int)scn->dev->nthreads); #pragma omp parallel for schedule(static) for(irealisation = 0; irealisation < (int64_t)nrealisations; ++irealisation) { + struct time t0, t1; const int ithread = omp_get_thread_num(); struct ssp_rng* rng = rngs[ithread]; struct accum* acc_temp = &acc_tp[ithread]; + struct accum* acc_time = &acc_ti[ithread]; struct accum* acc_flux = &acc_fl[ithread]; struct accum* acc_fcon = &acc_fc[ithread]; struct accum* acc_frad = &acc_fr[ithread]; @@ -489,9 +523,13 @@ XD(solve_boundary_flux) double time; int flux_mask = 0; res_T res_local = RES_OK; + res_T res_simul = RES_OK; if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */ + /* Begin time registration */ + time_current(&t0); + time = sample_time(rng, time_range); /* Sample a position onto the boundary */ @@ -545,9 +583,18 @@ XD(solve_boundary_flux) flux_mask = 0; if(hr > 0) flux_mask |= FLUX_FLAG_RADIATIVE; if(hc > 0) flux_mask |= FLUX_FLAG_CONVECTIVE; - res_local = XD(boundary_flux_realisation)(scn, rng, iprim, uv, time, + + res_simul = XD(boundary_flux_realisation)(scn, rng, iprim, uv, time, solid_side, fp_to_meter, Tarad, Tref, flux_mask, T_brf); - if(res_local == RES_OK) { + + /* Stop time registration */ + time_sub(&t0, time_current(&t1), &t0); + + if(res_simul != RES_OK && res_simul != RES_BAD_OP) { + ATOMIC_SET(&res, res_simul); + continue; + } else if(res_simul == RES_OK) { /* Update accumulators */ + const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; const double Tboundary = T_brf[0]; const double Tradiative = T_brf[1]; const double Tfluid = T_brf[2]; @@ -558,6 +605,10 @@ XD(solve_boundary_flux) acc_temp->sum += Tboundary; acc_temp->sum2 += Tboundary*Tboundary; ++acc_temp->count; + /* Time */ + acc_time->sum += usec; + acc_time->sum2 += usec*usec; + ++acc_time->count; /* Overwall flux */ acc_flux->sum += w_total; acc_flux->sum2 += w_total*w_total; @@ -570,25 +621,25 @@ XD(solve_boundary_flux) acc_frad->sum += w_rad; acc_frad->sum2 += w_rad*w_rad; ++acc_frad->count; - } else if(res_local != RES_BAD_OP) { - ATOMIC_SET(&res, res_local); - continue; } } if(res != RES_OK) goto error; /* Redux the per thread accumulators */ sum_accums(acc_tp, scn->dev->nthreads, &acc_tp[0]); + sum_accums(acc_ti, scn->dev->nthreads, &acc_ti[0]); sum_accums(acc_fc, scn->dev->nthreads, &acc_fc[0]); sum_accums(acc_fr, scn->dev->nthreads, &acc_fr[0]); sum_accums(acc_fl, scn->dev->nthreads, &acc_fl[0]); ASSERT(acc_tp[0].count == acc_fl[0].count); + ASSERT(acc_tp[0].count == acc_ti[0].count); ASSERT(acc_tp[0].count == acc_fr[0].count); ASSERT(acc_tp[0].count == acc_fc[0].count); /* Setup the estimated values */ estimator_setup_realisations_count(estimator, nrealisations, acc_tp[0].count); estimator_setup_temperature(estimator, acc_tp[0].sum, acc_tp[0].sum2); + estimator_setup_realisation_time(estimator, acc_ti[0].sum, acc_ti[0].sum2); estimator_setup_flux(estimator, FLUX_CONVECTIVE, acc_fc[0].sum, acc_fc[0].sum2); estimator_setup_flux(estimator, FLUX_RADIATIVE, acc_fr[0].sum, acc_fr[0].sum2); estimator_setup_flux(estimator, FLUX_TOTAL, acc_fl[0].sum, acc_fl[0].sum2); @@ -599,6 +650,7 @@ exit: MEM_RM(scn->dev->allocator, rngs); } if(acc_tp) MEM_RM(scn->dev->allocator, acc_tp); + if(acc_ti) MEM_RM(scn->dev->allocator, acc_ti); if(acc_fc) MEM_RM(scn->dev->allocator, acc_fc); if(acc_fr) MEM_RM(scn->dev->allocator, acc_fr); if(acc_fl) MEM_RM(scn->dev->allocator, acc_fl); diff --git a/src/sdis_solve_medium_Xd.h b/src/sdis_solve_medium_Xd.h @@ -20,6 +20,7 @@ #include "sdis_scene_c.h" #include <rsys/algorithm.h> +#include <rsys/clock_time.h> #include <rsys/dynamic_array.h> #include "sdis_Xd_begin.h" @@ -215,7 +216,8 @@ XD(solve_medium) struct ssp_rng_proxy* rng_proxy = NULL; struct ssp_rng** rngs = NULL; struct sdis_estimator* estimator = NULL; - struct accum* accums = NULL; + struct accum* acc_temps = NULL; + struct accum* acc_times = NULL; int64_t irealisation; int cumul_is_init = 0; size_t i; @@ -257,9 +259,13 @@ XD(solve_medium) if(res != RES_OK) goto error; } - /* Create the per thread accumulator */ - accums = MEM_CALLOC(scn->dev->allocator, scn->dev->nthreads, sizeof(*accums)); - if(!accums) { res = RES_MEM_ERR; goto error; } + /* Create the per thread accumulators */ + acc_temps = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_temps)); + if(!acc_temps) { res = RES_MEM_ERR; goto error; } + acc_times = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_times)); + if(!acc_times) { res = RES_MEM_ERR; goto error; } /* Compute the enclosure cumulative */ darray_enclosure_cumul_init(scn->dev->allocator, &cumul); @@ -286,9 +292,11 @@ XD(solve_medium) omp_set_num_threads((int)scn->dev->nthreads); #pragma omp parallel for schedule(static) for(irealisation = 0; irealisation < (int64_t)nrealisations; ++irealisation) { + struct time t0, t1; const int ithread = omp_get_thread_num(); struct ssp_rng* rng = rngs[ithread]; - struct accum* accum = accums + ithread; + struct accum* acc_temp = &acc_temps[ithread]; + struct accum* acc_time = &acc_times[ithread]; struct green_path_handle* pgreen_path = NULL; struct green_path_handle green_path = GREEN_PATH_HANDLE_NULL; const struct enclosure* enc = NULL; @@ -298,9 +306,12 @@ XD(solve_medium) double time; double pos[DIM]; res_T res_local = RES_OK; + res_T res_simul = RES_OK; if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */ + time_current(&t0); + if(!out_green) { /* Sample the time */ time = sample_time(rng, time_range); @@ -331,19 +342,17 @@ XD(solve_medium) } /* Run a probe realisation */ - res_local = XD(probe_realisation)((size_t)irealisation, scn, rng, mdm, pos, + res_simul = XD(probe_realisation)((size_t)irealisation, scn, rng, mdm, pos, time, fp_to_meter, Tarad, Tref, pgreen_path, pheat_path, &weight); - if(res_local != RES_OK) { - if(res_local != RES_BAD_OP) { ATOMIC_SET(&res, res_local); continue; } - } else { - accum->sum += weight; - accum->sum2 += weight*weight; - ++accum->count; + + if(res_simul != RES_OK && res_simul != RES_BAD_OP) { + ATOMIC_SET(&res, res_simul); + continue; } /* Finalize the registered path */ if(pheat_path) { - pheat_path->status = res_local == RES_OK + pheat_path->status = res_simul == RES_OK ? SDIS_HEAT_PATH_SUCCEED : SDIS_HEAT_PATH_FAILED; @@ -355,18 +364,36 @@ XD(solve_medium) if(res_local != RES_OK) { ATOMIC_SET(&res, res_local); continue; } } } + + /* Stop time registration */ + time_sub(&t0, time_current(&t1), &t0); + + /* Update accumulators */ + if(res_simul == RES_OK) { + const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; + acc_temp->sum += weight; acc_temp->sum2 += weight*weight; ++acc_temp->count; + acc_time->sum += usec; acc_time->sum2 += usec*usec; ++acc_time->count; + } } if(res != RES_OK) goto error; /* Setup the estimated temperature */ if(out_estimator) { - struct accum acc; - sum_accums(accums, scn->dev->nthreads, &acc); - estimator_setup_realisations_count(estimator, nrealisations, acc.count); - estimator_setup_temperature(estimator, acc.sum, acc.sum2); + struct accum acc_temp; + struct accum acc_time; + + sum_accums(acc_temps, scn->dev->nthreads, &acc_temp); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + ASSERT(acc_temp.count == acc_time.count); + + estimator_setup_realisations_count(estimator, nrealisations, acc_temp.count); + estimator_setup_temperature(estimator, acc_temp.sum, acc_temp.sum2); + estimator_setup_realisation_time(estimator, acc_time.sum, acc_time.sum2); } if(out_green) { + struct accum acc_time; + /* Redux the per thread green function into the green of the 1st thread */ green = greens[0]; /* Return the green of the 1st thread */ greens[0] = NULL; /* Make invalid the 1st green for 'on exit' clean up*/ @@ -374,7 +401,8 @@ XD(solve_medium) if(res != RES_OK) goto error; /* Finalize the estimated green */ - res = green_function_finalize(green, rng_proxy); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + res = green_function_finalize(green, rng_proxy, &acc_time); if(res != RES_OK) goto error; } @@ -391,7 +419,8 @@ exit: } MEM_RM(scn->dev->allocator, greens); } - if(accums) MEM_RM(scn->dev->allocator, accums); + if(acc_temps) MEM_RM(scn->dev->allocator, acc_temps); + if(acc_times) MEM_RM(scn->dev->allocator, acc_times); if(cumul_is_init) darray_enclosure_cumul_release(&cumul); if(rng_proxy) SSP(rng_proxy_ref_put(rng_proxy)); if(out_estimator) *out_estimator = estimator; diff --git a/src/sdis_solve_probe_Xd.h b/src/sdis_solve_probe_Xd.h @@ -20,6 +20,7 @@ #include "sdis_realisation.h" #include "sdis_scene_c.h" +#include <rsys/clock_time.h> #include <star/ssp.h> #include <omp.h> @@ -47,7 +48,8 @@ XD(solve_probe) struct sdis_green_function** greens = NULL; struct ssp_rng_proxy* rng_proxy = NULL; struct ssp_rng** rngs = NULL; - struct accum* accums = NULL; + struct accum* acc_temps = NULL; + struct accum* acc_times = NULL; int64_t irealisation = 0; size_t i; ATOMIC res = RES_OK; @@ -89,8 +91,12 @@ XD(solve_probe) } /* Create the per thread accumulator */ - accums = MEM_CALLOC(scn->dev->allocator, scn->dev->nthreads, sizeof(*accums)); - if(!accums) { res = RES_MEM_ERR; goto error; } + acc_temps = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_temps)); + if(!acc_temps) { res = RES_MEM_ERR; goto error; } + acc_times = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_times)); + if(!acc_times) { res = RES_MEM_ERR; goto error; } /* Retrieve the medium in which the submitted position lies */ res = scene_get_medium(scn, position, NULL, &medium); @@ -116,19 +122,25 @@ XD(solve_probe) omp_set_num_threads((int)scn->dev->nthreads); #pragma omp parallel for schedule(static) for(irealisation = 0; irealisation < (int64_t)nrealisations; ++irealisation) { + struct time t0, t1; const int ithread = omp_get_thread_num(); struct ssp_rng* rng = rngs[ithread]; - struct accum* accum = &accums[ithread]; + struct accum* acc_temp = &acc_temps[ithread]; + struct accum* acc_time = &acc_times[ithread]; struct green_path_handle* pgreen_path = NULL; struct green_path_handle green_path = GREEN_PATH_HANDLE_NULL; struct sdis_heat_path* pheat_path = NULL; struct sdis_heat_path heat_path; double w = NaN; double time; - res_T res_local; + res_T res_local = RES_OK; + res_T res_simul = RES_OK; if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */ + /* Begin time registration */ + time_current(&t0); + if(!out_green) { time = sample_time(rng, time_range); if(register_paths) { @@ -145,19 +157,17 @@ XD(solve_probe) pgreen_path = &green_path; } - res_local = XD(probe_realisation)((size_t)irealisation, scn, rng, medium, + res_simul = XD(probe_realisation)((size_t)irealisation, scn, rng, medium, position, time, fp_to_meter, Tarad, Tref, pgreen_path, pheat_path, &w); - if(res_local == RES_OK) { - accum->sum += w; - accum->sum2 += w*w; - ++accum->count; - } else if(res_local != RES_BAD_OP) { - ATOMIC_SET(&res, res_local); + + /* Handle fatal error */ + if(res_simul != RES_OK && res_simul != RES_BAD_OP) { + ATOMIC_SET(&res, res_simul); continue; } if(pheat_path) { - pheat_path->status = res_local == RES_OK + pheat_path->status = res_simul == RES_OK ? SDIS_HEAT_PATH_SUCCEED : SDIS_HEAT_PATH_FAILED; @@ -169,18 +179,36 @@ XD(solve_probe) if(res_local != RES_OK) { ATOMIC_SET(&res, res_local); continue; } } } + + /* Stop time registration */ + time_sub(&t0, time_current(&t1), &t0); + + /* Update accumulators */ + if(res_simul == RES_OK) { + const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; + acc_temp->sum += w; acc_temp->sum2 += w*w; ++acc_temp->count; + acc_time->sum += usec; acc_time->sum2 += usec*usec; ++acc_time->count; + } } if(res != RES_OK) goto error; - /* Setup the estimated temperature */ + /* Setup the estimated temperature and per realisation time */ if(out_estimator) { - struct accum acc; - sum_accums(accums, scn->dev->nthreads, &acc); - estimator_setup_realisations_count(estimator, nrealisations, acc.count); - estimator_setup_temperature(estimator, acc.sum, acc.sum2); + struct accum acc_temp; + struct accum acc_time; + + sum_accums(acc_temps, scn->dev->nthreads, &acc_temp); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + ASSERT(acc_temp.count == acc_time.count); + + estimator_setup_realisations_count(estimator, nrealisations, acc_temp.count); + estimator_setup_temperature(estimator, acc_temp.sum, acc_temp.sum2); + estimator_setup_realisation_time(estimator, acc_time.sum, acc_time.sum2); } if(out_green) { + struct accum acc_time; + /* Redux the per thread green function into the green of the 1st thread */ green = greens[0]; /* Return the green of the 1st thread */ greens[0] = NULL; /* Make invalid the 1st green for 'on exit' clean up*/ @@ -188,7 +216,8 @@ XD(solve_probe) if(res != RES_OK) goto error; /* Finalize the estimated green */ - res = green_function_finalize(green, rng_proxy); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + res = green_function_finalize(green, rng_proxy, &acc_time); if(res != RES_OK) goto error; } @@ -205,7 +234,8 @@ exit: } MEM_RM(scn->dev->allocator, greens); } - if(accums) MEM_RM(scn->dev->allocator, accums); + if(acc_temps) MEM_RM(scn->dev->allocator, acc_temps); + if(acc_times) MEM_RM(scn->dev->allocator, acc_times); if(rng_proxy) SSP(rng_proxy_ref_put(rng_proxy)); if(out_green) *out_green = green; if(out_estimator) *out_estimator = estimator; diff --git a/src/sdis_solve_probe_boundary_Xd.h b/src/sdis_solve_probe_boundary_Xd.h @@ -20,6 +20,7 @@ #include "sdis_realisation.h" #include "sdis_scene_c.h" +#include <rsys/clock_time.h> #include <star/ssp.h> #include <omp.h> @@ -48,7 +49,8 @@ XD(solve_probe_boundary) struct sdis_green_function** greens = NULL; struct ssp_rng_proxy* rng_proxy = NULL; struct ssp_rng** rngs = NULL; - struct accum* accums = NULL; + struct accum* acc_temps = NULL; + struct accum* acc_times = NULL; int64_t irealisation = 0; size_t i; ATOMIC res = RES_OK; @@ -131,8 +133,12 @@ XD(solve_probe_boundary) } /* Create the per thread accumulator */ - accums = MEM_CALLOC(scn->dev->allocator, scn->dev->nthreads, sizeof(*accums)); - if(!accums) { res = RES_MEM_ERR; goto error; } + acc_temps = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_temps)); + if(!acc_temps) { res = RES_MEM_ERR; goto error; } + acc_times = MEM_CALLOC + (scn->dev->allocator, scn->dev->nthreads, sizeof(*acc_times)); + if(!acc_times) { res = RES_MEM_ERR; goto error; } /* Create the per thread green function */ if(out_green) { @@ -154,19 +160,25 @@ XD(solve_probe_boundary) omp_set_num_threads((int)scn->dev->nthreads); #pragma omp parallel for schedule(static) for(irealisation = 0; irealisation < (int64_t)nrealisations; ++irealisation) { + struct time t0, t1; const int ithread = omp_get_thread_num(); struct ssp_rng* rng = rngs[ithread]; - struct accum* accum = &accums[ithread]; + struct accum* acc_temp = &acc_temps[ithread]; + struct accum* acc_time = &acc_times[ithread]; struct green_path_handle* pgreen_path = NULL; struct green_path_handle green_path = GREEN_PATH_HANDLE_NULL; struct sdis_heat_path* pheat_path = NULL; struct sdis_heat_path heat_path; double w = NaN; double time; - res_T res_local; + res_T res_local = RES_OK; + res_T res_simul = RES_OK; if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */ + /* Begin time registration */ + time_current(&t0); + if(!out_green) { time = sample_time(rng, time_range); if(register_paths) { @@ -182,19 +194,17 @@ XD(solve_probe_boundary) pgreen_path = &green_path; } - res_local = XD(boundary_realisation)(scn, rng, iprim, uv, time, side, + res_simul = XD(boundary_realisation)(scn, rng, iprim, uv, time, side, fp_to_meter, Tarad, Tref, pgreen_path, pheat_path, &w); - if(res_local == RES_OK) { - accum->sum += w; - accum->sum2 += w*w; - ++accum->count; - } else if(res_local != RES_BAD_OP) { - ATOMIC_SET(&res, res_local); + + /* Handle fatal error */ + if(res_simul != RES_OK && res_simul != RES_BAD_OP) { + ATOMIC_SET(&res, res_simul); continue; } if(pheat_path) { - pheat_path->status = res_local == RES_OK + pheat_path->status = res_simul == RES_OK ? SDIS_HEAT_PATH_SUCCEED : SDIS_HEAT_PATH_FAILED; @@ -206,18 +216,36 @@ XD(solve_probe_boundary) if(res_local != RES_OK) { ATOMIC_SET(&res, res_local); continue; } } } + + /* Stop time registration */ + time_sub(&t0, time_current(&t1), &t0); + + /* Update accumulators */ + if(res_simul == RES_OK) { + const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; + acc_temp->sum += w; acc_temp->sum2 += w*w; ++acc_temp->count; + acc_time->sum += usec; acc_time->sum2 += usec*usec; ++acc_time->count; + } } if(res != RES_OK) goto error; - /* Setup the estimated temperature */ + /* Setup the estimated temperature and per realisation time */ if(out_estimator) { - struct accum acc; - sum_accums(accums, scn->dev->nthreads, &acc); - estimator_setup_realisations_count(estimator, nrealisations, acc.count); - estimator_setup_temperature(estimator, acc.sum, acc.sum2); + struct accum acc_temp; + struct accum acc_time; + + sum_accums(acc_temps, scn->dev->nthreads, &acc_temp); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + ASSERT(acc_temp.count == acc_time.count); + + estimator_setup_realisations_count(estimator, nrealisations, acc_temp.count); + estimator_setup_temperature(estimator, acc_temp.sum, acc_temp.sum2); + estimator_setup_realisation_time(estimator, acc_time.sum, acc_time.sum2); } if(out_green) { + struct accum acc_time; + /* Redux the per thread green function into the green of the 1st thread */ green = greens[0]; /* Return the green of the 1st thread */ greens[0] = NULL; /* Make invalid the 1st green for 'on exit' clean up*/ @@ -225,7 +253,8 @@ XD(solve_probe_boundary) if(res != RES_OK) goto error; /* Finalize the estimated green */ - res = green_function_finalize(green, rng_proxy); + sum_accums(acc_times, scn->dev->nthreads, &acc_time); + res = green_function_finalize(green, rng_proxy, &acc_time); if(res != RES_OK) goto error; } @@ -242,7 +271,8 @@ exit: } MEM_RM(scn->dev->allocator, greens); } - if(accums) MEM_RM(scn->dev->allocator, accums); + if(acc_temps) MEM_RM(scn->dev->allocator, acc_temps); + if(acc_times) MEM_RM(scn->dev->allocator, acc_times); if(rng_proxy) SSP(rng_proxy_ref_put(rng_proxy)); if(out_green) *out_green = green; if(out_estimator) *out_estimator = estimator; @@ -279,6 +309,7 @@ XD(solve_probe_boundary_flux) enum sdis_side solid_side, fluid_side; struct sdis_interface_fragment frag; struct accum* acc_tp = NULL; /* Per thread temperature accumulator */ + struct accum* acc_ti = NULL; /* Per thread realisation time */ struct accum* acc_fl = NULL; /* Per thread flux accumulator */ struct accum* acc_fc = NULL; /* Per thread convective flux accumulator */ struct accum* acc_fr = NULL; /* Per thread radiative flux accumulator */ @@ -369,6 +400,7 @@ XD(solve_probe_boundary_flux) if(!Dst) { res = RES_MEM_ERR; goto error; } \ } (void)0 ALLOC_ACCUMS(acc_tp); + ALLOC_ACCUMS(acc_ti); ALLOC_ACCUMS(acc_fc); ALLOC_ACCUMS(acc_fl); ALLOC_ACCUMS(acc_fr); @@ -387,19 +419,24 @@ XD(solve_probe_boundary_flux) omp_set_num_threads((int)scn->dev->nthreads); #pragma omp parallel for schedule(static) for(irealisation = 0; irealisation < (int64_t)nrealisations; ++irealisation) { + struct time t0, t1; const int ithread = omp_get_thread_num(); struct ssp_rng* rng = rngs[ithread]; struct accum* acc_temp = &acc_tp[ithread]; + struct accum* acc_time = &acc_ti[ithread]; struct accum* acc_flux = &acc_fl[ithread]; struct accum* acc_fcon = &acc_fc[ithread]; struct accum* acc_frad = &acc_fr[ithread]; double time, epsilon, hc, hr; int flux_mask = 0; double T_brf[3] = { 0, 0, 0 }; - res_T res_local; + res_T res_simul = RES_OK; if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */ + /* Begin time registration */ + time_current(&t0); + time = sample_time(rng, time_range); /* Compute hr and hc */ @@ -412,9 +449,17 @@ XD(solve_probe_boundary_flux) flux_mask = 0; if(hr > 0) flux_mask |= FLUX_FLAG_RADIATIVE; if(hc > 0) flux_mask |= FLUX_FLAG_CONVECTIVE; - res_local = XD(boundary_flux_realisation)(scn, rng, iprim, uv, time, + res_simul = XD(boundary_flux_realisation)(scn, rng, iprim, uv, time, solid_side, fp_to_meter, Tarad, Tref, flux_mask, T_brf); - if(res_local == RES_OK) { + + /* Stop time registration */ + time_sub(&t0, time_current(&t1), &t0); + + if(res_simul != RES_OK && res_simul != RES_BAD_OP) { + ATOMIC_SET(&res, res_simul); + continue; + } else if(res_simul == RES_OK) { /* Update accumulators */ + const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; const double Tboundary = T_brf[0]; const double Tradiative = T_brf[1]; const double Tfluid = T_brf[2]; @@ -425,6 +470,10 @@ XD(solve_probe_boundary_flux) acc_temp->sum += Tboundary; acc_temp->sum2 += Tboundary*Tboundary; ++acc_temp->count; + /* Time */ + acc_time->sum += usec; + acc_time->sum2 += usec*usec; + ++acc_time->count; /* Overwall flux */ acc_flux->sum += w_total; acc_flux->sum2 += w_total*w_total; @@ -437,25 +486,25 @@ XD(solve_probe_boundary_flux) acc_frad->sum += w_rad; acc_frad->sum2 += w_rad*w_rad; ++acc_frad->count; - } else if(res_local != RES_BAD_OP) { - ATOMIC_SET(&res, res_local); - continue; } } if(res != RES_OK) goto error; /* Redux the per thread accumulators */ sum_accums(acc_tp, scn->dev->nthreads, &acc_tp[0]); + sum_accums(acc_ti, scn->dev->nthreads, &acc_ti[0]); sum_accums(acc_fc, scn->dev->nthreads, &acc_fc[0]); sum_accums(acc_fr, scn->dev->nthreads, &acc_fr[0]); sum_accums(acc_fl, scn->dev->nthreads, &acc_fl[0]); ASSERT(acc_tp[0].count == acc_fl[0].count); + ASSERT(acc_tp[0].count == acc_ti[0].count); ASSERT(acc_tp[0].count == acc_fr[0].count); ASSERT(acc_tp[0].count == acc_fc[0].count); /* Setup the estimated values */ estimator_setup_realisations_count(estimator, nrealisations, acc_tp[0].count); estimator_setup_temperature(estimator, acc_tp[0].sum, acc_tp[0].sum2); + estimator_setup_realisation_time(estimator, acc_ti[0].sum, acc_ti[0].sum2); estimator_setup_flux(estimator, FLUX_CONVECTIVE, acc_fc[0].sum, acc_fc[0].sum2); estimator_setup_flux(estimator, FLUX_RADIATIVE, acc_fr[0].sum, acc_fr[0].sum2); estimator_setup_flux(estimator, FLUX_TOTAL, acc_fl[0].sum, acc_fl[0].sum2); @@ -466,6 +515,7 @@ exit: MEM_RM(scn->dev->allocator, rngs); } if(acc_tp) MEM_RM(scn->dev->allocator, acc_tp); + if(acc_ti) MEM_RM(scn->dev->allocator, acc_ti); if(acc_fc) MEM_RM(scn->dev->allocator, acc_fc); if(acc_fr) MEM_RM(scn->dev->allocator, acc_fr); if(acc_fl) MEM_RM(scn->dev->allocator, acc_fl); diff --git a/src/test_sdis_conducto_radiative.c b/src/test_sdis_conducto_radiative.c @@ -384,6 +384,7 @@ main(int argc, char** argv) OK(ssp_rng_create(&allocator, &ssp_rng_kiss, &rng)); FOR_EACH(isimul, 0, nsimuls) { struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_estimator* estimator; struct sdis_estimator* estimator2; struct sdis_green_function* green; @@ -402,11 +403,13 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); u = (pos[0] + 1) / thickness; ref = u * Ts1 + (1-u) * Ts0; printf("Temperature at (%g, %g, %g) = %g ~ %g +/- %g\n", SPLIT3(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(nfails + nreals == N); @@ -426,6 +429,8 @@ main(int argc, char** argv) OK(sdis_solve_probe(scn, 10, pos, time_range, 1, -1, Tref, SDIS_HEAT_PATH_ALL, &estimator)); OK(sdis_estimator_ref_put(estimator)); + + printf("\n"); } /* Release memory */ diff --git a/src/test_sdis_conducto_radiative_2d.c b/src/test_sdis_conducto_radiative_2d.c @@ -390,6 +390,7 @@ main(int argc, char** argv) OK(ssp_rng_create(&allocator, &ssp_rng_kiss, &rng)); FOR_EACH(isimul, 0, nsimuls) { struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_estimator* estimator; struct sdis_estimator* estimator2; struct sdis_green_function* green; @@ -407,11 +408,13 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); u = (pos[0] + 1) / thickness; ref = u * Ts1 + (1-u) * Ts0; printf("Temperature at (%g, %g) = %g ~ %g +/- %g\n", SPLIT2(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(nfails + nreals == N); @@ -431,6 +434,8 @@ main(int argc, char** argv) OK(sdis_solve_probe (scn, 10, pos, time_range, 1, -1, Tref, SDIS_HEAT_PATH_ALL, &estimator)); OK(sdis_estimator_ref_put(estimator)); + + printf("\n"); } /* Release memory */ diff --git a/src/test_sdis_convection.c b/src/test_sdis_convection.c @@ -167,6 +167,7 @@ main(int argc, char** argv) { struct mem_allocator allocator; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc mc_time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_medium* fluid = NULL; struct sdis_medium* solid = NULL; @@ -266,7 +267,7 @@ main(int argc, char** argv) /* Test in 3D for various time values. */ nu = (6 * H) / (RHO*CP); Tinf = (H*(T0 + T1 + T2 + T3 + T4 + T5)) / (6 * H); - printf("Temperature of the box at (%g %g %g)\n", SPLIT3(pos)); + printf(">>> Temperature of the box at (%g %g %g)\n\n", SPLIT3(pos)); FOR_EACH(i, 0, 5) { double time = i ? (double)i / nu : INF; double time_range[2]; @@ -279,10 +280,12 @@ main(int argc, char** argv) /* Solve in 3D */ OK(sdis_solve_probe(box_scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &mc_time)); OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); CHK(nfails + nreals == N); - printf(" t=%g : %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + printf("Temperature at %g = %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", mc_time.E, mc_time.SE); if(nfails) printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(eq_eps(T.E, ref, T.SE * 3)); @@ -297,12 +300,13 @@ main(int argc, char** argv) } OK(sdis_estimator_ref_put(estimator)); + printf("\n"); } /* Test in 2D for various time values. */ nu = (4 * H) / (RHO*CP); Tinf = (H * (T0 + T1 + T2 + T3)) / (4 * H); - printf("Temperature of the square at (%g %g)\n", SPLIT2(pos)); + printf(">>> Temperature of the square at (%g %g)\n\n", SPLIT2(pos)); FOR_EACH(i, 0, 5) { double time = i ? (double)i / nu : INF; double time_range[2]; @@ -318,7 +322,9 @@ main(int argc, char** argv) OK(sdis_estimator_get_failure_count(estimator, &nfails)); CHK(nfails + nreals == N); OK(sdis_estimator_get_temperature(estimator, &T)); - printf(" t=%g : %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + OK(sdis_estimator_get_realisation_time(estimator, &mc_time)); + printf("Temperature at %g = %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", mc_time.E, mc_time.SE); if(nfails) printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(eq_eps(T.E, ref, T.SE * 3)); @@ -333,6 +339,7 @@ main(int argc, char** argv) } OK(sdis_estimator_ref_put(estimator)); + printf("\n"); } OK(sdis_scene_ref_put(box_scn)); diff --git a/src/test_sdis_convection_non_uniform.c b/src/test_sdis_convection_non_uniform.c @@ -177,6 +177,7 @@ main(int argc, char** argv) { struct mem_allocator allocator; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc mc_time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_medium* fluid = NULL; struct sdis_medium* solid = NULL; @@ -283,7 +284,7 @@ main(int argc, char** argv) nu = (HC0 + HC1 + HC2 + HC3 + HC4 + HC5) / (RHO * CP); Tinf = (HC0 * T0 + HC1 * T1 + HC2 * T2 + HC3 * T3 + HC4 * T4 + HC5 * T5) / (HC0 + HC1 + HC2 + HC3 + HC4 + HC5); - printf("Temperature of the box at (%g %g %g)\n", SPLIT3(pos)); + printf(">>> Temperature of the box at (%g %g %g)\n\n", SPLIT3(pos)); FOR_EACH(i, 0, 5) { double time = i ? (double)i / nu : INF; double time_range[2]; @@ -298,7 +299,9 @@ main(int argc, char** argv) OK(sdis_estimator_get_failure_count(estimator, &nfails)); CHK(nfails + nreals == N); OK(sdis_estimator_get_temperature(estimator, &T)); - printf(" t=%g : %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + OK(sdis_estimator_get_realisation_time(estimator, &mc_time)); + printf("Temperature at %g = %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", mc_time.E, mc_time.SE); if(nfails) printf("#failures = %lu/%lu\n", (unsigned long)nfails,(unsigned long)N); CHK(eq_eps(T.E, ref, T.SE * 3)); @@ -313,12 +316,13 @@ main(int argc, char** argv) } OK(sdis_estimator_ref_put(estimator)); + printf("\n"); } /* Test in 2D for various time values. */ nu = (HC0 + HC1 + HC2 + HC3) / (RHO * CP); Tinf = (HC0 * T0 + HC1 * T1 + HC2 * T2 + HC3 * T3) / (HC0 + HC1 + HC2 + HC3); - printf("Temperature of the square at (%g %g)\n", SPLIT2(pos)); + printf(">>> Temperature of the square at (%g %g)\n\n", SPLIT2(pos)); FOR_EACH(i, 0, 5) { double time = i ? (double)i / nu : INF; double time_range[2]; @@ -332,7 +336,9 @@ main(int argc, char** argv) OK(sdis_estimator_get_failure_count(estimator, &nfails)); CHK(nfails + nreals == N); OK(sdis_estimator_get_temperature(estimator, &T)); - printf(" t=%g : %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + OK(sdis_estimator_get_realisation_time(estimator, &mc_time)); + printf("Temperature at %g = %g ~ %g +/- %g\n", time, ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", mc_time.E, mc_time.SE); if(nfails) printf("#failures = %lu/%lu\n", (unsigned long)nfails,(unsigned long)N); CHK(eq_eps(T.E, ref, T.SE * 3)); @@ -347,6 +353,7 @@ main(int argc, char** argv) } OK(sdis_estimator_ref_put(estimator)); + printf("\n"); } OK(sdis_scene_ref_put(box_scn)); diff --git a/src/test_sdis_flux.c b/src/test_sdis_flux.c @@ -135,6 +135,7 @@ solve(struct sdis_scene* scn, const double pos[]) struct sdis_estimator* estimator2; struct sdis_green_function* green; struct sdis_mc T; + struct sdis_mc time; size_t nreals; size_t nfails; double ref; @@ -152,6 +153,7 @@ solve(struct sdis_scene* scn, const double pos[]) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); OK(sdis_scene_get_dimension(scn, &dim)); @@ -166,6 +168,7 @@ solve(struct sdis_scene* scn, const double pos[]) break; default: FATAL("Unreachable code.\n"); break; } + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); printf("Elapsed time = %s\n\n", dump); diff --git a/src/test_sdis_solve_boundary.c b/src/test_sdis_solve_boundary.c @@ -148,15 +148,18 @@ check_estimator const double ref) { struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; size_t nreals; size_t nfails; CHK(estimator && nrealisations); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); printf("%g ~ %g +/- %g\n", ref, T.E, T.SE); - printf("#failures = %lu/%lu\n", + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); + printf("#failures = %lu/%lu\n\n", (unsigned long)nfails, (unsigned long)nrealisations); CHK(nfails + nreals == nrealisations); CHK(nfails < N/1000); diff --git a/src/test_sdis_solve_medium.c b/src/test_sdis_solve_medium.c @@ -203,6 +203,7 @@ main(int argc, char** argv) struct s3dut_mesh* msh0 = NULL; struct s3dut_mesh* msh1 = NULL; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_medium* solid0 = NULL; struct sdis_medium* solid1 = NULL; @@ -362,8 +363,10 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); printf("Shape0 temperature = "STR(Tf0)" ~ %g +/- %g\n", T.E, T.SE); - printf("#failures = %lu/%lu\n", (unsigned long)nfails, N); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); + printf("#failures = %lu/%lu\n\n", (unsigned long)nfails, N); CHK(eq_eps(T.E, Tf0, T.SE)); CHK(nreals + nfails == N); OK(sdis_estimator_ref_put(estimator)); @@ -372,8 +375,10 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); printf("Shape1 temperature = "STR(Tf1)" ~ %g +/- %g\n", T.E, T.SE); - printf("#failures = %lu/%lu\n", (unsigned long)nfails, N); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); + printf("#failures = %lu/%lu\n\n", (unsigned long)nfails, N); CHK(eq_eps(T.E, Tf1, T.SE)); CHK(nreals + nfails == N); OK(sdis_estimator_ref_put(estimator)); @@ -398,10 +403,12 @@ main(int argc, char** argv) BA(sdis_solve_medium(scn, N, solid1, trange, 1.f, -1, 0, 0, &estimator)); OK(sdis_solve_medium(scn, Np, solid0, trange, 1.f, -1, 0, 0, &estimator)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); ref = Tf0 * v0/v + Tf1 * v1/v; printf("Shape0 + Shape1 temperature = %g ~ %g +/- %g\n", ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, Np); CHK(eq_eps(T.E, ref, T.SE*3)); diff --git a/src/test_sdis_solve_medium_2d.c b/src/test_sdis_solve_medium_2d.c @@ -188,6 +188,7 @@ main(int argc, char** argv) { struct mem_allocator allocator; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_medium* solid0 = NULL; struct sdis_medium* solid1 = NULL; @@ -336,10 +337,12 @@ main(int argc, char** argv) /* Estimate the temperature of the square */ OK(sdis_solve_medium(scn, N, solid0, trange, 1.f, -1, 0, 0, &estimator)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); printf("Square temperature = "STR(Tf0)" ~ %g +/- %g\n", T.E, T.SE); - printf("#failures = %lu / %lu\n", (unsigned long)nfails, N); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); + printf("#failures = %lu / %lu\n\n", (unsigned long)nfails, N); CHK(eq_eps(T.E, Tf0, T.SE)); CHK(nreals + nfails == N); OK(sdis_estimator_ref_put(estimator)); @@ -347,10 +350,12 @@ main(int argc, char** argv) /* Estimate the temperature of the disk */ OK(sdis_solve_medium(scn, N, solid1, trange, 1.f, -1, 0, 0, &estimator)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); printf("Disk temperature = "STR(Tf1)" ~ %g +/- %g\n", T.E, T.SE); - printf("#failures = %lu / %lu\n", (unsigned long)nfails, N); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); + printf("#failures = %lu / %lu\n\n", (unsigned long)nfails, N); CHK(eq_eps(T.E, Tf1, T.SE)); CHK(nreals + nfails == N); OK(sdis_estimator_ref_put(estimator)); @@ -369,10 +374,12 @@ main(int argc, char** argv) BA(sdis_solve_medium(scn, N, solid1, trange, 1.f, -1, 0, 0, &estimator)); OK(sdis_solve_medium(scn, Np, solid0, trange, 1.f, -1, 0, 0, &estimator)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); ref = Tf0 * a0/a + Tf1 * a1/a; printf("Square + Disk temperature = %g ~ %g +/- %g\n", ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu / %lu\n", (unsigned long)nfails, Np); CHK(eq_eps(T.E, ref, 3*T.SE)); CHK(nreals + nfails == Np); 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 F = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_medium* solid = NULL; struct sdis_medium* fluid = NULL; @@ -379,9 +380,14 @@ main(int argc, char** argv) BA(sdis_estimator_get_temperature(NULL, &T)); OK(sdis_estimator_get_temperature(estimator, &T)); + BA(sdis_estimator_get_realisation_time(estimator, NULL)); + BA(sdis_estimator_get_realisation_time(NULL, &time)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); + ref = 300; printf("Temperature at (%g, %g, %g) = %g ~ %g +/- %g\n", SPLIT3(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(nfails + nreals == N); 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 time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_data* data = NULL; struct sdis_estimator* estimator = NULL; @@ -246,11 +247,13 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); /* Print the estimation results */ ref = 350 * pos[2] + (1-pos[2]) * 300; printf("Temperature at (%g, %g, %g) = %g ~ %g +/- %g\n", SPLIT3(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); /* Check the results */ diff --git a/src/test_sdis_solve_probe2_2d.c b/src/test_sdis_solve_probe2_2d.c @@ -143,6 +143,7 @@ main(int argc, char** argv) { struct mem_allocator allocator; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_data* data = NULL; struct sdis_estimator* estimator = NULL; @@ -244,11 +245,14 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); + /* Print the estimation results */ ref = 350 * pos[0] + (1-pos[0]) * 300; printf("Temperature at (%g, %g) = %g ~ %g +/- %g\n", SPLIT2(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); /* Check the results */ diff --git a/src/test_sdis_solve_probe3.c b/src/test_sdis_solve_probe3.c @@ -168,6 +168,7 @@ main(int argc, char** argv) { struct mem_allocator allocator; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_data* data = NULL; struct sdis_estimator* estimator = NULL; @@ -302,11 +303,13 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); /* Print the estimation results */ ref = 350 * pos[2] + (1-pos[2]) * 300; printf("Temperature at (%g, %g, %g) = %g ~ %g +/- %g\n", SPLIT3(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); /* Check the results */ diff --git a/src/test_sdis_solve_probe3_2d.c b/src/test_sdis_solve_probe3_2d.c @@ -165,6 +165,7 @@ main(int argc, char** argv) { struct mem_allocator allocator; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_data* data = NULL; struct sdis_estimator* estimator = NULL; @@ -295,11 +296,13 @@ main(int argc, char** argv) OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); /* Print the estimation results */ ref = 350 * pos[0] + (1-pos[0]) * 300; printf("Temperature at (%g, %g) = %g ~ %g +/- %g\n", SPLIT2(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); /* Check the results */ diff --git a/src/test_sdis_solve_probe_2d.c b/src/test_sdis_solve_probe_2d.c @@ -135,6 +135,7 @@ main(int argc, char** argv) { struct mem_allocator allocator; struct sdis_mc T = SDIS_MC_NULL; + struct sdis_mc time = SDIS_MC_NULL; struct sdis_device* dev = NULL; struct sdis_medium* solid = NULL; struct sdis_medium* fluid = NULL; @@ -204,12 +205,13 @@ main(int argc, char** argv) OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator)); OK(sdis_estimator_get_realisation_count(estimator, &nreals)); OK(sdis_estimator_get_failure_count(estimator, &nfails)); - OK(sdis_estimator_get_temperature(estimator, &T)); + OK(sdis_estimator_get_realisation_time(estimator, &time)); ref = 300; printf("Temperature at (%g, %g) = %g ~ %g +/- %g\n", SPLIT2(pos), ref, T.E, T.SE); + printf("Time per realisation (in usec) = %g +/- %g\n", time.E, time.SE); printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(nfails + nreals == N); diff --git a/src/test_sdis_utils.c b/src/test_sdis_utils.c @@ -262,6 +262,10 @@ check_green_function(struct sdis_green_function* green) CHK(E + SE >= mc.E - mc.SE); CHK(E - SE <= mc.E + mc.SE); + OK(sdis_estimator_get_realisation_time(estimator, &mc)); + printf("Green per realisation time (in usec) = %g +/- %g\n", + mc.E, mc.SE); + OK(sdis_estimator_ref_put(estimator)); }