commit ebe4a68605dcdb88a90eaac73c66bc57f277cf2a
parent 38de701861bd9b5cbd3194a517ac06d401cdba6d
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Fri, 31 May 2024 16:27:35 +0200
Follow the enclosure rather than the medium
From now on, the path saves the enclosure in which it is located, rather
than the medium. This change has an impact on several parts of the code,
since it modifies the central data structures.
In particular, Robin's condition is no longer handled in the convective
path; the known fluid temperature is obtained at the boundary. This not
only makes sense from a physical point of view, since Robin's condition
is defined at a boundary, but is also a requirement, since in such a
situation, the enclosure can contain several media. Thus, Robin's
"fluid", i.e. the fluid whose temperature is known, cannot be recovered
from the enclosure. It must be recovered from the boundary.
Diffstat:
22 files changed, 574 insertions(+), 439 deletions(-)
diff --git a/src/sdis_heat_path.h b/src/sdis_heat_path.h
@@ -17,6 +17,7 @@
#define SDIS_HEAT_PATH_H
#include "sdis.h"
+#include "sdis_scene_c.h"
#include <rsys/dynamic_array.h>
#include <rsys/dynamic_array_size_t.h>
@@ -88,7 +89,7 @@ get_picard_order(const struct rwalk_context* ctx)
******************************************************************************/
struct rwalk {
struct sdis_rwalk_vertex vtx; /* Position and time of the Random walk */
- struct sdis_medium* mdm; /* Medium in which the random walk lies */
+ unsigned enc_id; /* Id of the enclosure in which the random walk lies */
struct s2d_hit hit_2d;
struct s3d_hit hit_3d;
@@ -100,7 +101,7 @@ struct rwalk {
};
#define RWALK_NULL__ { \
SDIS_RWALK_VERTEX_NULL__, \
- NULL, \
+ ENCLOSURE_ID_NULL, \
S2D_HIT_NULL__, \
S3D_HIT_NULL__, \
{0,0,0}, \
diff --git a/src/sdis_heat_path_boundary_Xd.h b/src/sdis_heat_path_boundary_Xd.h
@@ -25,6 +25,58 @@
#include "sdis_Xd_begin.h"
/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+/* This function checks whether the random walk is on a boundary and, if so,
+ * verifies that the temperature of the medium attached to the interface is
+ * known. This medium can be different from the medium of the enclosure. Indeed,
+ * the enclosure can contain several media used to set the temperatures of
+ * several boundary conditions. Hence this function, which queries the medium on
+ * the trajectory coming from a boundary */
+static res_T
+XD(handle_known_medium_temperature)
+ (struct sdis_scene* scn,
+ struct rwalk_context* ctx,
+ struct rwalk* rwalk,
+ struct temperature* T)
+{
+ struct sdis_interface* interf = NULL;
+ struct sdis_medium* mdm = NULL;
+ double temperature = SDIS_TEMPERATURE_NONE;
+ res_T res = RES_OK;
+ ASSERT(scn && ctx && rwalk && T);
+
+ /* Not at an interface */
+ if(SXD_HIT_NONE(&rwalk->XD(hit))) return RES_OK; /* Nothing to do */
+
+ interf = scene_get_interface(scn, rwalk->XD(hit).prim.prim_id);
+ mdm = rwalk->hit_side==SDIS_FRONT ? interf->medium_front: interf->medium_back;
+
+ temperature = medium_get_temperature(mdm, &rwalk->vtx);
+
+ /* Check if the temperature is known */
+ if(SDIS_TEMPERATURE_IS_UNKNOWN(temperature)) goto exit;
+
+ T->value += temperature;
+ T->done = 1;
+
+ if(ctx->green_path) {
+ res = green_path_set_limit_vertex
+ (ctx->green_path, mdm, &rwalk->vtx, rwalk->elapsed_time);
+ if(res != RES_OK) goto error;
+ }
+
+ if(ctx->heat_path) {
+ heat_path_get_last_vertex(ctx->heat_path)->weight = T->value;
+ }
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
+/*******************************************************************************
* Local functions
******************************************************************************/
res_T
@@ -42,7 +94,7 @@ XD(boundary_path)
double tmp;
res_T res = RES_OK;
ASSERT(scn && ctx && rwalk && rng && T);
- ASSERT(rwalk->mdm == NULL);
+ ASSERT(rwalk->enc_id == ENCLOSURE_ID_NULL);
ASSERT(!SXD_HIT_NONE(&rwalk->XD(hit)));
XD(setup_interface_fragment)
@@ -83,6 +135,32 @@ XD(boundary_path)
}
if(res != RES_OK) goto error;
+#if 1
+ if(T->done) goto exit;
+
+ /* Handling limit boundary condition, i.e. the trajectory originates from a
+ * boundary and the medium temperature is known (e.g. Robin's condition). To
+ * simplify data description, we allow in such a situation, to define several
+ * medium on the same enclosure, each with a fixed temperature, i.e. different
+ * conditions are defined for the different interfaces that detour the
+ * enclosure. As a result, no path can be sampled in this enclosure, which is
+ * beyond the system boundary. The boundary medium must therefore be
+ * interrogated from the interface.
+ *
+ * Note that we check this boundary condition with convective paths to handle
+ * Robin's boundary conditions. But we also make this check when passing
+ * through conduction when there's no reason why a solid should have a fixed
+ * temperature and not its boundary: it should be a Dirichlet condition.
+ * Although it's not physical, such systems can still be defined
+ * computationally, and in fact it's also a handy way of testing
+ * border cases */
+ if(T->func == XD(convective_path) || T->func == XD(conductive_path)) {
+ res = XD(handle_known_medium_temperature)(scn, ctx, rwalk, T);
+ if(res != RES_OK) goto error;
+ if(T->done) goto exit; /* That's all folks */
+ }
+#endif
+
exit:
return res;
error:
diff --git a/src/sdis_heat_path_boundary_Xd_c.h b/src/sdis_heat_path_boundary_Xd_c.h
@@ -26,18 +26,18 @@
#include "sdis_Xd_begin.h"
struct XD(find_reinjection_ray_args) {
- const struct sdis_medium* solid; /* Medium into which the reinjection occurs */
const struct rwalk* rwalk; /* Current random walk state */
float dir0[DIM]; /* Challenged ray direction */
float dir1[DIM]; /* Challenged ray direction */
double distance; /* Maximum reinjection distance */
+ unsigned solid_enc_id; /* Enclosure id into which the reinjection occurs */
/* Define if the random walk position can be moved or not to find a valid
* reinjection direction */
int can_move;
};
static const struct XD(find_reinjection_ray_args)
-XD(FIND_REINJECTION_RAY_ARGS_NULL) = { NULL, NULL, {0}, {0}, 0, 0 };
+XD(FIND_REINJECTION_RAY_ARGS_NULL) = { NULL, {0}, {0}, 0, ENCLOSURE_ID_NULL, 0 };
struct XD(reinjection_ray) {
double org[DIM]; /* Origin of the reinjection */
@@ -55,54 +55,112 @@ XD(REINJECTION_RAY_NULL) = { {0}, {0}, 0, SXD_HIT_NULL__, 0 };
/*******************************************************************************
* Helper functions
******************************************************************************/
-static INLINE int
+static INLINE res_T
XD(check_find_reinjection_ray_args)
- (const struct XD(find_reinjection_ray_args)* args)
+ (struct sdis_scene* scn,
+ const struct XD(find_reinjection_ray_args)* args)
{
- return args
- && args->solid
- && args->rwalk
- && args->distance > 0
- && fX(is_normalized)(args->dir0)
- && fX(is_normalized)(args->dir1);
+ const struct enclosure* enc = NULL;
+ struct sdis_medium* mdm = NULL;
+ res_T res = RES_OK;
+ ASSERT(scn);
+
+ /* Check pointers */
+ if(!args || !args->rwalk) return RES_BAD_ARG;
+
+ /* Check distance */
+ if(args->distance <= 0) return RES_BAD_ARG;
+
+ /* Check directions */
+ if(!fX(is_normalized)(args->dir0) || !fX(is_normalized)(args->dir1)) {
+ return RES_BAD_ARG;
+ }
+
+ /* Check enclosure id */
+ if(args->solid_enc_id == ENCLOSURE_ID_NULL) {
+ return RES_BAD_ARG;
+ }
+ enc = scene_get_enclosure(scn, args->solid_enc_id);
+
+ /* Check the enclosure */
+ enc = scene_get_enclosure(scn, args->solid_enc_id);
+ if(enc->medium_id != MEDIUM_ID_MULTI) {
+ if((res = scene_get_enclosure_medium(scn, enc, &mdm)) != RES_OK) return res;
+ if(sdis_medium_get_type(mdm) != SDIS_SOLID) {
+ res = RES_BAD_ARG;
+ }
+ }
+
+ return RES_OK;
}
-static INLINE int
+static INLINE res_T
XD(check_sample_reinjection_step_args)
- (const struct sample_reinjection_step_args* args)
+ (struct sdis_scene* scn,
+ const struct sample_reinjection_step_args* args)
{
- return args
- && args->rng
- && args->solid
- && args->solid->type == SDIS_SOLID
- && args->rwalk
- && args->distance > 0
- && (unsigned)args->side < SDIS_SIDE_NULL__;
+ const struct enclosure* enc = NULL;
+ struct sdis_medium* mdm = NULL;
+ res_T res = RES_OK;
+ ASSERT(scn);
+
+ /* Check pointers */
+ if(!args || !args->rng || !args->rwalk) return RES_BAD_ARG;
+
+ /* Check distance */
+ if(args->distance <= 0) return RES_BAD_ARG;
+
+ /* Check side */
+ if((unsigned)args->side >= SDIS_SIDE_NULL__) return RES_BAD_ARG;
+
+ /* Check enclosure id */
+ if(args->solid_enc_id == ENCLOSURE_ID_NULL) {
+ return RES_BAD_ARG;
+ }
+
+ /* Check the enclosure */
+ enc = scene_get_enclosure(scn, args->solid_enc_id);
+ if(enc->medium_id != MEDIUM_ID_MULTI) {
+ if((res = scene_get_enclosure_medium(scn, enc, &mdm)) != RES_OK) return res;
+ if(sdis_medium_get_type(mdm) != SDIS_SOLID) {
+ return RES_BAD_ARG;
+ }
+ }
+
+ return RES_OK;
}
-static INLINE int
+static INLINE res_T
XD(check_reinjection_step)(const struct reinjection_step* step)
{
- return step
- && fX(is_normalized)(step->direction)
- && step->distance > 0;
+ /* Check pointer */
+ if(!step) return RES_BAD_ARG;
+
+ /* Check direction */
+ if(!fX(is_normalized)(step->direction)) return RES_BAD_ARG;
+
+ /* Check distance */
+ if(step->distance <= 0) return RES_BAD_ARG;
+
+ return RES_OK;
}
-static INLINE int
+static INLINE res_T
XD(check_solid_reinjection_args)(const struct solid_reinjection_args* args)
{
- return args
- && XD(check_reinjection_step)(args->reinjection)
- && args->rng
- && args->rwalk
- && args->rwalk_ctx
- && args->T
- && args->fp_to_meter > 0;
+ /* Check pointers */
+ if(!args || !args->rng || !args->rwalk || !args->rwalk_ctx || !args->T)
+ return RES_BAD_ARG;
+
+ /* Check unit */
+ if(args->fp_to_meter <= 0) return RES_BAD_ARG;
+
+ return XD(check_reinjection_step)(args->reinjection);
}
/* Check that the interface fragment is consistent with the current state of
* the random walk */
-static INLINE int
+static INLINE res_T
XD(check_rwalk_fragment_consistency)
(const struct rwalk* rwalk,
const struct sdis_interface_fragment* frag)
@@ -110,20 +168,32 @@ XD(check_rwalk_fragment_consistency)
double N[DIM];
double uv[2] = {0, 0};
ASSERT(rwalk && frag);
+
+ /* Check intersection */
+ if(SXD_HIT_NONE(&rwalk->XD(hit))) return RES_BAD_ARG;
+
+ /* Check positions */
+ if(!dX(eq_eps)(rwalk->vtx.P, frag->P, 1.e-6)) return RES_BAD_ARG;
+
+ /* Check normals */
dX(normalize)(N, dX_set_fX(N, rwalk->XD(hit).normal));
- if( SXD_HIT_NONE(&rwalk->XD(hit))
- || !dX(eq_eps)(rwalk->vtx.P, frag->P, 1.e-6)
- || !dX(eq_eps)(N, frag->Ng, 1.e-6)
- || !( (IS_INF(rwalk->vtx.time) && IS_INF(frag->time))
- || eq_eps(rwalk->vtx.time, frag->time, 1.e-6))) {
- return 0;
+ if(!dX(eq_eps)(N, frag->Ng, 1.e-6)) return RES_BAD_ARG;
+
+ /* Check time */
+ if(!eq_eps(rwalk->vtx.time, frag->time, 1.e-6)
+ && !(IS_INF(rwalk->vtx.time) && IS_INF(frag->time))) {
+ return RES_BAD_ARG;
}
+
+ /* Check parametric coordinates */
#if (SDIS_XD_DIMENSION == 2)
uv[0] = rwalk->XD(hit).u;
#else
d2_set_f2(uv, rwalk->XD(hit).uv);
#endif
- return d2_eq_eps(uv, frag->uv, 1.e-6);
+ if(!d2_eq_eps(uv, frag->uv, 1.e-6)) return RES_BAD_ARG;
+
+ return RES_OK;
}
static void
@@ -146,7 +216,7 @@ XD(sample_reinjection_dir)
const uint64_t r = ssp_rng_uniform_uint64(rng, 0, 1);
ASSERT(rwalk && rng && dir);
ASSERT(!SXD_HIT_NONE(&rwalk->XD(hit)));
- ASSERT(!rwalk->mdm);
+ ASSERT(rwalk->enc_id == ENCLOSURE_ID_NULL);
if(r) {
dir[0] = rwalk->XD(hit).normal[0] - rwalk->XD(hit).normal[1];
@@ -163,7 +233,7 @@ XD(sample_reinjection_dir)
float frame[9];
ASSERT(rwalk && rng && dir);
ASSERT(!SXD_HIT_NONE(&rwalk->XD(hit)));
- ASSERT(!rwalk->mdm);
+ ASSERT(rwalk->enc_id == ENCLOSURE_ID_NULL);
ASSERT(fX(is_normalized)(rwalk->XD(hit).normal));
ssp_ran_circle_uniform_float(rng, dir, NULL);
@@ -315,10 +385,10 @@ XD(find_reinjection_ray)
/* # attempts to find a ray direction */
int MAX_ATTEMPTS = 1;
- /* Physical properties */
- struct sdis_interface* interf;
- struct sdis_medium* mdm0;
- struct sdis_medium* mdm1;
+ /* Enclosures */
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+ unsigned enc0_id = ENCLOSURE_ID_NULL;
+ unsigned enc1_id = ENCLOSURE_ID_NULL;
struct hit_filter_data filter_data = HIT_FILTER_DATA_NULL;
struct sXd(hit) hit;
@@ -334,12 +404,11 @@ XD(find_reinjection_ray)
float org[DIM];
const float range[2] = {0, FLT_MAX};
enum sdis_side side;
- unsigned enc_id = ENCLOSURE_ID_NULL;
int iattempt = 0;
res_T res = RES_OK;
ASSERT(scn && args && ray);
- ASSERT(XD(check_find_reinjection_ray_args)(args));
+ ASSERT(XD(check_find_reinjection_ray_args)(scn, args) == RES_OK);
*ray = XD(REINJECTION_RAY_NULL);
MAX_ATTEMPTS = args->can_move ? 2 : 1;
@@ -358,46 +427,32 @@ XD(find_reinjection_ray)
SXD(scene_view_trace_ray
(scn->sXd(view), org, args->dir1, range, &filter_data, &hit1));
- /* Retrieve the medium at the reinjection pos along dir0 */
+ /* Retrieve the enclosure at the reinjection pos along dir0 */
if(!SXD_HIT_NONE(&hit0)) {
- interf = scene_get_interface(scn, hit0.prim.prim_id);
+ scene_get_enclosure_ids(scn, hit0.prim.prim_id, enc_ids);
side = fX(dot)(args->dir0, hit0.normal) < 0 ? SDIS_FRONT : SDIS_BACK;
- mdm0 = interface_get_medium(interf, side);
+ enc0_id = enc_ids[side];
} else {
XD(move_pos)(dX(set)(tmp, ray->org), args->dir0, (float)args->distance);
- res = scene_get_enclosure_id_in_closed_boundaries(scn, tmp, &enc_id);
- if(res == RES_BAD_OP) { enc_id = ENCLOSURE_ID_NULL; res = RES_OK; }
+ res = scene_get_enclosure_id_in_closed_boundaries(scn, tmp, &enc0_id);
+ if(res == RES_BAD_OP) { enc0_id = ENCLOSURE_ID_NULL; res = RES_OK; }
if(res != RES_OK) goto error;
-
- mdm0 = NULL;
- if(enc_id != ENCLOSURE_ID_NULL) {
- const struct enclosure* enc = scene_get_enclosure(scn, enc_id);
- res = scene_get_enclosure_medium(scn, enc, &mdm0);
- if(res != RES_OK) goto error;
- }
}
- /* Retrieve the medium at the reinjection pos along dir1 */
+ /* Retrieve the enclosure at the reinjection pos along dir1 */
if(!SXD_HIT_NONE(&hit1)) {
- interf = scene_get_interface(scn, hit1.prim.prim_id);
+ scene_get_enclosure_ids(scn, hit1.prim.prim_id, enc_ids);
side = fX(dot)(args->dir1, hit1.normal) < 0 ? SDIS_FRONT : SDIS_BACK;
- mdm1 = interface_get_medium(interf, side);
+ enc1_id = enc_ids[side];
} else {
XD(move_pos)(dX(set)(tmp, ray->org), args->dir1, (float)args->distance);
- res = scene_get_enclosure_id_in_closed_boundaries(scn, tmp, &enc_id);
- if(res == RES_BAD_OP) { enc_id = ENCLOSURE_ID_NULL; res = RES_OK; }
+ res = scene_get_enclosure_id_in_closed_boundaries(scn, tmp, &enc1_id);
+ if(res == RES_BAD_OP) { enc1_id = ENCLOSURE_ID_NULL; res = RES_OK; }
if(res != RES_OK) goto error;
-
- mdm1 = NULL;
- if(enc_id != ENCLOSURE_ID_NULL) {
- const struct enclosure* enc = scene_get_enclosure(scn, enc_id);
- res = scene_get_enclosure_medium(scn, enc, &mdm1);
- if(res != RES_OK) goto error;
- }
}
dst0 = dst1 = -1;
- if(mdm0 == args->solid) { /* Check reinjection consistency */
+ if(enc0_id == args->solid_enc_id) { /* Check reinjection consistency */
if(hit0.distance <= dst_adjusted) {
dst0 = hit0.distance;
} else {
@@ -405,7 +460,7 @@ XD(find_reinjection_ray)
hit0 = SXD_HIT_NULL;
}
}
- if(mdm1 == args->solid) { /* Check reinjection consistency */
+ if(enc1_id == args->solid_enc_id) { /* Check reinjection consistency */
if(hit1.distance <= dst_adjusted) {
dst1 = hit1.distance;
} else {
@@ -518,30 +573,25 @@ XD(find_reinjection_ray_and_check_validity)
struct XD(reinjection_ray)* ray)
{
double pos[DIM];
- struct sdis_medium* reinject_mdm;
res_T res = RES_OK;
ASSERT(scn && args && ray);
- ASSERT(XD(check_find_reinjection_ray_args)(args));
+ ASSERT(XD(check_find_reinjection_ray_args)(scn, args) == RES_OK);
/* Select a reinjection direction */
res = XD(find_reinjection_ray)(scn, args, ray);
if(res != RES_OK) goto error;
if(SXD_HIT_NONE(&ray->hit)) {
- const struct enclosure* enc = NULL;
unsigned enc_id = ENCLOSURE_ID_NULL;
/* Obtain the enclosure in which the reinjection position lies */
XD(move_pos)(dX(set)(pos, ray->org), ray->dir, (float)ray->dst);
res = scene_get_enclosure_id_in_closed_boundaries(scn, pos, &enc_id);
if(res != RES_OK) goto error;
- enc = scene_get_enclosure(scn, enc_id);
- /* Check medium consistency at the reinjection position */
- res = scene_get_enclosure_medium(scn, enc, &reinject_mdm);
- if(res != RES_OK) goto error;
- if(reinject_mdm != args->solid) {
+ /* Check enclosure consistency at the reinjection position */
+ if(enc_id != args->solid_enc_id) {
res = RES_BAD_OP;
goto error;
}
@@ -622,14 +672,14 @@ XD(sample_reinjection_step_solid_fluid)
XD(FIND_REINJECTION_RAY_ARGS_NULL);
struct XD(reinjection_ray) ray = XD(REINJECTION_RAY_NULL);
+ /* Enclosures */
+ const struct enclosure* solid_enc = NULL;
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+
/* In 2D it is useless to try to resample a reinjection direction since there
* is only one possible direction */
const int MAX_ATTEMPTS = DIM == 2 ? 1 : 10;
- /* Enclosure */
- const struct enclosure* solid_enclosure = NULL;
- unsigned enc_ids[2];
-
/* Miscellaneous variables */
float dir0[DIM]; /* Sampled direction */
float dir1[DIM]; /* Sampled direction reflected */
@@ -638,7 +688,7 @@ XD(sample_reinjection_step_solid_fluid)
/* Pre-conditions */
ASSERT(scn && args && step);
- ASSERT(XD(check_sample_reinjection_step_args)(args));
+ ASSERT(XD(check_sample_reinjection_step_args)(scn, args) == RES_OK);
/* Initialise the reinjection step */
*step = REINJECTION_STEP_NULL;
@@ -648,11 +698,11 @@ XD(sample_reinjection_step_solid_fluid)
* enclosure has no geometrical existence and it is sufficient to return a
* valid reinjection step which will be used to select the next step. Note
* that if the trajectory passes through the solid enclosure, it will stop,
- * i.e. the temperature of the solid should be fixed. If it doesn't, it's an
+ * i.e. the temperature of the solid should be fixed. If it doesn't, it's a
* error */
scene_get_enclosure_ids(scn, args->rwalk->XD(hit).prim.prim_id, enc_ids);
- solid_enclosure = scene_get_enclosure(scn, enc_ids[args->side]);
- if(solid_enclosure->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
+ solid_enc = scene_get_enclosure(scn, enc_ids[args->side]);
+ if(solid_enc->medium_id == MEDIUM_ID_MULTI) {
step->XD(hit) = SXD_HIT_NULL;
fX(normalize)(step->direction, args->rwalk->XD(hit).normal);
if(args->side == SDIS_BACK) fX(minus)(step->direction, step->direction);
@@ -675,7 +725,7 @@ XD(sample_reinjection_step_solid_fluid)
}
/* Find the reinjection step */
- find_reinject_ray_args.solid = args->solid;
+ find_reinject_ray_args.solid_enc_id = args->solid_enc_id;
find_reinject_ray_args.rwalk = args->rwalk;
find_reinject_ray_args.distance = args->distance;
find_reinject_ray_args.can_move = 1;
@@ -709,7 +759,7 @@ XD(sample_reinjection_step_solid_fluid)
/* Post-conditions */
ASSERT(dX(eq)(args->rwalk->vtx.P, ray.org));
- ASSERT(XD(check_reinjection_step)(step));
+ ASSERT(XD(check_reinjection_step)(step) == RES_OK);
exit:
return res;
@@ -745,23 +795,23 @@ XD(sample_reinjection_step_solid_solid)
const int MAX_ATTEMPTS = DIM == 2 ? 1 : 10;
/* Enclosure */
- const struct enclosure* enclosure_bck = NULL;
- const struct enclosure* enclosure_frt = NULL;
- int multiple_media_bck = 0;
- int multiple_media_frt = 0;
unsigned enc_ids[2];
+ const struct enclosure* enc_frt = NULL;
+ const struct enclosure* enc_bck = NULL;
float dir_frt_samp[DIM]; /* Sampled direction */
float dir_frt_refl[DIM]; /* Sampled direction reflected */
float dir_bck_samp[DIM]; /* Negated sampled direction */
float dir_bck_refl[DIM]; /* Negated sampled direction reflected */
+ int multi_frt = 0;
+ int multi_bck = 0;
int iattempt = 0; /* #attempts to find a reinjection dir */
res_T res = RES_OK;
/* Pre-conditions */
ASSERT(scn && args_frt && args_bck && step_frt && step_bck);
- ASSERT(XD(check_sample_reinjection_step_args)(args_frt));
- ASSERT(XD(check_sample_reinjection_step_args)(args_bck));
+ ASSERT(XD(check_sample_reinjection_step_args)(scn, args_frt) == RES_OK);
+ ASSERT(XD(check_sample_reinjection_step_args)(scn, args_bck) == RES_OK);
ASSERT(args_frt->side == SDIS_FRONT);
ASSERT(args_bck->side == SDIS_BACK);
ASSERT(SXD_PRIMITIVE_EQ(&args_frt->rwalk->XD(hit).prim, &args_bck->rwalk->XD(hit).prim));
@@ -783,24 +833,22 @@ XD(sample_reinjection_step_solid_solid)
* i.e. the temperature of the solid should be fixed. If it doesn't, it's an
* error */
scene_get_enclosure_ids(scn, args_frt->rwalk->XD(hit).prim.prim_id, enc_ids);
- enclosure_frt = scene_get_enclosure(scn, enc_ids[SDIS_FRONT]);
- enclosure_bck = scene_get_enclosure(scn, enc_ids[SDIS_BACK]);
- if(enclosure_frt->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
+ enc_frt = scene_get_enclosure(scn, enc_ids[SDIS_FRONT]);
+ enc_bck = scene_get_enclosure(scn, enc_ids[SDIS_BACK]);
+ if(enc_frt->medium_id == MEDIUM_ID_MULTI) {
step_frt->XD(hit) = SXD_HIT_NULL;
fX(normalize)(step_frt->direction, args_frt->rwalk->XD(hit).normal);
step_frt->distance = (float)args_frt->distance;
- multiple_media_frt = 1;
+ multi_frt = 1;
}
- if(enclosure_bck->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
+ if(enc_bck->medium_id == MEDIUM_ID_MULTI) {
step_bck->XD(hit) = SXD_HIT_NULL;
fX(normalize)(step_bck->direction, args_bck->rwalk->XD(hit).normal);
- fX(minus)(step_bck->direction, step_bck->direction);
step_bck->distance = (float)args_bck->distance;
- multiple_media_bck = 1;
+ multi_bck = 1;
}
- if(multiple_media_frt && multiple_media_bck)
- goto exit; /* That's all folks! */
+ if(multi_frt && multi_bck) goto exit; /* That's all folks */
dX(set)(rwalk_pos_backup, rwalk->vtx.P);
iattempt = 0;
@@ -817,10 +865,9 @@ XD(sample_reinjection_step_solid_solid)
/* Reject the sampling of the re-injection step if it has already been
* defined, i.e. if the enclosure is a limit condition */
- if(!multiple_media_frt) {
-
+ if(!multi_frt) {
/* Find the reinjection ray for the front side */
- find_reinject_ray_frt_args.solid = args_frt->solid;
+ find_reinject_ray_frt_args.solid_enc_id = args_frt->solid_enc_id;
find_reinject_ray_frt_args.rwalk = args_frt->rwalk;
find_reinject_ray_frt_args.distance = args_frt->distance;
find_reinject_ray_frt_args.can_move = 1;
@@ -837,10 +884,9 @@ XD(sample_reinjection_step_solid_solid)
/* Reject the sampling of the re-injection step if it has already been
* defined, i.e. if the enclosure is a limit condition */
- if(!multiple_media_bck) {
-
+ if(!multi_bck) {
/* Select the reinjection direction and distance for the back side */
- find_reinject_ray_bck_args.solid = args_bck->solid;
+ find_reinject_ray_bck_args.solid_enc_id = args_bck->solid_enc_id;
find_reinject_ray_bck_args.rwalk = args_bck->rwalk;
find_reinject_ray_bck_args.distance = args_bck->distance;
find_reinject_ray_bck_args.can_move = 1;
@@ -856,7 +902,7 @@ XD(sample_reinjection_step_solid_solid)
/* If random walk was moved to find a valid rinjection ray on back side,
* one has to find a valid reinjection on front side from the new pos */
- if(ray_bck.position_was_moved && !multiple_media_frt) {
+ if(ray_bck.position_was_moved) {
find_reinject_ray_frt_args.can_move = 0;
res = XD(find_reinjection_ray_and_check_validity)
(scn, &find_reinject_ray_frt_args, &ray_frt);
@@ -880,17 +926,17 @@ XD(sample_reinjection_step_solid_solid)
}
/* Setup the front and back reinjection steps */
- if(!multiple_media_frt) {
+ if(!multi_frt) {
step_frt->XD(hit) = ray_frt.hit;
step_frt->distance = ray_frt.dst;
fX(set)(step_frt->direction, ray_frt.dir);
- ASSERT(XD(check_reinjection_step)(step_frt)); /* Post-condition */
+ ASSERT(XD(check_reinjection_step)(step_frt) == RES_OK); /* Post-condition */
}
- if(!multiple_media_bck) {
+ if(!multi_bck) {
step_bck->XD(hit) = ray_bck.hit;
step_bck->distance = ray_bck.dst;
fX(set)(step_bck->direction, ray_bck.dir);
- ASSERT(XD(check_reinjection_step)(step_bck)); /* Post-condition */
+ ASSERT(XD(check_reinjection_step)(step_bck) == RES_OK); /* Post-condition */
}
exit:
@@ -901,17 +947,28 @@ error:
res_T
XD(solid_reinjection)
- (struct sdis_medium* solid,
+ (struct sdis_scene* scn,
+ const unsigned solid_enc_id,
struct solid_reinjection_args* args)
{
+ /* Properties */
struct solid_props props = SOLID_PROPS_NULL;
+ struct sdis_medium* solid = NULL;
+ const struct enclosure* enc = NULL;
+
double reinject_dst_m; /* Reinjection distance in meters */
double mu;
res_T res = RES_OK;
- ASSERT(solid && XD(check_solid_reinjection_args)(args));
+ ASSERT(XD(check_solid_reinjection_args)(args) == RES_OK);
+ ASSERT(solid_enc_id != ENCLOSURE_ID_NULL);
reinject_dst_m = args->reinjection->distance * args->fp_to_meter;
+ /* Get the enclosure medium properties */
+ enc = scene_get_enclosure(scn, solid_enc_id);
+ res = scene_get_enclosure_medium(scn, enc, &solid);
+ if(res != RES_OK) goto error;
+ ASSERT(sdis_medium_get_type(solid) == SDIS_SOLID);
res = solid_get_properties(solid, &args->rwalk->vtx, &props);
if(res != RES_OK) goto error;
@@ -921,10 +978,10 @@ XD(solid_reinjection)
if(res != RES_OK) goto error;
/* Time rewind */
- args->rwalk->mdm = solid; /* Medium into which the time is rewind */
+ args->rwalk->enc_id = solid_enc_id; /* Enclosure into which the time is rewind */
mu = (2*DIM*props.lambda)/(props.rho*props.cp*reinject_dst_m*reinject_dst_m);
res = time_rewind
- (mu, props.t0, args->rng, args->rwalk, args->rwalk_ctx, args->T);
+ (scn, 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 */
@@ -939,14 +996,14 @@ XD(solid_reinjection)
/* The random walk is in the solid */
if(args->reinjection->XD(hit).distance != args->reinjection->distance) {
args->T->func = XD(conductive_path);
- args->rwalk->mdm = solid;
+ args->rwalk->enc_id = solid_enc_id;
args->rwalk->XD(hit) = SXD_HIT_NULL;
args->rwalk->hit_side = SDIS_SIDE_NULL__;
/* The random walk is at a boundary */
} else {
args->T->func = XD(boundary_path);
- args->rwalk->mdm = NULL;
+ args->rwalk->enc_id = ENCLOSURE_ID_NULL;
args->rwalk->XD(hit) = args->reinjection->XD(hit);
if(fX(dot)(args->reinjection->XD(hit).normal, args->reinjection->direction) < 0) {
args->rwalk->hit_side = SDIS_FRONT;
@@ -1035,7 +1092,7 @@ XD(check_Tref)
if((Bound) < 0) { \
log_err(scn->dev, \
"%s: the "Name" temperature cannot be negative " \
- "to sample a radiative path (T"Name" = %g K).\n", \
+ "to sample a radiative path -- T"Name" = %g K\n", \
func_name, (Bound)); \
return RES_BAD_OP_IRRECOVERABLE; \
} \
@@ -1054,7 +1111,7 @@ XD(check_Tref)
if(SDIS_TEMPERATURE_IS_UNKNOWN(Tref)) {
log_err(scn->dev,
- "%s: the reference temperature is unknown at `"FORMAT_VECX". "
+ "%s: the reference temperature is unknown at ("FORMAT_VECX"). "
"Sampling a radiative path requires a valid reference temperature field.\n",
func_name, SPLITX(pos));
return RES_BAD_OP_IRRECOVERABLE;
@@ -1062,7 +1119,7 @@ XD(check_Tref)
if(Tref < 0) {
log_err(scn->dev,
- "%s: the reference temperature is negative at `"FORMAT_VECX" (Tref = %g K). "
+ "%s: the reference temperature is negative at ("FORMAT_VECX") and Tref = %g K. "
"Sampling a radiative path requires a known, positive reference "
"temperature field.\n",
func_name, SPLITX(pos), Tref);
@@ -1071,7 +1128,7 @@ XD(check_Tref)
if(Tref < scn->tmin || scn->tmax < Tref) {
log_err(scn->dev,
- "%s: invalid reference temperature at `"FORMAT_VECX"' (Tref = %g K). "
+ "%s: invalid reference temperature at ("FORMAT_VECX") and Tref=%g K. "
"It must be included in the provided temperature range "
"(Tmin = %g K; Tmax = %g K)\n",
func_name, SPLITX(pos), Tref, scn->tmin, scn->tmax);
diff --git a/src/sdis_heat_path_boundary_Xd_solid_fluid_picard1.h b/src/sdis_heat_path_boundary_Xd_solid_fluid_picard1.h
@@ -108,6 +108,9 @@ XD(solid_fluid_boundary_picard1_path)
/* Fragment on the fluid side of the boundary */
struct sdis_interface_fragment frag_fluid;
+ /* The enclosures split by the boundary */
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+
/* Data attached to the boundary */
struct sdis_interface* interf = NULL;
struct sdis_medium* solid = NULL;
@@ -135,7 +138,7 @@ XD(solid_fluid_boundary_picard1_path)
res_T res = RES_OK;
ASSERT(scn && rwalk && rng && T && ctx);
- ASSERT(XD(check_rwalk_fragment_consistency)(rwalk, frag));
+ ASSERT(XD(check_rwalk_fragment_consistency)(rwalk, frag) == RES_OK);
/* Retrieve the solid and the fluid split by the boundary */
interf = scene_get_interface(scn, rwalk->XD(hit).prim.prim_id);
@@ -149,6 +152,9 @@ XD(solid_fluid_boundary_picard1_path)
ASSERT(fluid->type == SDIS_FLUID);
}
+ /* Get the enclosures split by the boundary */
+ scene_get_enclosure_ids(scn, rwalk->XD(hit).prim.prim_id, enc_ids);
+
/* Setup a fragment for the fluid side */
frag_fluid = *frag;
frag_fluid.side = fluid_side;
@@ -177,7 +183,7 @@ XD(solid_fluid_boundary_picard1_path)
/* Sample a reinjection step */
samp_reinject_step_args.rng = rng;
- samp_reinject_step_args.solid = solid;
+ samp_reinject_step_args.solid_enc_id = enc_ids[solid_side];
samp_reinject_step_args.rwalk = rwalk;
samp_reinject_step_args.distance = delta_boundary;
samp_reinject_step_args.side = solid_side;
@@ -247,7 +253,7 @@ XD(solid_fluid_boundary_picard1_path)
/* Switch in convective path */
if(r < p_conv) {
T->func = XD(convective_path);
- rwalk->mdm = fluid;
+ rwalk->enc_id = enc_ids[fluid_side];
rwalk->hit_side = fluid_side;
break;
}
@@ -264,7 +270,7 @@ XD(solid_fluid_boundary_picard1_path)
solid_reinject_args.rng = rng;
solid_reinject_args.T = T;
solid_reinject_args.fp_to_meter = scn->fp_to_meter;
- res = XD(solid_reinjection)(solid, &solid_reinject_args);
+ res = XD(solid_reinjection)(scn, enc_ids[solid_side], &solid_reinject_args);
if(res != RES_OK) goto error;
break;
}
@@ -281,7 +287,7 @@ XD(solid_fluid_boundary_picard1_path)
/* Sample a radiative path and get the Tref at its end. */
T_s = *T;
rwalk_s = *rwalk;
- rwalk_s.mdm = fluid;
+ rwalk_s.enc_id = enc_ids[fluid_side];
rwalk_s.hit_side = fluid_side;
res = XD(radiative_path)(scn, ctx, &rwalk_s, rng, &T_s);
if(res != RES_OK) goto error;
diff --git a/src/sdis_heat_path_boundary_Xd_solid_fluid_picardN.h b/src/sdis_heat_path_boundary_Xd_solid_fluid_picardN.h
@@ -81,7 +81,7 @@ XD(sample_path)
/* Init the random walk */
rwalk.vtx = rwalk_from->vtx;
- rwalk.mdm = rwalk_from->mdm;
+ rwalk.enc_id = rwalk_from->enc_id;
rwalk.XD(hit) = rwalk_from->XD(hit);
rwalk.hit_side = rwalk_from->hit_side;
@@ -144,6 +144,9 @@ XD(solid_fluid_boundary_picardN_path)
struct sdis_heat_vertex hvtx = SDIS_HEAT_VERTEX_NULL;
struct sdis_heat_vertex hvtx_s = SDIS_HEAT_VERTEX_NULL;
+ /* The enclosures split by the boundary */
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+
/* Data attached to the boundary */
struct sdis_interface* interf = NULL;
struct sdis_medium* solid = NULL;
@@ -171,8 +174,7 @@ XD(solid_fluid_boundary_picardN_path)
res_T res = RES_OK;
ASSERT(scn && rwalk && rng && T && ctx);
- ASSERT(XD(check_rwalk_fragment_consistency)(rwalk, frag));
-
+ ASSERT(XD(check_rwalk_fragment_consistency)(rwalk, frag) == RES_OK);
/* Fetch the Min/max temperature */
Tmin = ctx->Tmin;
@@ -194,6 +196,9 @@ XD(solid_fluid_boundary_picardN_path)
ASSERT(fluid->type == SDIS_FLUID);
}
+ /* Get the enclosures split by the boundary */
+ scene_get_enclosure_ids(scn, rwalk->XD(hit).prim.prim_id, enc_ids);
+
/* Check that no net flux is set for this interface since the provided
* picardN algorithm does not handle it */
res = check_net_flux(scn, interf, frag, get_picard_order(ctx));
@@ -218,7 +223,7 @@ XD(solid_fluid_boundary_picardN_path)
/* Sample a reinjection step */
samp_reinject_step_args.rng = rng;
- samp_reinject_step_args.solid = solid;
+ samp_reinject_step_args.solid_enc_id = enc_ids[solid_side];
samp_reinject_step_args.rwalk = rwalk;
samp_reinject_step_args.distance = delta_boundary;
samp_reinject_step_args.side = solid_side;
@@ -271,7 +276,7 @@ XD(solid_fluid_boundary_picardN_path)
/* Switch in convective path */
if(r < p_conv) {
T->func = XD(convective_path);
- rwalk->mdm = fluid;
+ rwalk->enc_id = enc_ids[fluid_side];
rwalk->hit_side = fluid_side;
break;
}
@@ -288,7 +293,7 @@ XD(solid_fluid_boundary_picardN_path)
solid_reinject_args.rng = rng;
solid_reinject_args.T = T;
solid_reinject_args.fp_to_meter = scn->fp_to_meter;
- res = XD(solid_reinjection)(solid, &solid_reinject_args);
+ res = XD(solid_reinjection)(scn, enc_ids[solid_side], &solid_reinject_args);
if(res != RES_OK) goto error;
break;
}
@@ -302,7 +307,7 @@ XD(solid_fluid_boundary_picardN_path)
/* Sample a radiative path */
T_s = *T;
rwalk_s = *rwalk;
- rwalk_s.mdm = fluid;
+ rwalk_s.enc_id = enc_ids[fluid_side];
rwalk_s.hit_side = fluid_side;
res = XD(radiative_path)(scn, ctx, &rwalk_s, rng, &T_s);
if(res != RES_OK) goto error;
diff --git a/src/sdis_heat_path_boundary_Xd_solid_solid.h b/src/sdis_heat_path_boundary_Xd_solid_solid.h
@@ -51,10 +51,11 @@ XD(solid_solid_boundary_path)
SOLID_REINJECTION_ARGS_NULL;
/* Data attached to the boundary */
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
struct sdis_interface* interf = NULL;
struct sdis_medium* solid_frt = NULL;
struct sdis_medium* solid_bck = NULL;
- struct sdis_medium* solid = NULL;
+ unsigned solid_enc_id = ENCLOSURE_ID_NULL; /* Solid to re-inject */
double lambda_frt;
double lambda_bck;
@@ -67,10 +68,11 @@ XD(solid_solid_boundary_path)
res_T res = RES_OK;
ASSERT(scn && ctx && frag && rwalk && rng && T);
- ASSERT(XD(check_rwalk_fragment_consistency)(rwalk, frag));
+ ASSERT(XD(check_rwalk_fragment_consistency)(rwalk, frag) == RES_OK);
(void)frag, (void)ctx;
- /* Retrieve the two solids split by the boundary */
+ /* Retrieve the two enclosures and associated media split by the boundary */
+ scene_get_enclosure_ids(scn, rwalk->XD(hit).prim.prim_id, enc_ids);
interf = scene_get_interface(scn, rwalk->XD(hit).prim.prim_id);
solid_frt = interface_get_medium(interf, SDIS_FRONT);
solid_bck = interface_get_medium(interf, SDIS_BACK);
@@ -93,8 +95,8 @@ XD(solid_solid_boundary_path)
/* Sample a front/back reinjection steps */
samp_reinject_step_frt_args.rng = rng;
samp_reinject_step_bck_args.rng = rng;
- samp_reinject_step_frt_args.solid = solid_frt;
- samp_reinject_step_bck_args.solid = solid_bck;
+ samp_reinject_step_frt_args.solid_enc_id = enc_ids[SDIS_FRONT];
+ samp_reinject_step_bck_args.solid_enc_id = enc_ids[SDIS_BACK];
samp_reinject_step_frt_args.rwalk = rwalk;
samp_reinject_step_bck_args.rwalk = rwalk;
samp_reinject_step_frt_args.distance = delta_boundary_frt;
@@ -150,10 +152,10 @@ XD(solid_solid_boundary_path)
if(r < proba) { /* Reinject in front */
reinject_step = &reinject_step_frt;
- solid = solid_frt;
+ solid_enc_id = enc_ids[SDIS_FRONT];
} else { /* Reinject in back */
reinject_step = &reinject_step_bck;
- solid = solid_bck;
+ solid_enc_id = enc_ids[SDIS_BACK];
}
/* Perform the reinjection into the solid */
@@ -163,7 +165,7 @@ XD(solid_solid_boundary_path)
solid_reinject_args.rwalk_ctx = ctx;
solid_reinject_args.T = T;
solid_reinject_args.fp_to_meter = scn->fp_to_meter;
- res = XD(solid_reinjection)(solid, &solid_reinject_args);
+ res = XD(solid_reinjection)(scn, solid_enc_id, &solid_reinject_args);
if(res != RES_OK) goto error;
exit:
diff --git a/src/sdis_heat_path_boundary_c.h b/src/sdis_heat_path_boundary_c.h
@@ -30,14 +30,14 @@ struct sdis_medium;
******************************************************************************/
struct sample_reinjection_step_args {
struct ssp_rng* rng; /* Random number generator to use */
- const struct sdis_medium* solid; /* Solid in which to reinject */
struct rwalk* rwalk; /* Current state of the random walk */
double distance; /* Maximum Reinjection distance */
+ unsigned solid_enc_id; /* Enclosured Id of the solid in which to reinject */
enum sdis_side side; /* Side of the boundary to re-inject */
};
#define SAMPLE_REINJECTION_STEP_ARGS_NULL__ \
- {NULL, NULL, NULL, -1, SDIS_SIDE_NULL__}
+ {NULL, NULL, -1, ENCLOSURE_ID_NULL, SDIS_SIDE_NULL__}
static const struct sample_reinjection_step_args
SAMPLE_REINJECTION_STEP_ARGS_NULL = SAMPLE_REINJECTION_STEP_ARGS_NULL__;
@@ -98,12 +98,14 @@ static const struct solid_reinjection_args SOLID_REINJECTION_ARGS_NULL =
extern LOCAL_SYM res_T
solid_reinjection_2d
- (struct sdis_medium* solid,
+ (struct sdis_scene* scn,
+ const unsigned solid_enc_id,
struct solid_reinjection_args* args);
extern LOCAL_SYM res_T
solid_reinjection_3d
- (struct sdis_medium* solid,
+ (struct sdis_scene* scn,
+ const unsigned solid_enc_id,
struct solid_reinjection_args* args);
/*******************************************************************************
diff --git a/src/sdis_heat_path_conductive_delta_sphere_Xd.h b/src/sdis_heat_path_conductive_delta_sphere_Xd.h
@@ -129,7 +129,6 @@ XD(sample_next_step_robust)
res_T res = RES_OK;
ASSERT(scn && rng && pos && delta_solid > 0);
ASSERT(current_enc_id != ENCLOSURE_ID_NULL);
- ASSERT(current_enc_id != ENCLOSURE_ID_MULTI_MEDIA);
ASSERT(dir0 && dir1 && hit0 && hit1 && out_delta);
fX_set_dX(org, pos);
@@ -155,30 +154,18 @@ XD(sample_next_step_robust)
/* Check medium consistency */
if(current_enc_id != enc_id) {
#if 0
-#if DIM == 2
- log_warn(scn->dev,
- "%s: inconsistent medium during the solid random walk at {%g, %g}.\n",
- FUNC_NAME, SPLIT2(pos));
-#else
log_warn(scn->dev,
- "%s: inconsistent medium during the solid random walk at {%g, %g, %g}.\n",
- FUNC_NAME, SPLIT3(pos));
-#endif
+ "%s: inconsistent medium during the solid random walk -- "
+ "pos=("FORMAT_VECX")\n", FUNC_NAME, SPLITX(pos));
#endif
}
} while(current_enc_id != enc_id && ++iattempt < MAX_ATTEMPTS);
/* Handle error */
if(iattempt >= MAX_ATTEMPTS) {
-#if DIM == 2
- log_warn(scn->dev,
- "%s: could not find a next valid conductive step at {%g, %g}.\n",
- FUNC_NAME, SPLIT2(pos));
-#else
log_warn(scn->dev,
- "%s: could not find a next valid conductive step at {%g, %g, %g}.\n",
- FUNC_NAME, SPLIT3(pos));
-#endif
+ "%s: could not find a next valid conductive -- pos=("FORMAT_VECX")\n",
+ FUNC_NAME, SPLITX(pos));
res = RES_BAD_OP;
goto error;
}
@@ -351,7 +338,6 @@ XD(conductive_path_delta_sphere)
/* Check pre-conditions */
ASSERT(scn && rwalk && rng && T);
- ASSERT(rwalk->mdm->type == SDIS_SOLID);
(void)ctx, (void)istep; /* Avoid "unsued variable" warnings */
@@ -359,11 +345,13 @@ XD(conductive_path_delta_sphere)
if(res != RES_OK) goto error;
res = scene_get_enclosure_medium(scn, scene_get_enclosure(scn, enc_id), &mdm);
if(res != RES_OK) goto error;
+ ASSERT(sdis_medium_get_type(mdm) == SDIS_SOLID);
/* Check the random walk consistency */
- if(mdm != rwalk->mdm) {
+ if(enc_id != rwalk->enc_id) {
log_err(scn->dev, "%s: invalid solid random walk. "
- "Unexpected medium at {%g, %g, %g}.\n", FUNC_NAME, SPLIT3(rwalk->vtx.P));
+ "Unexpected enclosure -- pos=("FORMAT_VECX")\n",
+ FUNC_NAME, SPLITX(rwalk->vtx.P));
res = RES_BAD_OP_IRRECOVERABLE;
goto error;
}
@@ -410,7 +398,7 @@ XD(conductive_path_delta_sphere)
if(ctx->green_path) {
res = green_path_set_limit_vertex
- (ctx->green_path, rwalk->mdm, &rwalk->vtx, rwalk->elapsed_time);
+ (ctx->green_path, mdm, &rwalk->vtx, rwalk->elapsed_time);
if(res != RES_OK) goto error;
}
@@ -450,7 +438,7 @@ XD(conductive_path_delta_sphere)
/* Rewind the time */
delta_m = delta * scn->fp_to_meter;
mu = (2*DIM*props.lambda)/(props.rho*props.cp*delta_m*delta_m);
- res = time_rewind(mu, props.t0, rng, rwalk, ctx, T);
+ res = time_rewind(scn, mu, props.t0, rng, rwalk, ctx, T);
if(res != RES_OK) goto error;
if(T->done) break; /* Limit condition was reached */
@@ -479,12 +467,12 @@ XD(conductive_path_delta_sphere)
/* Register the power term for the green function */
if(ctx->green_path && props_ref.power != SDIS_VOLUMIC_POWER_NONE) {
res = green_path_add_power_term
- (ctx->green_path, rwalk->mdm, &rwalk->vtx, green_power_term);
+ (ctx->green_path, mdm, &rwalk->vtx, green_power_term);
if(res != RES_OK) goto error;
}
T->func = XD(boundary_path);
- rwalk->mdm = NULL; /* The random walk is at an interface between 2 media */
+ rwalk->enc_id = ENCLOSURE_ID_NULL; /* At the interface between 2 media */
exit:
return res;
diff --git a/src/sdis_heat_path_conductive_wos_Xd.h b/src/sdis_heat_path_conductive_wos_Xd.h
@@ -69,24 +69,21 @@ error:
* Helper function
******************************************************************************/
static res_T
-XD(check_medium_consistency)
+XD(check_enclosure_consistency)
(struct sdis_scene* scn,
const struct rwalk* rwalk)
{
unsigned enc_id = ENCLOSURE_ID_NULL;
- struct sdis_medium* mdm = NULL;
res_T res = RES_OK;
ASSERT(rwalk);
res = scene_get_enclosure_id_in_closed_boundaries(scn, rwalk->vtx.P, &enc_id);
if(res != RES_OK) goto error;
- res = scene_get_enclosure_medium(scn, scene_get_enclosure(scn, enc_id), &mdm);
- if(res != RES_OK) goto error;
- /* Check medium consistency */
- if(mdm != rwalk->mdm) {
+ /* Check enclosure consistency */
+ if(enc_id != rwalk->enc_id) {
log_err(scn->dev,
- "%s:%s: invalid solid walk. Unexpected medium (position: "FORMAT_VECX").\n",
+ "%s:%s: invalid solid walk. Unexpected enclosure -- pos=("FORMAT_VECX")\n",
__FILE__, FUNC_NAME, SPLITX(rwalk->vtx.P));
res = RES_BAD_OP_IRRECOVERABLE;
goto error;
@@ -103,6 +100,7 @@ XD(time_travel)
(struct sdis_scene* scn,
struct rwalk* rwalk,
struct ssp_rng* rng,
+ struct sdis_medium* mdm,
const double alpha, /* Diffusivity, i.e. lambda/(rho*cp) */
const double t0, /* Initial time [s] */
double* distance, /* Displacement [m/fp_to_meter] */
@@ -166,14 +164,14 @@ XD(time_travel)
dX(add)(rwalk->vtx.P, rwalk->vtx.P, dir);
/* Fetch the initial temperature */
- temperature = medium_get_temperature(rwalk->mdm, &rwalk->vtx);
+ temperature = medium_get_temperature(mdm, &rwalk->vtx);
if(SDIS_TEMPERATURE_IS_UNKNOWN(temperature)) {
log_err(scn->dev,
"%s:%s: the path reaches the initial condition but the "
- "%s temperature remains unknown -- position=%g, %g, %g\n",
+ "%s temperature remains unknown -- pos=("FORMAT_VECX")\n",
__FILE__, FUNC_NAME,
- medium_type_to_string(sdis_medium_get_type(rwalk->mdm)),
- SPLIT3(rwalk->vtx.P));
+ medium_type_to_string(sdis_medium_get_type(mdm)),
+ SPLITX(rwalk->vtx.P));
res = RES_BAD_ARG;
goto error;
}
@@ -275,16 +273,15 @@ compute_hit_side_3d
}
#endif
-/* Verify that the submitted position is in the expected medium */
+/* Verify that the submitted position is in the expected enclosure */
static res_T
XD(check_diffusion_position)
(struct sdis_scene* scn,
- const struct sdis_medium* expected_medium,
+ const unsigned expected_enc_id,
const double pos[DIM])
{
- struct sdis_interface* interf = NULL;
- struct sdis_medium* mdm = NULL;
- enum sdis_side side;
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+ enum sdis_side side = SDIS_SIDE_NULL__;
struct sXd(hit) hit = SXD_HIT_NULL;
float wos_pos[DIM] = {0};
@@ -292,7 +289,8 @@ XD(check_diffusion_position)
res_T res = RES_OK;
/* Check pre-conditions */
- ASSERT(scn && expected_medium && pos);
+ ASSERT(scn && pos);
+ ASSERT(expected_enc_id != ENCLOSURE_ID_NULL);
/* Look for the nearest surface within 1 mm of the position to be checked. By
* limiting the search radius we speed up the closest point query. If no
@@ -311,11 +309,10 @@ XD(check_diffusion_position)
SXD(scene_view_closest_point(scn->sXd(view), wos_pos, wos_radius, NULL, &hit));
if(SXD_HIT_NONE(&hit)) goto exit;
- /* Fetch interface properties and check path consistency */
- interf = scene_get_interface(scn, hit.prim.prim_id);
+ /* Check path consistency */
+ scene_get_enclosure_ids(scn, hit.prim.prim_id, enc_ids);
side = XD(compute_hit_side)(&hit, pos);
- mdm = side == SDIS_FRONT ? interf->medium_front : interf->medium_back;
- if(mdm != expected_medium) {
+ if(enc_ids[side] != expected_enc_id) {
res = RES_BAD_ARG;
goto error;
}
@@ -337,8 +334,7 @@ XD(setup_hit_wos)
struct sXd(attrib) attr;
/* Properties */
- struct sdis_interface* interf = NULL;
- struct sdis_medium* mdm = NULL;
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
enum sdis_side side = SDIS_SIDE_NULL__;
/* Miscellaneous */
@@ -360,13 +356,12 @@ XD(setup_hit_wos)
dX_set_fX(tgt, attr.value);
side = XD(compute_hit_side)(hit, rwalk->vtx.P);
- /* Fetch interface properties and check path consistency */
- interf = scene_get_interface(scn, hit->prim.prim_id);
- mdm = side == SDIS_FRONT ? interf->medium_front : interf->medium_back;
- if(mdm != rwalk->mdm) {
+ /* Check path consistency */
+ scene_get_enclosure_ids(scn, hit->prim.prim_id, enc_ids);
+ if(enc_ids[side] != rwalk->enc_id) {
log_err(scn->dev,
- "%s:%s: the conductive path has reached an invalid interface; "
- "unexpected medium (position: "FORMAT_VECX"; side: %s).\n",
+ "%s:%s: the conductive path has reached an invalid interface. "
+ "Unexpected enclosure -- pos=("FORMAT_VECX"), side=%s\n",
__FILE__, FUNC_NAME, SPLITX(tgt), side == SDIS_FRONT ? "front" : "back");
res = RES_BAD_OP_IRRECOVERABLE;
goto error;
@@ -395,8 +390,7 @@ XD(setup_hit_rt)
struct rwalk* rwalk)
{
/* Properties */
- struct sdis_interface* interf = NULL;
- struct sdis_medium* mdm = NULL;
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
enum sdis_side side = SDIS_SIDE_NULL__;
/* Miscellaneous */
@@ -416,12 +410,11 @@ XD(setup_hit_rt)
side = dX(dot)(N, dir) > 0 ? SDIS_BACK : SDIS_FRONT;
/* Fetch interface properties and check path consistency */
- interf = scene_get_interface(scn, hit->prim.prim_id);
- mdm = side == SDIS_FRONT ? interf->medium_front : interf->medium_back;
- if(mdm != rwalk->mdm) {
+ scene_get_enclosure_ids(scn, hit->prim.prim_id, enc_ids);
+ if(enc_ids[side] != rwalk->enc_id) {
log_err(scn->dev,
- "%s:%s: the conductive path has reached an invalid interface; "
- "unexpected medium (position: "FORMAT_VECX"; side: %s).\n",
+ "%s:%s: the conductive path has reached an invalid interface. "
+ "Unexpected enclosure -- pos=("FORMAT_VECX"), side=%s\n",
__FILE__, FUNC_NAME, SPLITX(tgt), side == SDIS_FRONT ? "front" : "back");
res = RES_BAD_OP_IRRECOVERABLE;
goto error;
@@ -497,7 +490,7 @@ XD(sample_next_position)
* The next diffusion step would then detect an error. This is why we use a
* new function based on the same geometric operator used in the present
* algorithm. */
- res = XD(check_diffusion_position)(scn, rwalk->mdm, pos);
+ res = XD(check_diffusion_position)(scn, rwalk->enc_id, pos);
/* Diffusion position is valid => move the path to the new position */
if(res == RES_OK) {
@@ -531,8 +524,8 @@ XD(sample_next_position)
* we don't care to save it. */
if(SXD_HIT_NONE(&hit)) {
log_err(scn->dev,
- "%s:%s: unable to find the next diffusion position "
- "(position: "FORMAT_VECX"; direction: "FORMAT_VECX"; distance: %g\n",
+ "%s:%s: unable to find the next diffusion position -- "
+ "position=("FORMAT_VECX"), direction=("FORMAT_VECX"), distance=%g\n",
__FILE__, FUNC_NAME, SPLITX(pos), SPLITX(dir), wos_distance);
res = RES_BAD_OP_IRRECOVERABLE;
goto error;
@@ -562,6 +555,7 @@ XD(conductive_path_wos)
struct temperature* T)
{
/* Properties */
+ const struct enclosure* enc = NULL;
struct sdis_medium* mdm = NULL;
struct solid_props props_ref = SOLID_PROPS_NULL;
struct solid_props props = SOLID_PROPS_NULL;
@@ -577,18 +571,18 @@ XD(conductive_path_wos)
/* Check pre-conditions */
ASSERT(scn && ctx && rwalk && rng && T);
- ASSERT(sdis_medium_get_type(rwalk->mdm) == SDIS_SOLID);
/* Is green evaluated evaluated */
green = ctx->green_path != NULL;
- /* Keep track of the solid. After conduction, a boundary may have been
- * reached, so the random walk medium is NULL. However, this medium is still
- * needed to update the green path. Hence this backup */
- mdm = rwalk->mdm;
+ res = XD(check_enclosure_consistency)(scn, rwalk);
+ if(res != RES_OK) goto error;
- res = XD(check_medium_consistency)(scn, rwalk);
+ /* Get the enclosure medium */
+ enc = scene_get_enclosure(scn, rwalk->enc_id);
+ res = scene_get_enclosure_medium(scn, enc, &mdm);
if(res != RES_OK) goto error;
+ ASSERT(sdis_medium_get_type(mdm) == SDIS_SOLID);
/* Retrieve the solid properties at the current position. Use them to verify
* that those that are supposed to be constant by the conductive random walk
@@ -599,7 +593,7 @@ XD(conductive_path_wos)
* position. By comparing them to the properties along the random walk, we
* thus verify that the properties are constant throughout the random walk
* with respect to the properties of the reinjected position. */
- solid_get_properties(rwalk->mdm, &rwalk->vtx, &props_ref);
+ solid_get_properties(mdm, &rwalk->vtx, &props_ref);
props = props_ref;
/* The algorithm assumes that lambda, rho and cp are constants. The
@@ -631,7 +625,7 @@ XD(conductive_path_wos)
if(res != RES_OK) goto error;
/* Going back in time */
- res = XD(time_travel)(scn, rwalk, rng, alpha, props.t0, &dst, T);
+ res = XD(time_travel)(scn, rwalk, rng, mdm, alpha, props.t0, &dst, T);
if(res != RES_OK) goto error;
/* Add the volumic power density */
@@ -652,14 +646,14 @@ XD(conductive_path_wos)
/* The path reaches a boundary */
if(!SXD_HIT_NONE(&rwalk->XD(hit))) {
T->func = XD(boundary_path);
- rwalk->mdm = NULL;
+ rwalk->enc_id = ENCLOSURE_ID_NULL;
break;
}
#undef REGISTER_VERTEX
/* Retreive and check solid properties at the new position */
- res = solid_get_properties(rwalk->mdm, &rwalk->vtx, &props);
+ res = solid_get_properties(mdm, &rwalk->vtx, &props);
if(res != RES_OK) goto error;
res = check_solid_constant_properties(scn->dev, green, wos, &props_ref, &props);
if(res != RES_OK) goto error;
diff --git a/src/sdis_heat_path_convective_Xd.h b/src/sdis_heat_path_convective_Xd.h
@@ -66,57 +66,20 @@ error:
* Helper functions
******************************************************************************/
static res_T
-XD(register_heat_vertex_in_fluid)
- (struct sdis_scene* scn,
- struct rwalk_context* ctx,
- struct rwalk* rwalk,
- const double weight)
-{
- struct sdis_rwalk_vertex vtx = SDIS_RWALK_VERTEX_NULL;
- struct hit_filter_data filter_data = HIT_FILTER_DATA_NULL;
- const float empirical_dst = 0.1f;
- const float range[2] = {0, FLT_MAX};
- float org[DIM];
- float dir[DIM];
- float pos[DIM];
- float dst;
- struct sXd(hit) hit;
-
- if(!ctx->heat_path) return RES_OK;
-
- ASSERT(!SXD_HIT_NONE(&rwalk->XD(hit)));
-
- fX_set_dX(org, rwalk->vtx.P);
- fX(set)(dir, rwalk->XD(hit).normal);
- if(rwalk->hit_side == SDIS_BACK) fX(minus)(dir, dir);
-
- filter_data.XD(hit) = rwalk->XD(hit);
- filter_data.epsilon = 1.e-6;
- SXD(scene_view_trace_ray(scn->sXd(view), org, dir, range, &filter_data, &hit));
- dst = SXD_HIT_NONE(&hit) ? empirical_dst : hit.distance * 0.5f;
-
- vtx = rwalk->vtx;
- fX(add)(pos, org, fX(mulf)(dir, dir, dst));
- dX_set_fX(vtx.P, pos);
-
- return register_heat_vertex(ctx->heat_path, &vtx, weight,
- SDIS_HEAT_VERTEX_CONVECTION, (int)ctx->nbranchings);
-}
-
-static res_T
XD(handle_known_fluid_temperature)
(struct sdis_scene* scn,
struct rwalk_context* ctx,
struct rwalk* rwalk,
+ struct sdis_medium* mdm,
struct temperature* T)
{
double temperature;
int known_temperature;
res_T res = RES_OK;
ASSERT(scn && ctx && rwalk && T);
- ASSERT(sdis_medium_get_type(rwalk->mdm) == SDIS_FLUID);
+ ASSERT(sdis_medium_get_type(mdm) == SDIS_FLUID);
- temperature = fluid_get_temperature(rwalk->mdm, &rwalk->vtx);
+ temperature = fluid_get_temperature(mdm, &rwalk->vtx);
/* Check if the temperature is known */
known_temperature = SDIS_TEMPERATURE_IS_KNOWN(temperature);
@@ -127,12 +90,13 @@ XD(handle_known_fluid_temperature)
if(ctx->green_path) {
res = green_path_set_limit_vertex
- (ctx->green_path, rwalk->mdm, &rwalk->vtx, rwalk->elapsed_time);
+ (ctx->green_path, mdm, &rwalk->vtx, rwalk->elapsed_time);
if(res != RES_OK) goto error;
}
- res = XD(register_heat_vertex_in_fluid)(scn, ctx, rwalk, T->value);
- if(res != RES_OK) goto error;
+ if(ctx->heat_path) {
+ heat_path_get_last_vertex(ctx->heat_path)->weight = T->value;
+ }
exit:
return res;
@@ -151,7 +115,6 @@ XD(handle_convective_path_startup)
float org[DIM] = {0};
res_T res = RES_OK;
ASSERT(scn && rwalk && path_starts_in_fluid);
- ASSERT(sdis_medium_get_type(rwalk->mdm) == SDIS_FLUID);
*path_starts_in_fluid = SXD_HIT_NONE(&rwalk->XD(hit));
if(*path_starts_in_fluid == 0) goto exit; /* Nothing to do */
@@ -180,54 +143,28 @@ error:
}
static res_T
-XD(fetch_fluid_enclosure)
+XD(check_enclosure)
(struct sdis_scene* scn,
- struct rwalk* rwalk,
- const struct enclosure** out_enclosure)
+ const struct rwalk* rwalk,
+ const struct enclosure* enc)
{
- const struct sdis_interface* interf;
- const struct enclosure* enc;
- unsigned enc_ids[2];
- unsigned enc_id;
res_T res = RES_OK;
- ASSERT(scn && rwalk && out_enclosure);
- ASSERT(sdis_medium_get_type(rwalk->mdm) == SDIS_FLUID);
- ASSERT(!SXD_HIT_NONE(&rwalk->XD(hit)));
-
- /* Fetch the current interface and its associated enclosures */
- interf = scene_get_interface(scn, rwalk->XD(hit).prim.prim_id);
- scene_get_enclosure_ids(scn, rwalk->XD(hit).prim.prim_id, enc_ids);
-
- /* Find the enclosure identifier of the current medium */
- ASSERT(interf->medium_front != interf->medium_back);
- if(rwalk->mdm == interf->medium_front) {
- enc_id = enc_ids[0];
- ASSERT(rwalk->hit_side == SDIS_FRONT);
- } else {
- ASSERT(rwalk->mdm == interf->medium_back);
- enc_id = enc_ids[1];
- ASSERT(rwalk->hit_side == SDIS_BACK);
- }
+ ASSERT(scn && rwalk && enc);
- /* Fetch the enclosure data */
- enc = scene_get_enclosure(scn, enc_id);
- ASSERT(enc != NULL);
- if(enc->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
+ if(enc->medium_id == MEDIUM_ID_MULTI) {
/* The enclosures with multiple media are used to describe limit
* conditions and therefore they cannot be fetched */
log_err(scn->dev,
- "%s: enclosure with multiple media at {%g, %g, %g}. "
- "Path should be reached a limit condition before.\n",
- FUNC_NAME, rwalk->vtx.P[0], rwalk->vtx.P[1], DIM==3 ? rwalk->vtx.P[2]:0);
+ "%s: enclosure with multiple media at ("FORMAT_VECX"). "
+ "The path should be reached a limit condition before.\n",
+ FUNC_NAME, SPLITX(rwalk->vtx.P));
res = RES_BAD_ARG;
goto error;
}
exit:
- *out_enclosure = enc;
return res;
error:
- enc = NULL;
goto exit;
}
@@ -242,10 +179,17 @@ XD(convective_path)
struct ssp_rng* rng,
struct temperature* T)
{
- struct sXd(attrib) attr_P, attr_N;
+ /* Properties */
struct fluid_props props_ref = FLUID_PROPS_NULL;
const struct sdis_interface* interf = NULL;
+ struct sdis_medium* mdm = NULL;
+
+ /* Enclosure */
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
const struct enclosure* enc = NULL;
+
+ /* Miscellaneous */
+ struct sXd(attrib) attr_P, attr_N;
struct sXd(hit)* rwalk_hit = NULL;
double r;
#if SDIS_XD_DIMENSION == 2
@@ -255,13 +199,19 @@ XD(convective_path)
#endif
int path_starts_in_fluid;
res_T res = RES_OK;
- (void)rng, (void)ctx;
+
ASSERT(scn && ctx && rwalk && rng && T);
- ASSERT(rwalk->mdm->type == SDIS_FLUID);
+ (void)rng, (void)ctx; /* Avoid "unsued variable" warnings */
rwalk_hit = &rwalk->XD(hit);
- res = XD(handle_known_fluid_temperature)(scn, ctx, rwalk, T);
+ /* Get the enclosure medium */
+ enc = scene_get_enclosure(scn, rwalk->enc_id);
+ if((res = XD(check_enclosure)(scn, rwalk, enc)) != RES_OK) goto error;
+ if((res = scene_get_enclosure_medium(scn, enc, &mdm)) != RES_OK) goto error;
+ ASSERT(sdis_medium_get_type(mdm) == SDIS_FLUID);
+
+ res = XD(handle_known_fluid_temperature)(scn, ctx, rwalk, mdm, T);
if(res != RES_OK) goto error;
if(T->done) goto exit; /* The fluid temperature is known */
@@ -270,13 +220,10 @@ XD(convective_path)
res = XD(handle_convective_path_startup)(scn, rwalk, &path_starts_in_fluid);
if(res != RES_OK) goto error;
- res = XD(fetch_fluid_enclosure)(scn, rwalk, &enc);
- if(res != RES_OK) goto error;
-
/* Retrieve the fluid properties at the current position. Use them to verify
* that those that are supposed to be constant by the convective random walk
* remain the same. */
- res = fluid_get_properties(rwalk->mdm, &rwalk->vtx, &props_ref);
+ res = fluid_get_properties(mdm, &rwalk->vtx, &props_ref);
if(res != RES_OK) goto error;
/* The hc upper bound can be 0 if h is uniformly 0. In that case the result
@@ -284,7 +231,7 @@ XD(convective_path)
if(enc->hc_upper_bound == 0) {
ASSERT(path_starts_in_fluid); /* Cannot be in the fluid without starting there. */
rwalk->vtx.time = props_ref.t0;
- res = XD(handle_known_fluid_temperature)(scn, ctx, rwalk, T);
+ res = XD(handle_known_fluid_temperature)(scn, ctx, rwalk, mdm, T);
if(res != RES_OK) goto error;
if(T->done) {
goto exit; /* Stop the random walk */
@@ -304,7 +251,7 @@ XD(convective_path)
double mu;
/* Fetch fluid properties */
- res = fluid_get_properties(rwalk->mdm, &rwalk->vtx, &props);
+ res = fluid_get_properties(mdm, &rwalk->vtx, &props);
if(res != RES_OK) goto error;
res = check_fluid_constant_properties(scn->dev, &props_ref, &props);
@@ -312,7 +259,7 @@ XD(convective_path)
/* Sample the time using the upper bound. */
mu = enc->hc_upper_bound / (props.rho * props.cp) * enc->S_over_V;
- res = time_rewind(mu, props.t0, rng, rwalk, ctx, T);
+ res = time_rewind(scn, mu, props.t0, rng, rwalk, ctx, T);
if(res != RES_OK) goto error;
if(T->done) break; /* Limit condition was reached */
@@ -343,16 +290,19 @@ XD(convective_path)
dX_set_fX(rwalk->vtx.P, attr_P.value);
fX(set)(rwalk_hit->normal, attr_N.value);
- /* Fetch the interface of the sampled point. */
- interf = scene_get_interface(scn, rwalk_hit->prim.prim_id);
- if(rwalk->mdm == interf->medium_front) {
- rwalk->hit_side = SDIS_FRONT;
- } else if(rwalk->mdm == interf->medium_back) {
+ /* Define the interface side */
+ scene_get_enclosure_ids(scn, rwalk_hit->prim.prim_id, enc_ids);
+ if(rwalk->enc_id == enc_ids[SDIS_BACK]) {
rwalk->hit_side = SDIS_BACK;
+ } else if(rwalk->enc_id == enc_ids[SDIS_FRONT]) {
+ rwalk->hit_side = SDIS_FRONT;
} else {
- FATAL("Unexpected fluid interface.\n");
+ FATAL("Unexpected fluid interface\n");
}
+ /* Get the interface properties */
+ interf = scene_get_interface(scn, rwalk_hit->prim.prim_id);
+
/* Register the new vertex against the heat path */
res = register_heat_vertex(ctx->heat_path, &rwalk->vtx, T->value,
SDIS_HEAT_VERTEX_CONVECTION, (int)ctx->nbranchings);
@@ -380,7 +330,7 @@ XD(convective_path)
rwalk_hit->distance = 0;
T->func = XD(boundary_path);
- rwalk->mdm = NULL; /* The random walk is at an interface between 2 media */
+ rwalk->enc_id = ENCLOSURE_ID_NULL; /* Interface between 2 enclosures */
exit:
return res;
diff --git a/src/sdis_heat_path_radiative_Xd.h b/src/sdis_heat_path_radiative_Xd.h
@@ -58,7 +58,8 @@ XD(trace_radiative_path)
const struct sdis_interface* interf = NULL;
struct hit_filter_data filter_data = HIT_FILTER_DATA_NULL;
struct sdis_interface_fragment frag = SDIS_INTERFACE_FRAGMENT_NULL;
- struct sdis_medium* chk_mdm = NULL;
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+ unsigned chk_enc_id = ENCLOSURE_ID_NULL;
double alpha;
double epsilon;
double r;
@@ -156,35 +157,27 @@ XD(trace_radiative_path)
r = ssp_rng_canonical(rng);
if(r < epsilon) {
T->func = XD(boundary_path);
- rwalk->mdm = NULL; /* The random walk is at an interface between 2 media */
+ rwalk->enc_id = ENCLOSURE_ID_NULL; /* Interface between 2 enclosures */
break;
}
/* Normalize the normal of the interface and ensure that it points toward the
* current medium */
fX(normalize)(N, rwalk->XD(hit).normal);
- if(rwalk->hit_side == SDIS_BACK){
- chk_mdm = interf->medium_back;
- fX(minus)(N, N);
- } else {
- chk_mdm = interf->medium_front;
- }
+ if(rwalk->hit_side == SDIS_BACK) fX(minus)(N, N);
- if(chk_mdm != rwalk->mdm) {
- /* To ease the setting of models, the external enclosure is allowed to be
- * incoherent regarding media. Here a radiative path is allowed to join
- * 2 different fluids. */
- const int outside = scene_is_outside
- (scn, rwalk->hit_side, rwalk->XD(hit).prim.prim_id);
- if(outside && chk_mdm->type == SDIS_FLUID) {
- rwalk->mdm = chk_mdm;
- } else {
- log_warn(scn->dev, "%s: inconsistent medium definition at `%g %g %g'.\n",
- FUNC_NAME, SPLIT3(rwalk->vtx.P));
- res = RES_BAD_OP;
- goto error;
- }
+ /* Check that the radiative path still lies in the same enclosure */
+ scene_get_enclosure_ids(scn, rwalk->XD(hit).prim.prim_id, enc_ids);
+ chk_enc_id = rwalk->hit_side == SDIS_FRONT ? enc_ids[0] : enc_ids[1];
+
+ if(chk_enc_id != rwalk->enc_id) {
+ log_warn(scn->dev,
+ "%s: the radiative path has escaped from its cavity -- pos=(%g, %g, %g)\n",
+ FUNC_NAME, SPLIT3(rwalk->vtx.P));
+ res = RES_BAD_OP;
+ goto error;
}
+
alpha = interface_side_get_specular_fraction(interf, SDIS_INTERN_SOURCE_ID, &frag);
r = ssp_rng_canonical(rng);
if(r < alpha) { /* Sample specular part */
diff --git a/src/sdis_misc.c b/src/sdis_misc.c
@@ -24,17 +24,25 @@
res_T
time_rewind
- (const double mu,
+ (struct sdis_scene* scn,
+ const double mu,
const double t0,
struct ssp_rng* rng,
struct rwalk* rwalk,
const struct rwalk_context* ctx,
struct temperature* T)
{
+ const struct enclosure* enc = NULL;
+ struct sdis_medium* mdm = NULL;
double temperature;
double tau;
res_T res = RES_OK;
- ASSERT(rwalk && rng && T);
+ ASSERT(scn && rwalk && ctx && rng && T);
+
+ /* Get the current medium */
+ enc = scene_get_enclosure(scn, rwalk->enc_id);
+ res = scene_get_enclosure_medium(scn, enc, &mdm);
+ if(res != RES_OK) goto error;
/* Sample the time using the upper bound. */
tau = ssp_ran_exp(rng, mu);
@@ -52,11 +60,11 @@ time_rewind
if(rwalk->vtx.time > t0) goto exit;
/* Fetch the initial temperature */
- temperature = medium_get_temperature(rwalk->mdm, &rwalk->vtx);
+ temperature = medium_get_temperature(mdm, &rwalk->vtx);
if(SDIS_TEMPERATURE_IS_UNKNOWN(temperature)) {
- log_err(rwalk->mdm->dev, "the path reaches the limit condition but the "
+ log_err(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)),
+ medium_type_to_string(sdis_medium_get_type(mdm)),
SPLIT3(rwalk->vtx.P));
res = RES_BAD_ARG;
goto error;
@@ -75,8 +83,8 @@ time_rewind
}
if(ctx->green_path) {
- res = green_path_set_limit_vertex(ctx->green_path, rwalk->mdm,
- &rwalk->vtx, rwalk->elapsed_time);
+ res = green_path_set_limit_vertex
+ (ctx->green_path, mdm, &rwalk->vtx, rwalk->elapsed_time);
if(res != RES_OK) goto error;
}
diff --git a/src/sdis_misc.h b/src/sdis_misc.h
@@ -149,7 +149,8 @@ register_heat_vertex
extern LOCAL_SYM res_T
time_rewind
- (const double mu,
+ (struct sdis_scene* scn,
+ const double mu,
const double t0, /* Initial time */
struct ssp_rng* rng,
struct rwalk* rwalk,
diff --git a/src/sdis_realisation.c b/src/sdis_realisation.c
@@ -19,16 +19,37 @@
/*******************************************************************************
* Helper functions
******************************************************************************/
-static INLINE int
-check_ray_realisation_args(const struct ray_realisation_args* args)
+static INLINE res_T
+check_ray_realisation_args
+ (struct sdis_scene* scn,
+ const struct ray_realisation_args* args)
{
- return args
- && args->rng
- && args->medium
- && args->medium->type == SDIS_FLUID
- && args->time >= 0
- && args->picard_order > 0
- && (unsigned)args->diff_algo < SDIS_DIFFUSION_ALGORITHMS_COUNT__;
+ const struct enclosure* enc = NULL;
+ struct sdis_medium* mdm = NULL;
+ res_T res = RES_OK;
+ ASSERT(scn);
+
+ /* Check pointers */
+ if(!args || !args->rng) return RES_BAD_ARG;
+
+ if(args->time < 0) return RES_BAD_ARG;
+ if(args->picard_order <= 0) return RES_BAD_ARG;
+
+ if((unsigned)args->diff_algo >= SDIS_DIFFUSION_ALGORITHMS_COUNT__) {
+ return RES_BAD_ARG;
+ }
+
+ /* Check the enclosure id */
+ if(args->enc_id == ENCLOSURE_ID_NULL) {
+ return RES_BAD_ARG;
+ }
+ enc = scene_get_enclosure(scn, args->enc_id);
+ if((res = scene_get_enclosure_medium(scn, enc, &mdm)) != RES_OK) return res;
+ if(sdis_medium_get_type(mdm) != SDIS_FLUID) {
+ return RES_BAD_ARG;
+ }
+
+ return RES_OK;
}
/*******************************************************************************
@@ -51,13 +72,13 @@ ray_realisation_3d
struct temperature T = TEMPERATURE_NULL;
float dir[3];
res_T res = RES_OK;
- ASSERT(scn && weight && check_ray_realisation_args(args));
+ ASSERT(scn && weight && check_ray_realisation_args(scn, args) == RES_OK);
d3_set(rwalk.vtx.P, args->position);
rwalk.vtx.time = args->time;
rwalk.hit_3d = S3D_HIT_NULL;
rwalk.hit_side = SDIS_SIDE_NULL__;
- rwalk.mdm = args->medium;
+ rwalk.enc_id = args->enc_id;
ctx.heat_path = args->heat_path;
ctx.Tmin = scn->tmin;
@@ -69,7 +90,7 @@ ray_realisation_3d
ctx.max_branchings = args->picard_order - 1;
ctx.irealisation = args->irealisation;
ctx.diff_algo = args->diff_algo;
-
+
f3_set_d3(dir, args->direction);
/* Register the starting position against the heat path */
diff --git a/src/sdis_realisation.h b/src/sdis_realisation.h
@@ -60,7 +60,7 @@ sample_coupled_path_3d
******************************************************************************/
struct probe_realisation_args {
struct ssp_rng* rng;
- struct sdis_medium* medium; /* Medium into which the realisation starts */
+ unsigned enc_id; /* Enclosure into which the realisation starts */
double position[3]; /* Probe position */
double time; /* Observation time */
size_t picard_order; /* Picard order to estimate radiative temperature */
@@ -71,7 +71,7 @@ struct probe_realisation_args {
};
#define PROBE_REALISATION_ARGS_NULL__ { \
NULL, /* RNG */ \
- NULL, /* Medium */ \
+ ENCLOSURE_ID_NULL, /* Enclosure */ \
{0,0,0}, /* Position */ \
-1, /* Observation time */ \
0, /* Picard order */ \
@@ -182,7 +182,7 @@ boundary_flux_realisation_3d
******************************************************************************/
struct ray_realisation_args {
struct ssp_rng* rng;
- struct sdis_medium* medium; /* Medium into which the realisation starts */
+ unsigned enc_id; /* Enclosure into which the realisation starts */
double position[3]; /* Ray position */
double direction[3]; /* Ray direction */
double time; /* Observation time */
@@ -193,7 +193,7 @@ struct ray_realisation_args {
};
#define RAY_REALISATION_ARGS_NULL__ { \
NULL, /* RNG */ \
- NULL, /* Medium */ \
+ ENCLOSURE_ID_NULL, /* Enclosure */ \
{0,0,0}, /* Position */ \
{0,0,0}, /* Direction */ \
-1, /* Observation time */ \
diff --git a/src/sdis_realisation_Xd.h b/src/sdis_realisation_Xd.h
@@ -32,18 +32,19 @@
#ifndef SDIS_REALISATION_XD_H
#define SDIS_REALISATION_XD_H
-static INLINE int
+static INLINE res_T
check_probe_realisation_args(const struct probe_realisation_args* args)
{
return args
&& args->rng
- && args->medium
+ && args->enc_id != ENCLOSURE_ID_NULL
&& args->time >= 0
&& args->picard_order > 0
- && (unsigned)args->diff_algo < SDIS_DIFFUSION_ALGORITHMS_COUNT__;
+ && (unsigned)args->diff_algo < SDIS_DIFFUSION_ALGORITHMS_COUNT__
+ ? RES_OK : RES_BAD_ARG;
}
-static INLINE int
+static INLINE res_T
check_boundary_realisation_args(const struct boundary_realisation_args* args)
{
return args
@@ -55,10 +56,11 @@ check_boundary_realisation_args(const struct boundary_realisation_args* args)
&& args->time >= 0
&& args->picard_order > 0
&& (args->side == SDIS_FRONT || args->side == SDIS_BACK)
- && (unsigned)args->diff_algo < SDIS_DIFFUSION_ALGORITHMS_COUNT__;
+ && (unsigned)args->diff_algo < SDIS_DIFFUSION_ALGORITHMS_COUNT__
+ ? RES_OK : RES_BAD_ARG;
}
-static INLINE int
+static INLINE res_T
check_boundary_flux_realisation_args
(const struct boundary_flux_realisation_args* args)
{
@@ -71,7 +73,8 @@ check_boundary_flux_realisation_args
&& args->time >= 0
&& args->picard_order > 0
&& (args->solid_side == SDIS_FRONT || args->solid_side == SDIS_BACK)
- && (unsigned)args->diff_algo < SDIS_DIFFUSION_ALGORITHMS_COUNT__;
+ && (unsigned)args->diff_algo < SDIS_DIFFUSION_ALGORITHMS_COUNT__
+ ? RES_OK : RES_BAD_ARG;
}
#endif /* SDIS_REALISATION_XD_H */
@@ -171,27 +174,39 @@ XD(probe_realisation)
struct probe_realisation_args* args,
double* weight)
{
+ /* Starting enclosure/medium */
+ const struct enclosure* enc = NULL;
+ struct sdis_medium* mdm = NULL;
+
+ /* Random walk */
struct rwalk_context ctx = RWALK_CONTEXT_NULL;
struct rwalk rwalk = RWALK_NULL;
struct temperature T = TEMPERATURE_NULL;
+
+ /* Miscellaneous */
enum sdis_heat_vertex_type type;
double t0;
double (*get_initial_temperature)
(const struct sdis_medium* mdm,
const struct sdis_rwalk_vertex* vtx);
res_T res = RES_OK;
- ASSERT(scn && weight && check_probe_realisation_args(args));
+ ASSERT(scn && weight && check_probe_realisation_args(args) == RES_OK);
+
+ /* Get the enclosure medium */
+ enc = scene_get_enclosure(scn, args->enc_id);
+ res = scene_get_enclosure_medium(scn, enc, &mdm);
+ if(res != RES_OK) goto error;
- switch(args->medium->type) {
+ switch(sdis_medium_get_type(mdm)) {
case SDIS_FLUID:
T.func = XD(convective_path);
get_initial_temperature = fluid_get_temperature;
- t0 = fluid_get_t0(args->medium);
+ t0 = fluid_get_t0(mdm);
break;
case SDIS_SOLID:
T.func = XD(conductive_path);
get_initial_temperature = solid_get_temperature;
- t0 = solid_get_t0(args->medium);
+ t0 = solid_get_t0(mdm);
break;
default: FATAL("Unreachable code\n"); break;
}
@@ -200,7 +215,7 @@ XD(probe_realisation)
rwalk.vtx.time = args->time;
/* Register the starting position against the heat path */
- type = args->medium->type == SDIS_SOLID
+ type = sdis_medium_get_type(mdm) == SDIS_SOLID
? SDIS_HEAT_VERTEX_CONDUCTION
: SDIS_HEAT_VERTEX_CONVECTION;
res = register_heat_vertex(args->heat_path, &rwalk.vtx, 0, type, 0);
@@ -210,7 +225,7 @@ XD(probe_realisation)
double tmp;
/* Check the initial condition. */
rwalk.vtx.time = t0;
- tmp = get_initial_temperature(args->medium, &rwalk.vtx);
+ tmp = get_initial_temperature(mdm, &rwalk.vtx);
if(SDIS_TEMPERATURE_IS_KNOWN(tmp)) {
*weight = tmp;
goto exit;
@@ -225,7 +240,7 @@ XD(probe_realisation)
}
rwalk.XD(hit) = SXD_HIT_NULL;
- rwalk.mdm = args->medium;
+ rwalk.enc_id = args->enc_id;
ctx.green_path = args->green_path;
ctx.heat_path = args->heat_path;
@@ -267,13 +282,13 @@ XD(boundary_realisation)
float st[2];
#endif
res_T res = RES_OK;
- ASSERT(scn && weight && check_boundary_realisation_args(args));
+ ASSERT(scn && weight && check_boundary_realisation_args(args) == RES_OK);
T.func = XD(boundary_path);
rwalk.hit_side = args->side;
rwalk.XD(hit).distance = 0;
rwalk.vtx.time = args->time;
- rwalk.mdm = NULL; /* The random walk is at an interface between 2 media */
+ rwalk.enc_id = ENCLOSURE_ID_NULL; /* At an interface between 2 enclosures */
#if SDIS_XD_DIMENSION == 2
st = (float)args->uv[0];
@@ -332,14 +347,14 @@ XD(boundary_flux_realisation)
struct boundary_flux_realisation_args* args,
struct bound_flux_result* result)
{
+ /* Random walk */
struct rwalk_context ctx = RWALK_CONTEXT_NULL;
- struct rwalk rwalk;
- struct temperature T;
+ struct rwalk rwalk = RWALK_NULL;
+ struct temperature T = TEMPERATURE_NULL;
+
+ /* Boundary */
struct sXd(attrib) attr;
struct sXd(primitive) prim;
- struct sdis_interface* interf = NULL;
- struct sdis_medium* fluid_mdm = NULL;
-
#if SDIS_XD_DIMENSION == 2
float st;
#else
@@ -347,13 +362,17 @@ XD(boundary_flux_realisation)
#endif
double P[SDIS_XD_DIMENSION];
float N[SDIS_XD_DIMENSION];
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+ enum sdis_side fluid_side;
+
+ /* Miscellaneous */
double Tmin, Tmin2, Tmin3;
double That, That2, That3;
- enum sdis_side fluid_side;
res_T res = RES_OK;
char compute_radiative;
char compute_convective;
- ASSERT(scn && result && check_boundary_flux_realisation_args(args));
+
+ ASSERT(scn && result && check_boundary_flux_realisation_args(args) == RES_OK);
#if SDIS_XD_DIMENSION == 2
#define SET_PARAM(Dest, Src) (Dest).u = (Src);
@@ -386,12 +405,12 @@ XD(boundary_flux_realisation)
SXD(primitive_get_attrib(&prim, SXD_GEOMETRY_NORMAL, st, &attr));
fX(set)(N, attr.value);
- #define RESET_WALK(Side, Mdm) { \
+ #define RESET_WALK(Side, EncId) { \
rwalk = RWALK_NULL; \
rwalk.hit_side = (Side); \
rwalk.XD(hit).distance = 0; \
rwalk.vtx.time = args->time; \
- rwalk.mdm = (Mdm); \
+ rwalk.enc_id = (EncId); \
rwalk.XD(hit).prim = prim; \
SET_PARAM(rwalk.XD(hit), st); \
ctx.Tmin = Tmin; \
@@ -408,19 +427,18 @@ XD(boundary_flux_realisation)
} (void)0
/* Compute boundary temperature */
- RESET_WALK(args->solid_side, NULL);
+ RESET_WALK(args->solid_side, ENCLOSURE_ID_NULL);
T.func = XD(boundary_path);
res = XD(sample_coupled_path)(scn, &ctx, &rwalk, args->rng, &T);
if(res != RES_OK) return res;
result->Tboundary = T.value;
- /* Fetch the fluid medium */
- interf = scene_get_interface(scn, (unsigned)args->iprim);
- fluid_mdm = interface_get_medium(interf, fluid_side);
+ /* Get the enclosures */
+ scene_get_enclosure_ids(scn, (unsigned)args->iprim, enc_ids);
/* Compute radiative temperature */
if(compute_radiative) {
- RESET_WALK(fluid_side, fluid_mdm);
+ RESET_WALK(fluid_side, enc_ids[fluid_side]);
T.func = XD(radiative_path);
res = XD(sample_coupled_path)(scn, &ctx, &rwalk, args->rng, &T);
if(res != RES_OK) return res;
@@ -430,7 +448,7 @@ XD(boundary_flux_realisation)
/* Compute fluid temperature */
if(compute_convective) {
- RESET_WALK(fluid_side, fluid_mdm);
+ RESET_WALK(fluid_side, enc_ids[fluid_side]);
T.func = XD(convective_path);
res = XD(sample_coupled_path)(scn, &ctx, &rwalk, args->rng, &T);
if(res != RES_OK) return res;
diff --git a/src/sdis_scene.c b/src/sdis_scene.c
@@ -460,7 +460,7 @@ scene_get_enclosure_medium
ASSERT(scn && enc && out_mdm);
/* Check that the enclosure doesn't surround multiple media */
- if(enc->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
+ if(enc->medium_id == MEDIUM_ID_MULTI) {
log_warn(scn->dev,
"%s: invalid medium request. The enclosure includes several media.\n",
FUNC_NAME);
diff --git a/src/sdis_scene_Xd.h b/src/sdis_scene_Xd.h
@@ -823,7 +823,7 @@ XD(register_enclosure)(struct sdis_scene* scn, struct sencXd(enclosure)* enc)
/* Setup the medium id of the enclosure */
if(header.enclosed_media_count > 1) {
- enc_data->medium_id = ENCLOSURE_ID_MULTI_MEDIA;
+ enc_data->medium_id = MEDIUM_ID_MULTI;
} else {
SENCXD(enclosure_get_medium(enc, 0, &enc_data->medium_id));
}
diff --git a/src/sdis_scene_c.h b/src/sdis_scene_c.h
@@ -26,8 +26,8 @@
#include <limits.h>
-#define ENCLOSURE_ID_NULL UINT_MAX-1
-#define ENCLOSURE_ID_MULTI_MEDIA (UINT_MAX-1)
+#define MEDIUM_ID_MULTI UINT_MAX
+#define ENCLOSURE_ID_NULL UINT_MAX
struct prim_prop {
struct sdis_interface* interf;
@@ -98,7 +98,7 @@ enclosure_init(struct mem_allocator* allocator, struct enclosure* enc)
enc->S_over_V = 0;
enc->V = 0;
enc->hc_upper_bound = 0;
- enc->medium_id = ENCLOSURE_ID_MULTI_MEDIA;
+ enc->medium_id = MEDIUM_ID_MULTI;
}
static INLINE void
diff --git a/src/sdis_solve_camera.c b/src/sdis_solve_camera.c
@@ -84,7 +84,7 @@ static res_T
solve_pixel
(struct sdis_scene* scn,
struct ssp_rng* rng,
- struct sdis_medium* mdm,
+ const unsigned enc_id,
const struct sdis_camera* cam,
const double time_range[2], /* Observation time */
const size_t ipix[2], /* Pixel coordinate in the image space */
@@ -99,7 +99,7 @@ solve_pixel
struct sdis_heat_path* pheat_path = NULL;
size_t irealisation;
res_T res = RES_OK;
- ASSERT(scn && mdm && rng && cam && ipix && nrealisations);
+ ASSERT(scn && rng && cam && ipix && nrealisations);
ASSERT(pix_sz && pix_sz[0] > 0 && pix_sz[1] > 0);
ASSERT(pixel && time_range);
@@ -136,7 +136,7 @@ solve_pixel
/* Launch the realisation */
realis_args.rng = rng;
- realis_args.medium = mdm;
+ realis_args.enc_id = enc_id;
realis_args.time = time;
realis_args.picard_order = picard_order;
realis_args.heat_path = pheat_path;
@@ -195,7 +195,7 @@ static res_T
solve_tile
(struct sdis_scene* scn,
struct ssp_rng* rng,
- struct sdis_medium* mdm,
+ const unsigned enc_id,
const struct sdis_camera* cam,
const double time_range[2],
const size_t tile_org[2], /* Origin of the tile in pixel space */
@@ -211,7 +211,7 @@ solve_tile
size_t mcode; /* Morton code of the tile pixel */
size_t npixels;
res_T res = RES_OK;
- ASSERT(scn && rng && mdm && cam && spp);
+ ASSERT(scn && rng && cam && spp);
ASSERT(tile_size && tile_size[0] && tile_size[1]);
ASSERT(pix_sz && pix_sz[0] > 0 && pix_sz[1] > 0 && time_range);
@@ -241,8 +241,8 @@ solve_tile
estimator = estimator_buffer_grab(buf, ipix_image[0], ipix_image[1]);
}
res = solve_pixel
- (scn, rng, mdm, cam, time_range, ipix_image, spp, register_paths, pix_sz,
- picard_order, diff_algo, estimator, pixel);
+ (scn, rng, enc_id, cam, time_range, ipix_image, spp, register_paths,
+ pix_sz, picard_order, diff_algo, estimator, pixel);
if(res != RES_OK) goto error;
}
@@ -504,13 +504,13 @@ sdis_solve_camera
/* Stardis variables */
struct sdis_estimator_buffer* buf = NULL;
- struct sdis_medium* mdm = NULL;
/* Random number generators */
struct ssp_rng_proxy* rng_proxy = NULL;
struct ssp_rng** per_thread_rng = NULL;
- /* Enclosure in which the probe lies */
+ /* Enclosure & medium in which the probe lies */
+ struct sdis_medium* mdm = NULL;
unsigned enc_id = ENCLOSURE_ID_NULL;
/* Miscellaneous */
@@ -541,12 +541,13 @@ sdis_solve_camera
/* Retrieve the medium in which the submitted position lies */
res = scene_get_enclosure_id(scn, args->cam->position, &enc_id);
+ if(res != RES_OK) goto error;
res = scene_get_enclosure_medium(scn, scene_get_enclosure(scn, enc_id), &mdm);
if(res != RES_OK) goto error;
- if(mdm->type != SDIS_FLUID) {
+ if(sdis_medium_get_type(mdm) != SDIS_FLUID) {
log_err(scn->dev,
- "%s: the camera position `%g %g %g' must be in a fluid medium.\n",
+ "%s: the camera position (%g, %g, %g) must be in a fluid medium.\n",
FUNC_NAME, SPLIT3(args->cam->position));
res = RES_BAD_ARG;
goto error;
@@ -662,7 +663,7 @@ sdis_solve_camera
/* Draw the tile */
res_local = solve_tile
- (scn, rng, mdm, args->cam, args->time_range, tile_org, tile_sz,
+ (scn, rng, enc_id, args->cam, args->time_range, tile_org, tile_sz,
args->spp, register_paths, pix_sz, args->picard_order, args->diff_algo,
buf, tile);
if(res_local != RES_OK) {
diff --git a/src/sdis_solve_medium_Xd.h b/src/sdis_solve_medium_Xd.h
@@ -39,7 +39,7 @@
*/
struct enclosure_cumul {
- const struct enclosure* enc;
+ unsigned enc_id;
double cumul;
};
@@ -79,12 +79,13 @@ compute_medium_enclosure_cumulative
while(!htable_enclosure_iterator_eq(&it, &end)) {
struct enclosure_cumul enc_cumul;
const struct enclosure* enc = htable_enclosure_iterator_data_get(&it);
+ const unsigned* enc_id = htable_enclosure_iterator_key_get(&it);
htable_enclosure_iterator_next(&it);
if(sdis_medium_get_id(mdm) != enc->medium_id) continue;
accum += enc->V;
- enc_cumul.enc = enc;
+ enc_cumul.enc_id = *enc_id;
enc_cumul.cumul = accum;
res = darray_enclosure_cumul_push_back(cumul, &enc_cumul);
if(res != RES_OK) goto error;
@@ -105,7 +106,7 @@ error:
goto exit;
}
-static const struct enclosure*
+static unsigned
sample_medium_enclosure
(const struct darray_enclosure_cumul* cumul, struct ssp_rng* rng)
{
@@ -137,7 +138,7 @@ sample_medium_enclosure
enc_cumul_found = enc_cumuls + i;
}
- return enc_cumul_found->enc;
+ return enc_cumul_found->enc_id;
}
static INLINE res_T
@@ -217,17 +218,22 @@ check_compute_power_args(const struct sdis_compute_power_args* args)
******************************************************************************/
static res_T
XD(sample_enclosure_position)
- (const struct enclosure* enc,
+ (struct sdis_scene* scn,
+ const unsigned enc_id,
struct ssp_rng* rng,
double pos[DIM])
{
const size_t MAX_NCHALLENGES = 1000;
+
+ const struct enclosure* enc = NULL;
float lower[DIM], upper[DIM];
size_t ichallenge;
size_t i;
res_T res = RES_OK;
- ASSERT(enc && rng && pos);
+ ASSERT(scn && rng && pos);
+ ASSERT(enc_id != ENCLOSURE_ID_NULL);
+ enc = scene_get_enclosure(scn, enc_id);
SXD(scene_view_get_aabb(enc->sXd(view), lower, upper));
FOR_EACH(i, 0, DIM) {
@@ -396,13 +402,13 @@ XD(solve_medium)
struct accum* acc_time = &per_thread_acc_time[ithread];
struct green_path_handle* pgreen_path = NULL;
struct green_path_handle green_path = GREEN_PATH_HANDLE_NULL;
- const struct enclosure* enc = NULL;
struct sdis_heat_path* pheat_path = NULL;
struct sdis_heat_path heat_path;
double weight;
double time;
double pos[DIM];
size_t n;
+ unsigned enc_id = ENCLOSURE_ID_NULL;
int pcent;
res_T res_local = RES_OK;
res_T res_simul = RES_OK;
@@ -425,8 +431,8 @@ XD(solve_medium)
/* Uniformly Sample an enclosure that surround the submitted medium and
* uniformly sample a position into it */
- enc = sample_medium_enclosure(&cumul, rng);
- res_local = XD(sample_enclosure_position)(enc, rng, pos);
+ enc_id = sample_medium_enclosure(&cumul, rng);
+ res_local = XD(sample_enclosure_position)(scn, enc_id, rng, pos);
if(res_local != RES_OK) {
log_err(scn->dev, "%s: could not sample a medium position.\n", FUNC_NAME);
ATOMIC_SET(&res, res_local);
@@ -435,7 +441,7 @@ XD(solve_medium)
/* Run a probe realisation */
realis_args.rng = rng;
- realis_args.medium = args->medium;
+ realis_args.enc_id = enc_id;
realis_args.time = time;
realis_args.picard_order = args->picard_order;
realis_args.green_path = pgreen_path;
@@ -685,10 +691,10 @@ XD(compute_power)
struct ssp_rng* rng = per_thread_rng[ithread];
struct accum* acc_mpow = &per_thread_acc_mpow[ithread];
struct accum* acc_time = &per_thread_acc_time[ithread];
- const struct enclosure* enc = NULL;
double power = 0;
double usec = 0;
size_t n = 0;
+ unsigned enc_id = ENCLOSURE_ID_NULL;
int pcent = 0;
res_T res_local = RES_OK;
@@ -702,8 +708,8 @@ XD(compute_power)
/* Uniformly Sample an enclosure that surround the submitted medium and
* uniformly sample a position into it */
- enc = sample_medium_enclosure(&cumul, rng);
- res_local = XD(sample_enclosure_position)(enc, rng, vtx.P);
+ enc_id = sample_medium_enclosure(&cumul, rng);
+ res_local = XD(sample_enclosure_position)(scn, enc_id, rng, vtx.P);
if(res_local != RES_OK) {
log_err(scn->dev, "%s: could not sample a medium position.\n", FUNC_NAME);
ATOMIC_SET(&res, res_local);
diff --git a/src/sdis_solve_probe_Xd.h b/src/sdis_solve_probe_Xd.h
@@ -118,7 +118,6 @@ XD(solve_one_probe)
{
/* Enclosure in which the probe lies */
unsigned enc_id = ENCLOSURE_ID_NULL;
- struct sdis_medium* mdm = NULL;
size_t irealisation = 0;
res_T res = RES_OK;
@@ -131,8 +130,6 @@ XD(solve_one_probe)
/* Retrieve the medium in which the submitted position lies */
res = scene_get_enclosure_id(scn, args->position, &enc_id);
if(res != RES_OK) goto error;
- res = scene_get_enclosure_medium(scn, scene_get_enclosure(scn, enc_id), &mdm);
- if(res != RES_OK) goto error;
FOR_EACH(irealisation, 0, args->nrealisations) {
struct probe_realisation_args realis_args = PROBE_REALISATION_ARGS_NULL;
@@ -149,7 +146,7 @@ XD(solve_one_probe)
/* Run a realisation */
realis_args.rng = rng;
- realis_args.medium = mdm;
+ realis_args.enc_id = enc_id;
realis_args.time = time;
realis_args.picard_order = args->picard_order;
realis_args.irealisation = irealisation;
@@ -196,7 +193,6 @@ XD(solve_probe)
size_t nthreads = 0;
/* Stardis variables */
- struct sdis_medium* mdm = NULL;
struct sdis_estimator* estimator = NULL;
struct sdis_green_function* green = NULL;
struct sdis_green_function** per_thread_green = NULL;
@@ -206,6 +202,7 @@ XD(solve_probe)
struct ssp_rng** per_thread_rng = NULL;
/* Enclosure in which the probe lies */
+ const struct enclosure* enc = NULL;
unsigned enc_id = ENCLOSURE_ID_NULL;
/* Miscellaneous */
@@ -263,11 +260,18 @@ XD(solve_probe)
if(!per_thread_acc_temp) { res = RES_MEM_ERR; goto error; }
if(!per_thread_acc_time) { res = RES_MEM_ERR; goto error; }
- /* Retrieve the medium in which the submitted position lies */
+ /* Retrieve the enclosure in which the submitted position lies */
res = scene_get_enclosure_id(scn, args->position, &enc_id);
if(res != RES_OK) goto error;
- res = scene_get_enclosure_medium(scn, scene_get_enclosure(scn, enc_id), &mdm);
- if(res != RES_OK) goto error;
+
+ /* Check that the enclosure does not contain multiple materials */
+ enc = scene_get_enclosure(scn, enc_id);
+ if(enc->medium_id == MEDIUM_ID_MULTI) {
+ log_err(scn->dev, "%s: probe is in an enclosure with several media "
+ "-- pos=("FORMAT_VECX")\n", FUNC_NAME, SPLITX(args->position));
+ res = RES_BAD_OP;
+ goto error;
+ }
/* Create the per thread green function */
if(out_green) {
@@ -337,7 +341,7 @@ XD(solve_probe)
/* Invoke the probe realisation */
realis_args.rng = rng;
- realis_args.medium = mdm;
+ realis_args.enc_id = enc_id;
realis_args.time = time;
realis_args.picard_order = args->picard_order;
realis_args.green_path = pgreen_path;
