commit b80e9f194fd95ae330570ed440cc621120c78337
parent 34a3152c88013b25186945e52a8678cb41bc5b7c
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Wed, 29 May 2024 19:01:09 +0200
Start tracking enclosures instead of media
The sampled path (conductive, convective or radiative) must remain in
the same enclosure, except when it reaches a boundary where it is
coupled to other physics. This assumption is purely geometric and
can be verified by controlling the enclosures during the
random walk. Media are physical properties that can be shared between
enclosures. Their verification is therefore less robust than that of the
enclosures.
At the current stage, the random walk still keeps track of the current
medium, but this commit replaces the function that retrieves a medium
from a position with its counterpart that returns an enclosure. We can
then retrieve the enclosure medium as required by the solver.
Enclosure tracking could also be the solution to enable radiative path
sampling in an enclosure that surrounds several media, such as those
used to describe fluids with a known temperature for different dirichlet
boundary conditions. While convective paths cannot be sampled in this
type of enclosure, radiative paths can, although an error is currently
being reported. This currently prevents the rendering of scenes whose
"visible" boundaries have different boundary conditions.
Diffstat:
9 files changed, 195 insertions(+), 104 deletions(-)
diff --git a/src/sdis_heat_path_boundary_Xd_c.h b/src/sdis_heat_path_boundary_Xd_c.h
@@ -334,6 +334,7 @@ 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;
@@ -358,27 +359,41 @@ XD(find_reinjection_ray)
(scn->sXd(view), org, args->dir1, range, &filter_data, &hit1));
/* Retrieve the medium at the reinjection pos along dir0 */
- if(SXD_HIT_NONE(&hit0)) {
- XD(move_pos)(dX(set)(tmp, ray->org), args->dir0, (float)args->distance);
- res = scene_get_medium_in_closed_boundaries(scn, tmp, &mdm0);
- if(res == RES_BAD_OP) { mdm0 = NULL; res = RES_OK; }
- if(res != RES_OK) goto error;
- } else {
+ if(!SXD_HIT_NONE(&hit0)) {
interf = scene_get_interface(scn, hit0.prim.prim_id);
side = fX(dot)(args->dir0, hit0.normal) < 0 ? SDIS_FRONT : SDIS_BACK;
mdm0 = interface_get_medium(interf, 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; }
+ 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 */
- if(SXD_HIT_NONE(&hit1)) {
- XD(move_pos)(dX(set)(tmp, ray->org), args->dir1, (float)args->distance);
- res = scene_get_medium_in_closed_boundaries(scn, tmp, &mdm1);
- if(res == RES_BAD_OP) { mdm1 = NULL; res = RES_OK; }
- if(res != RES_OK) goto error;
- } else {
+ if(!SXD_HIT_NONE(&hit1)) {
interf = scene_get_interface(scn, hit1.prim.prim_id);
side = fX(dot)(args->dir1, hit1.normal) < 0 ? SDIS_FRONT : SDIS_BACK;
mdm1 = interface_get_medium(interf, 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; }
+ 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;
@@ -514,11 +529,18 @@ XD(find_reinjection_ray_and_check_validity)
if(res != RES_OK) goto error;
if(SXD_HIT_NONE(&ray->hit)) {
- /* Check medium consistency at the reinjection position */
+ 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_medium_in_closed_boundaries(scn, pos, &reinject_mdm);
+ 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) {
res = RES_BAD_OP;
goto error;
@@ -630,7 +652,7 @@ XD(sample_reinjection_step_solid_fluid)
* error */
scene_get_enclosure_ids(scn, args->rwalk->hit.prim.prim_id, enc_ids);
solid_enclosure = scene_get_enclosure(scn, enc_ids[args->side]);
- if(solid_enclosure->medium_id == ENCLOSURE_MULTI_MEDIA) {
+ if(solid_enclosure->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
step->hit = SXD_HIT_NULL;
fX(normalize)(step->direction, args->rwalk->hit.normal);
if(args->side == SDIS_BACK) fX(minus)(step->direction, step->direction);
@@ -763,13 +785,13 @@ XD(sample_reinjection_step_solid_solid)
scene_get_enclosure_ids(scn, args_frt->rwalk->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_MULTI_MEDIA) {
+ if(enclosure_frt->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
step_frt->hit = SXD_HIT_NULL;
fX(normalize)(step_frt->direction, args_frt->rwalk->hit.normal);
step_frt->distance = (float)args_frt->distance;
multiple_media_frt = 1;
}
- if(enclosure_bck->medium_id == ENCLOSURE_MULTI_MEDIA) {
+ if(enclosure_bck->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
step_bck->hit = SXD_HIT_NULL;
fX(normalize)(step_bck->direction, args_bck->rwalk->hit.normal);
fX(minus)(step_bck->direction, step_bck->direction);
diff --git a/src/sdis_heat_path_conductive_delta_sphere_Xd.h b/src/sdis_heat_path_conductive_delta_sphere_Xd.h
@@ -111,7 +111,7 @@ XD(sample_next_step)
static res_T
XD(sample_next_step_robust)
(struct sdis_scene* scn,
- struct sdis_medium* current_mdm,
+ const unsigned current_enc_id,
struct ssp_rng* rng,
const double pos[DIM],
const float delta_solid,
@@ -121,13 +121,15 @@ XD(sample_next_step_robust)
struct sXd(hit)* hit1, /* Hit used to adjust delta */
float* out_delta)
{
- struct sdis_medium* mdm;
+ unsigned enc_id = ENCLOSURE_ID_NULL;
float delta;
float org[DIM];
const size_t MAX_ATTEMPTS = 100;
size_t iattempt = 0;
res_T res = RES_OK;
- ASSERT(scn && current_mdm && rng && pos && delta_solid > 0);
+ 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);
@@ -141,19 +143,17 @@ XD(sample_next_step_robust)
/* Retrieve the medium of the next step */
if(hit0->distance > delta) {
XD(move_pos)(dX(set)(pos_next, pos), dir0, delta);
- res = scene_get_medium_in_closed_boundaries(scn, pos_next, &mdm);
- if(res == RES_BAD_OP) { mdm = NULL; res = RES_OK; }
+ res = scene_get_enclosure_id_in_closed_boundaries(scn, pos_next, &enc_id);
+ if(res == RES_BAD_OP) { enc_id = ENCLOSURE_ID_NULL; res = RES_OK; }
if(res != RES_OK) goto error;
} else {
- struct sdis_interface* interf;
- enum sdis_side side;
- interf = scene_get_interface(scn, hit0->prim.prim_id);
- side = fX(dot)(dir0, hit0->normal) < 0 ? SDIS_FRONT : SDIS_BACK;
- mdm = interface_get_medium(interf, side);
+ unsigned enc_ids[2] = {ENCLOSURE_ID_NULL, ENCLOSURE_ID_NULL};
+ scene_get_enclosure_ids(scn, hit0->prim.prim_id, enc_ids);
+ enc_id = fX(dot)(dir0, hit0->normal) < 0 ? enc_ids[0] : enc_ids[1];
}
/* Check medium consistency */
- if(current_mdm != mdm) {
+ if(current_enc_id != enc_id) {
#if 0
#if DIM == 2
log_warn(scn->dev,
@@ -166,7 +166,7 @@ XD(sample_next_step_robust)
#endif
#endif
}
- } while(current_mdm != mdm && ++iattempt < MAX_ATTEMPTS);
+ } while(current_enc_id != enc_id && ++iattempt < MAX_ATTEMPTS);
/* Handle error */
if(iattempt >= MAX_ATTEMPTS) {
@@ -336,24 +336,38 @@ XD(conductive_path_delta_sphere)
struct ssp_rng* rng,
struct XD(temperature)* T)
{
- double position_start[DIM];
+ /* Enclosure/medium in which the conductive path starts */
+ struct sdis_medium* mdm = NULL;
+ unsigned enc_id = ENCLOSURE_ID_NULL;
+
+ /* Physical properties */
struct solid_props props_ref = SOLID_PROPS_NULL;
double green_power_term = 0;
- struct sdis_medium* mdm;
+
+ /* Miscellaneous */
+ double position_start[DIM] = {0};
size_t istep = 0; /* Help for debug */
res_T res = RES_OK;
+
+ /* Check pre-conditions */
ASSERT(scn && rwalk && rng && T);
ASSERT(rwalk->mdm->type == SDIS_SOLID);
- (void)ctx, (void)istep;
+
+ (void)ctx, (void)istep; /* Avoid "unsued variable" warnings */
+
+ 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 the random walk consistency */
- res = scene_get_medium_in_closed_boundaries(scn, rwalk->vtx.P, &mdm);
- if(res != RES_OK || mdm != rwalk->mdm) {
+ if(mdm != rwalk->mdm) {
log_err(scn->dev, "%s: invalid solid random walk. "
"Unexpected medium at {%g, %g, %g}.\n", FUNC_NAME, SPLIT3(rwalk->vtx.P));
res = RES_BAD_OP_IRRECOVERABLE;
goto error;
}
+
/* Save the submitted position */
dX(set)(position_start, rwalk->vtx.P);
@@ -410,7 +424,7 @@ XD(conductive_path_delta_sphere)
fX_set_dX(org, rwalk->vtx.P);
/* Sample the direction to walk toward and compute the distance to travel */
- res = XD(sample_next_step_robust)(scn, mdm, rng, rwalk->vtx.P,
+ res = XD(sample_next_step_robust)(scn, enc_id, rng, rwalk->vtx.P,
(float)props.delta, dir0, dir1, &hit0, &hit1, &delta);
if(res != RES_OK) goto error;
diff --git a/src/sdis_heat_path_conductive_wos_Xd.h b/src/sdis_heat_path_conductive_wos_Xd.h
@@ -30,11 +30,14 @@ XD(check_medium_consistency)
(struct sdis_scene* scn,
const struct XD(rwalk)* rwalk)
{
+ unsigned enc_id = ENCLOSURE_ID_NULL;
struct sdis_medium* mdm = NULL;
res_T res = RES_OK;
ASSERT(rwalk);
- res = scene_get_medium_in_closed_boundaries(scn, rwalk->vtx.P, &mdm);
+ 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 */
diff --git a/src/sdis_heat_path_convective_Xd.h b/src/sdis_heat_path_convective_Xd.h
@@ -211,7 +211,7 @@ XD(fetch_fluid_enclosure)
/* Fetch the enclosure data */
enc = scene_get_enclosure(scn, enc_id);
ASSERT(enc != NULL);
- if(enc->medium_id == ENCLOSURE_MULTI_MEDIA) {
+ if(enc->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
/* The enclosures with multiple media are used to describe limit
* conditions and therefore they cannot be fetched */
log_err(scn->dev,
diff --git a/src/sdis_scene.c b/src/sdis_scene.c
@@ -427,25 +427,57 @@ scene_get_interface(const struct sdis_scene* scn, const unsigned iprim)
}
res_T
-scene_get_medium
+scene_get_enclosure_id
(struct sdis_scene* scn,
const double pos[],
- struct sdis_medium** out_medium)
+ unsigned* enc_id)
{
return scene_is_2d(scn)
- ? scene_get_medium_2d(scn, pos, out_medium)
- : scene_get_medium_3d(scn, pos, out_medium);
+ ? scene_get_enclosure_id_2d(scn, pos, enc_id)
+ : scene_get_enclosure_id_3d(scn, pos, enc_id);
}
res_T
-scene_get_medium_in_closed_boundaries
+scene_get_enclosure_id_in_closed_boundaries
(struct sdis_scene* scn,
const double pos[],
- struct sdis_medium** out_medium)
+ unsigned* enc_id)
{
return scene_is_2d(scn)
- ? scene_get_medium_in_closed_boundaries_2d(scn, pos, out_medium)
- : scene_get_medium_in_closed_boundaries_3d(scn, pos, out_medium);
+ ? scene_get_enclosure_id_in_closed_boundaries_2d(scn, pos, enc_id)
+ : scene_get_enclosure_id_in_closed_boundaries_3d(scn, pos, enc_id);
+}
+
+res_T
+scene_get_enclosure_medium
+ (struct sdis_scene* scn,
+ const struct enclosure* enc,
+ struct sdis_medium** out_mdm)
+{
+ struct sdis_medium* mdm = NULL;
+ res_T res = RES_OK;
+
+ ASSERT(scn && enc && out_mdm);
+
+ /* Check that the enclosure doesn't surround multiple media */
+ if(enc->medium_id == ENCLOSURE_ID_MULTI_MEDIA) {
+ log_warn(scn->dev,
+ "%s: invalid medium request. The enclosure includes several media.\n",
+ FUNC_NAME);
+ res = RES_BAD_OP;
+ goto error;
+ }
+
+ /* Obtain enclosure medium */
+ ASSERT(enc->medium_id < darray_medium_size_get(&scn->media));
+ mdm = darray_medium_data_get(&scn->media)[enc->medium_id];
+
+error:
+ *out_mdm = mdm;
+ goto exit;
+exit:
+ mdm = NULL;
+ return res;
}
res_T
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_MULTI_MEDIA;
+ enc_data->medium_id = ENCLOSURE_ID_MULTI_MEDIA;
} else {
SENCXD(enclosure_get_medium(enc, 0, &enc_data->medium_id));
}
@@ -1094,13 +1094,12 @@ error:
/*******************************************************************************
* Local functions
******************************************************************************/
-static INLINE res_T
-XD(scene_get_medium)
+static res_T
+XD(scene_get_enclosure_id)
(struct sdis_scene* scn,
const double pos[DIM],
- struct sdis_medium** out_medium)
+ unsigned* out_enc_id)
{
- struct sdis_medium* medium = NULL;
size_t iprim, nprims;
float P[DIM];
/* Range of the parametric coordinate into which positions are challenged */
@@ -1110,6 +1109,7 @@ XD(scene_get_medium)
float st[3][2];
#endif
size_t nsteps = 3;
+ unsigned enc_id = ENCLOSURE_ID_NULL;
res_T res = RES_OK;
ASSERT(scn && pos);
@@ -1162,8 +1162,8 @@ XD(scene_get_medium)
} while((SXD_HIT_NONE(&hit) || HIT_ON_BOUNDARY(&hit, P, dir))
&& ++istep < nsteps);
- /* No valid intersection is found on the current primitive. Challenge
- * another. */
+ /* No valid intersection is found on the current primitive.
+ * Challenge another. */
if(istep >= nsteps) continue;
fX(normalize)(N, hit.normal);
@@ -1171,18 +1171,17 @@ XD(scene_get_medium)
/* Not too close and not roughly orthognonal */
if(hit.distance > 1.e-6 && absf(cos_N_dir) > 1.e-2f) {
- const struct enclosure* enclosure = NULL;
unsigned enc_ids[2];
- const struct sdis_interface* interf;
-
- interf = scene_get_interface(scn, hit.prim.prim_id);
scene_get_enclosure_ids(scn, hit.prim.prim_id, enc_ids);
- break;
+ enc_id = cos_N_dir < 0 ? enc_ids[0] : enc_ids[1];
+
+ break; /* That's all folks */
}
}
if(iprim >= nprims) {
- log_warn(scn->dev, "%s: could not retrieve the medium at {%g, %g, %g}.\n",
+ log_warn(scn->dev,
+ "%s: cannot retrieve current enclosure at {%g, %g, %g}.\n",
FUNC_NAME, P[0], P[1], DIM == 3 ? P[2] : 0);
res = RES_BAD_OP;
goto error;
@@ -1190,31 +1189,32 @@ XD(scene_get_medium)
if(iprim > 10 && iprim > (size_t)((double)nprims * 0.05)) {
log_warn(scn->dev,
- "%s: performance issue. Up to %lu primitives were tested to define the "
- "current medium at {%g, %g, %g}.\n",
+ "%s: performance issue. Up to %lu primitives were tested to find "
+ "current enclosure at {%g, %g, %g}.\n",
FUNC_NAME, (unsigned long)iprim, P[0], P[1], DIM == 3 ? P[2] : 0);
}
exit:
- *out_medium = medium;
+ *out_enc_id = enc_id;
return res;
error:
+ enc_id = ENCLOSURE_ID_NULL;
goto exit;
}
-static INLINE res_T
-XD(scene_get_medium_in_closed_boundaries)
+static res_T
+XD(scene_get_enclosure_id_in_closed_boundaries)
(struct sdis_scene* scn,
const double pos[DIM],
- struct sdis_medium** out_medium)
+ unsigned* out_enc_id)
{
- struct sdis_medium* medium = NULL;
- float P[DIM];
- float frame[DIM*DIM];
float dirs[6][3] = {{1,0,0},{-1,0,0},{0,1,0},{0,-1,0},{0,0,1},{0,0,-1}};
+ float frame[DIM*DIM];
+ float P[DIM];
+ unsigned enc_id = ENCLOSURE_ID_NULL;
int idir;
res_T res = RES_OK;
- ASSERT(scn && pos);
+ ASSERT(scn && pos && out_enc_id);
/* Build a frame that will be used to rotate the main axis by PI/4 around
* each axis. This can avoid numerical issues when geometry is discretized
@@ -1222,8 +1222,6 @@ XD(scene_get_medium_in_closed_boundaries)
#if DIM == 2
f22_rotation(frame, (float)PI/4);
#else
-/* N[0] = N[1] = N[2] = (float)(1.0 / sqrt(3.0));*/
-/* f33_basis(frame, N);*/
f33_rotation(frame, (float)PI/4, (float)PI/4, (float)PI/4);
#endif
@@ -1252,22 +1250,23 @@ XD(scene_get_medium_in_closed_boundaries)
/* Not too close and not roughly orthogonal */
if(hit.distance > 1.e-6 && absf(cos_N_dir) > 1.e-2f) {
- const struct sdis_interface* interf;
- interf = scene_get_interface(scn, hit.prim.prim_id);
- medium = interface_get_medium
- (interf, cos_N_dir < 0 ? SDIS_FRONT : SDIS_BACK);
- break;
+ unsigned enc_ids[2];
+ scene_get_enclosure_ids(scn, hit.prim.prim_id, enc_ids);
+ enc_id = cos_N_dir < 0 ? enc_ids[0] : enc_ids[1];
+
+ break; /* That's all folks */
}
}
- if(idir >= 2*DIM) {
- res = XD(scene_get_medium)(scn, pos, &medium);
+ if(idir >= 2*DIM) { /* Fallback to scene_get_enclosure_id function */
+ res = XD(scene_get_enclosure_id)(scn, pos, &enc_id);
if(res != RES_OK) goto error;
}
exit:
- *out_medium = medium;
+ *out_enc_id = enc_id;
return res;
error:
+ enc_id = ENCLOSURE_ID_NULL;
goto exit;
}
diff --git a/src/sdis_scene_c.h b/src/sdis_scene_c.h
@@ -26,6 +26,9 @@
#include <limits.h>
+#define ENCLOSURE_ID_NULL UINT_MAX-1
+#define ENCLOSURE_ID_MULTI_MEDIA (UINT_MAX-1)
+
struct prim_prop {
struct sdis_interface* interf;
unsigned front_enclosure; /* Id of the front facing enclosure */
@@ -71,8 +74,6 @@ medium_init(struct mem_allocator* allocator, struct sdis_medium** medium)
*medium = NULL;
}
-#define ENCLOSURE_MULTI_MEDIA UINT_MAX
-
struct enclosure {
struct s2d_scene_view* s2d_view;
struct s3d_scene_view* s3d_view;
@@ -97,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_MULTI_MEDIA;
+ enc->medium_id = ENCLOSURE_ID_MULTI_MEDIA;
}
static INLINE void
@@ -229,25 +230,31 @@ scene_get_interface
const unsigned iprim);
extern LOCAL_SYM res_T
-scene_get_medium
+scene_get_enclosure_id
(struct sdis_scene* scene,
const double position[],
- struct sdis_medium** medium);
+ unsigned* enclosure_id);
-/* This function assumes that the tested position lies into finite enclosure.
- * The medium into which it lies is thus retrieved by tracing a random ray
- * around the current position. For possible infinite enclosure, one has to use
- * the `scene_get_medium' function instead that, in counterpart, can be more
- * time consuming.
+/* This function assumes that the position under test lies within a finite
+ * enclosure. The enclosure in which it is located is therefore retrieved by
+ * tracing a random ray around the current position. For infinite enclosures,
+ * you need to use the `scene_get_enclosure_id' function, which in turn may take
+ * longer.
*
- * Note that actually, the function internally calls scene_get_medium if no
- * valid medium is found with the regular procedure. This may be due to
- * numerical issues or wrong assumptions on the current medium (its boundaries
- * are opened to infinity). */
+ * Note that the function actually calls scene_get_enclosure internally if no
+ * valid enclosure is found with the normal procedure. This may be due to
+ * numerical problems or incorrect assumptions about the current enclosure (its
+ * limits are open to infinity). */
extern LOCAL_SYM res_T
-scene_get_medium_in_closed_boundaries
- (struct sdis_scene* scn,
+scene_get_enclosure_id_in_closed_boundaries
+ (struct sdis_scene* scene,
const double position[],
+ unsigned* enclosure_id);
+
+extern LOCAL_SYM res_T
+scene_get_enclosure_medium
+ (struct sdis_scene* scene,
+ const struct enclosure* enclosure,
struct sdis_medium** medium);
extern LOCAL_SYM res_T
diff --git a/src/sdis_solve_camera.c b/src/sdis_solve_camera.c
@@ -504,12 +504,15 @@ sdis_solve_camera
/* Stardis variables */
struct sdis_estimator_buffer* buf = NULL;
- struct sdis_medium* medium = 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 */
+ unsigned enc_id = ENCLOSURE_ID_NULL;
+
/* Miscellaneous */
size_t ntiles_x, ntiles_y, ntiles, ntiles_adjusted;
size_t ntiles_proc; /* #tiles for the current proc */
@@ -537,10 +540,11 @@ sdis_solve_camera
}
/* Retrieve the medium in which the submitted position lies */
- res = scene_get_medium(scn, args->cam->position, &medium);
+ res = scene_get_enclosure_id(scn, args->cam->position, &enc_id);
+ res = scene_get_enclosure_medium(scn, scene_get_enclosure(scn, enc_id), &mdm);
if(res != RES_OK) goto error;
- if(medium->type != SDIS_FLUID) {
+ if(mdm->type != SDIS_FLUID) {
log_err(scn->dev,
"%s: the camera position `%g %g %g' must be in a fluid medium.\n",
FUNC_NAME, SPLIT3(args->cam->position));
@@ -658,7 +662,7 @@ sdis_solve_camera
/* Draw the tile */
res_local = solve_tile
- (scn, rng, medium, args->cam, args->time_range, tile_org, tile_sz,
+ (scn, rng, mdm, 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_probe_Xd.h b/src/sdis_solve_probe_Xd.h
@@ -116,7 +116,10 @@ XD(solve_one_probe)
struct accum* acc_temp,
struct accum* acc_time)
{
- struct sdis_medium* medium = NULL; /* Medium in which the probe lies */
+ /* 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;
ASSERT(scn && rng && check_solve_probe_args(args) == RES_OK);
@@ -126,7 +129,9 @@ XD(solve_one_probe)
*acc_time = ACCUM_NULL;
/* Retrieve the medium in which the submitted position lies */
- res = scene_get_medium(scn, args->position, &medium);
+ 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) {
@@ -144,7 +149,7 @@ XD(solve_one_probe)
/* Run a realisation */
realis_args.rng = rng;
- realis_args.medium = medium;
+ realis_args.medium = mdm;
realis_args.time = time;
realis_args.picard_order = args->picard_order;
realis_args.irealisation = irealisation;
@@ -191,7 +196,7 @@ XD(solve_probe)
size_t nthreads = 0;
/* Stardis variables */
- struct sdis_medium* medium = NULL;
+ struct sdis_medium* mdm = NULL;
struct sdis_estimator* estimator = NULL;
struct sdis_green_function* green = NULL;
struct sdis_green_function** per_thread_green = NULL;
@@ -200,6 +205,9 @@ XD(solve_probe)
struct ssp_rng_proxy* rng_proxy = NULL;
struct ssp_rng** per_thread_rng = NULL;
+ /* Enclosure in which the probe lies */
+ unsigned enc_id = ENCLOSURE_ID_NULL;
+
/* Miscellaneous */
struct accum* per_thread_acc_temp = NULL;
struct accum* per_thread_acc_time = NULL;
@@ -256,7 +264,9 @@ XD(solve_probe)
if(!per_thread_acc_time) { res = RES_MEM_ERR; goto error; }
/* Retrieve the medium in which the submitted position lies */
- res = scene_get_medium(scn, args->position, &medium);
+ 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;
/* Create the per thread green function */
@@ -327,7 +337,7 @@ XD(solve_probe)
/* Invoke the probe realisation */
realis_args.rng = rng;
- realis_args.medium = medium;
+ realis_args.medium = mdm;
realis_args.time = time;
realis_args.picard_order = args->picard_order;
realis_args.green_path = pgreen_path;