stardis-solver

Solve coupled heat transfers
git clone git://git.meso-star.fr/stardis-solver.git
Log | Files | Refs | README | LICENSE

commit b456a9e6d4b6ebc154e308673ef91b5de6f9b4fd
parent 67d8279567a65f9b094b7d5b21eb61a163e73d56
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Wed, 12 Jan 2022 17:00:11 +0100

Make the time_rewind function independent to the medium type

Diffstat:
Msrc/sdis_heat_path_boundary_Xd_c.h | 11+++++++++--
Msrc/sdis_heat_path_conductive_Xd.h | 7+++++--
Msrc/sdis_heat_path_convective_Xd.h | 52++++++----------------------------------------------
Msrc/sdis_medium_c.h | 31++++++++++++++++++++++---------
Msrc/sdis_misc.h | 16++++++++--------
Msrc/sdis_misc_Xd.h | 44+++++++++++++++++++-------------------------
6 files changed, 69 insertions(+), 92 deletions(-)

diff --git a/src/sdis_heat_path_boundary_Xd_c.h b/src/sdis_heat_path_boundary_Xd_c.h @@ -806,20 +806,27 @@ XD(solid_reinjection) (struct sdis_medium* solid, struct XD(solid_reinjection_args)* args) { + struct solid_props props = SOLID_PROPS_NULL; double reinject_dst_m; /* Reinjection distance in meters */ + double mu; res_T res = RES_OK; ASSERT(solid && XD(check_solid_reinjection_args)(args)); reinject_dst_m = args->reinjection->distance * args->fp_to_meter; + res = solid_get_properties(solid, &args->rwalk->vtx, &props); + if(res != RES_OK) goto error; + /* Manage the volumic power */ res = XD(handle_volumic_power) (solid, args->rwalk_ctx, args->rwalk, reinject_dst_m, args->T); if(res != RES_OK) goto error; /* Time rewind */ - res = XD(solid_time_rewind) - (solid, args->rng, reinject_dst_m, args->rwalk_ctx, args->rwalk, args->T); + args->rwalk->mdm = solid; /* Medium into which the time is rewind */ + mu = (2*DIM*props.lambda)/(props.rho*props.cp*reinject_dst_m*reinject_dst_m); + res = XD(time_rewind) + (mu, props.t0, args->rng, args->rwalk, args->rwalk_ctx, args->T); if(res != RES_OK) goto error; /* Test if a limit condition was reached */ diff --git a/src/sdis_heat_path_conductive_Xd.h b/src/sdis_heat_path_conductive_Xd.h @@ -434,7 +434,9 @@ XD(conductive_path) struct sXd(hit) hit0, hit1; struct solid_props props = SOLID_PROPS_NULL; double power_term = 0; + double mu; float delta; /* Random walk numerical parameter */ + double delta_m; float dir0[DIM], dir1[DIM]; float org[DIM]; @@ -492,8 +494,9 @@ XD(conductive_path) } /* Rewind the time */ - res = XD(solid_time_rewind) - (rwalk->mdm, rng, delta * scn->fp_to_meter, ctx, rwalk, T); + delta_m = delta * scn->fp_to_meter; + mu = (2*DIM*props.lambda)/(props.rho*props.cp*delta_m*delta_m); + res = XD(time_rewind)(mu, props.t0, rng, rwalk, ctx, T); if(res != RES_OK) goto error; if(T->done) break; /* Limit condition was reached */ diff --git a/src/sdis_heat_path_convective_Xd.h b/src/sdis_heat_path_convective_Xd.h @@ -172,8 +172,8 @@ XD(convective_path) if(SXD_HIT_NONE(&rwalk->hit)) { log_err(scn->dev, -"%s: the position %g %g %g lies in the surrounding fluid whose temperature must \n" -"be known.\n", FUNC_NAME, SPLIT3(rwalk->vtx.P)); + "%s: the position %g %g %g lies in the surrounding fluid whose " + "temperature must be known.\n", FUNC_NAME, SPLIT3(rwalk->vtx.P)); res = RES_BAD_OP; goto error; } @@ -238,7 +238,7 @@ XD(convective_path) struct sXd(primitive) prim; struct fluid_props props = FLUID_PROPS_NULL; double hc; - double mu, tau; + double mu; /* Fetch fluid properties */ res = fluid_get_properties(rwalk->mdm, &rwalk->vtx, &props); @@ -249,49 +249,9 @@ XD(convective_path) /* Sample the time using the upper bound. */ mu = enc->hc_upper_bound / (props.rho * props.cp) * enc->S_over_V; - tau = ssp_ran_exp(rng, mu); - - /* Increment the elapsed time */ - ASSERT(rwalk->vtx.time > props.t0); - rwalk->elapsed_time += MMIN(tau, rwalk->vtx.time - props.t0); - - if(rwalk->vtx.time != INF) { - rwalk->vtx.time = MMAX(rwalk->vtx.time - tau, props.t0); /* Time rewind */ - - /* Register the new vertex against the heat path */ - res = XD(register_heat_vertex_in_fluid)(scn, ctx, rwalk, T->value); - if(res != RES_OK) goto error; - - if(rwalk->vtx.time == props.t0) { - /* Check the initial condition. */ - tmp = fluid_get_temperature(rwalk->mdm, &rwalk->vtx); - if(tmp >= 0) { - T->value += tmp; - T->done = 1; - if(ctx->heat_path) { - /* Update the registered vertex data */ - struct sdis_heat_vertex* vtx; - vtx = heat_path_get_last_vertex(ctx->heat_path); - vtx->time = rwalk->vtx.time; - vtx->weight = T->value; - } - - if(ctx->green_path) { - res = green_path_set_limit_vertex(ctx->green_path, rwalk->mdm, - &rwalk->vtx, rwalk->elapsed_time); - if(res != RES_OK) goto error; - } - goto exit; - } - /* The initial condition should have been reached. */ - log_err(scn->dev, - "%s: undefined initial condition. " - "Time is %g but the temperature remains unknown.\n", - FUNC_NAME, props.t0); - res = RES_BAD_OP; - goto error; - } - } + res = XD(time_rewind)(mu, props.t0, rng, rwalk, ctx, T); + if(res != RES_OK) goto error; + if(T->done) break; /* Limit condition was reached */ /* Uniformly sample the enclosure. */ #if DIM == 2 diff --git a/src/sdis_medium_c.h b/src/sdis_medium_c.h @@ -104,13 +104,6 @@ static const struct solid_props SOLID_PROPS_NULL = SOLID_PROPS_NULL__; return RES_OK; \ } -static FINLINE unsigned -medium_get_id(const struct sdis_medium* mdm) -{ - ASSERT(mdm); - return mdm->id.index; -} - /******************************************************************************* * Fluid local functions ******************************************************************************/ @@ -275,7 +268,7 @@ solid_get_properties /******************************************************************************* * Generic functions ******************************************************************************/ -static FINLINE double +static INLINE double medium_get_temperature (const struct sdis_medium* mdm, const struct sdis_rwalk_vertex* vtx) { @@ -289,7 +282,7 @@ medium_get_temperature return temp; } -static FINLINE double +static INLINE double medium_get_t0(const struct sdis_medium* mdm) { double t0; @@ -301,6 +294,26 @@ medium_get_t0(const struct sdis_medium* mdm) } return t0; } + +static INLINE unsigned +medium_get_id(const struct sdis_medium* mdm) +{ + ASSERT(mdm); + return mdm->id.index; +} + +static INLINE const char* +medium_type_to_string(const enum sdis_medium_type type) +{ + const char* str = "none"; + switch(type) { + case SDIS_FLUID: str = "fluid"; break; + case SDIS_SOLID: str = "solid"; break; + default: FATAL("Unreachable code.\n"); break; + } + return str; +} + #undef MDM_TYPE #undef MDM_TYPE_solid #undef MDM_TYPE_fluid diff --git a/src/sdis_misc.h b/src/sdis_misc.h @@ -148,21 +148,21 @@ register_heat_vertex } extern LOCAL_SYM res_T -solid_time_rewind_2d - (struct sdis_medium* mdm, /* Solid into which the time is rewinded */ +time_rewind_2d + (const double mu, + const double t0, /* Initial time */ struct ssp_rng* rng, - const double dist_in_meter, - const struct rwalk_context* ctx, struct rwalk_2d* rwalk, + const struct rwalk_context* ctx, struct temperature_2d* T); extern LOCAL_SYM res_T -solid_time_rewind_3d - (struct sdis_medium* mdm, /* Solid into which the time is rewinded */ +time_rewind_3d + (const double mu, + const double t0, /* Initial time */ struct ssp_rng* rng, - const double dist_in_meter, - const struct rwalk_context* ctx, struct rwalk_3d* rwalk, + const struct rwalk_context* ctx, struct temperature_3d* T); /* Check the validity of the parametric coordinate onto a 2D primitive. If it diff --git a/src/sdis_misc_Xd.h b/src/sdis_misc_Xd.h @@ -24,47 +24,41 @@ #include "sdis_Xd_begin.h" res_T -XD(solid_time_rewind) - (struct sdis_medium* mdm, +XD(time_rewind) + (const double mu, + const double t0, struct ssp_rng* rng, - const double dist_in_meter, - const struct rwalk_context* ctx, struct XD(rwalk)* rwalk, + const struct rwalk_context* ctx, struct XD(temperature)* T) { - struct solid_props props = SOLID_PROPS_NULL; - double tau, mu; double temperature; + double tau; res_T res = RES_OK; - ASSERT(mdm && rng && ctx && rwalk && dist_in_meter > 0); - ASSERT(sdis_medium_get_type(mdm) == SDIS_SOLID); - ASSERT(T->done == 0); - - /* Fetch physical properties */ - res = solid_get_properties(mdm, &rwalk->vtx, &props); - if(res != RES_OK) goto error; + ASSERT(rwalk && rng && T); - /* Sample the time to reroll */ - mu = (2*DIM*props.lambda)/(props.rho*props.cp*dist_in_meter*dist_in_meter); + /* Sample the time using the upper bound. */ tau = ssp_ran_exp(rng, mu); /* Increment the elapsed time */ - ASSERT(rwalk->vtx.time >= props.t0); - rwalk->elapsed_time += MMIN(tau, rwalk->vtx.time - props.t0); + ASSERT(rwalk->vtx.time >= t0); + rwalk->elapsed_time += MMIN(tau, rwalk->vtx.time - t0); if(IS_INF(rwalk->vtx.time)) goto exit; /* Steady computation */ /* Time rewind */ - rwalk->vtx.time = MMAX(rwalk->vtx.time - tau, props.t0); + rwalk->vtx.time = MMAX(rwalk->vtx.time - tau, t0); /* Time rewind */ /* The path does not reach the limit condition */ - if(rwalk->vtx.time > props.t0) goto exit; + if(rwalk->vtx.time > t0) goto exit; - /* Fetch initial temperature */ - temperature = solid_get_temperature(mdm, &rwalk->vtx); + /* Fetch the initial temperature */ + temperature = medium_get_temperature(rwalk->mdm, &rwalk->vtx); if(temperature < 0) { - log_err(mdm->dev, "%s: the path reaches the limit condition but the " - "temperature remains unknown.\n", FUNC_NAME); + log_err(rwalk->mdm->dev, "the path reaches the limit condition but the " + "%s temperature remains unknown -- position=%g, %g, %g\n", + medium_type_to_string(sdis_medium_get_type(rwalk->mdm)), + SPLIT3(rwalk->vtx.P)); res = RES_BAD_ARG; goto error; } @@ -82,8 +76,8 @@ XD(solid_time_rewind) } if(ctx->green_path) { - res = green_path_set_limit_vertex(ctx->green_path, mdm, &rwalk->vtx, - rwalk->elapsed_time); + res = green_path_set_limit_vertex(ctx->green_path, rwalk->mdm, + &rwalk->vtx, rwalk->elapsed_time); if(res != RES_OK) goto error; }