commit fa403f13a85230402ffb965f0968bb4821aa937c
parent 38603e8f489d1e9c3398a682d0d7979ddc5170c0
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Thu, 12 Apr 2018 15:38:20 +0200
Implement the computation of unknown fluid temperature
The computation actually assumes that the convection coefficient of a
fluid enclosure is constant for the whole enclosure.
Diffstat:
3 files changed, 174 insertions(+), 19 deletions(-)
diff --git a/src/sdis_scene_Xd.h b/src/sdis_scene_Xd.h
@@ -518,6 +518,7 @@ XD(setup_enclosure_geometry)(struct sdis_scene* scn, struct sencXd(enclosure)* e
struct sXd(vertex_data) vdata = SXD_VERTEX_DATA_NULL;
struct enclosure enc_dummy;
struct enclosure* enc_data;
+ float S, V;
unsigned iprim, nprims, nverts;
#if DIM == 2
const struct enclosure2d_header* header;
@@ -571,10 +572,22 @@ XD(setup_enclosure_geometry)(struct sdis_scene* scn, struct sencXd(enclosure)* e
CALL(sXd(mesh_setup_indexed_vertices)(sXd_shape, nprims, XD(enclosure_indices),
nverts, &vdata, 1, enc));
#endif
-
CALL(sXd(scene_create)(sXd_dev, &sXd_scn));
CALL(sXd(scene_attach_shape)(sXd_scn, sXd_shape));
CALL(sXd(scene_view_create)(sXd_scn, SXD_SAMPLE, &enc_data->sXd(view)));
+
+ /* Compute the S/V ratio */
+#if DIM == 2
+ CALL(sXd(scene_view_compute_contour_length)(enc_data->sXd(view), &S));
+ CALL(sXd(scene_view_compute_area)(enc_data->sXd(view), &V));
+#else
+ CALL(sXd(scene_view_compute_area)(enc_data->sXd(view), &S));
+ CALL(sXd(scene_view_compute_volume)(enc_data->sXd(view), &V));
+#endif
+ /* The volume of the enclosure is actually negative since Star-Enc ensures
+ * that the normal of its primitives point outward the enclosure. Take its
+ * absolute value in order to ensure a postive value. */
+ enc_data->S_over_V = S/absf(V);
#undef CALL
/* Define the identifier of the enclosure primitives in the whole scene */
@@ -600,8 +613,8 @@ error:
goto exit;
}
-/* Build the Star-XD scene view and define its associated data of the fluid
- * enclosures */
+/* Build the Star-XD scene view and define its associated data of the finite
+ * fluid enclosures */
static res_T
XD(setup_enclosures)(struct sdis_scene* scn, struct sencXd(descriptor)* desc)
{
@@ -626,8 +639,8 @@ XD(setup_enclosures)(struct sdis_scene* scn, struct sencXd(descriptor)* desc)
mdm = darray_medium_cdata_get(&scn->media)[header->enclosed_medium];
ASSERT(mdm);
- /* Silently discard the solid enclosures */
- if(mdm->type == SDIS_MEDIUM_FLUID) {
+ /* Silently discard the solid and infinite enclosures */
+ if(mdm->type == SDIS_MEDIUM_FLUID && !header->is_infinite) {
res = XD(setup_enclosure_geometry)(scn, enc);
if(res != RES_OK) goto error;
}
diff --git a/src/sdis_scene_c.h b/src/sdis_scene_c.h
@@ -60,6 +60,8 @@ struct enclosure {
/* Map the id of the enclosure primitives to their primitive id into the
* whole scene */
struct darray_uint local2global;
+
+ double S_over_V; /* in 3D = surface/volume; in 2D = perimeter/area */
};
static INLINE void
@@ -69,6 +71,7 @@ enclosure_init(struct mem_allocator* allocator, struct enclosure* enc)
enc->s2d_view = NULL;
enc->s3d_view = NULL;
darray_uint_init(allocator, &enc->local2global);
+ enc->S_over_V = 0;
}
static INLINE void
@@ -90,6 +93,7 @@ enclosure_copy(struct enclosure* dst, const struct enclosure* src)
S2D(scene_view_ref_get(src->s2d_view));
dst->s2d_view = src->s2d_view;
}
+ dst->S_over_V = src->S_over_V;
return darray_uint_copy(&dst->local2global, &src->local2global);
}
@@ -109,6 +113,7 @@ enclosure_copy_and_release(struct enclosure* dst, struct enclosure* src)
dst->s2d_view = src->s2d_view;
src->s2d_view = NULL;
}
+ dst->S_over_V = src->S_over_V;
return RES_OK;
}
@@ -173,6 +178,27 @@ scene_get_medium
const double position[],
const struct sdis_medium** medium);
+static INLINE void
+scene_get_enclosure_ids
+ (const struct sdis_scene* scn,
+ const unsigned iprim,
+ unsigned encs[2]) /* Front and Back enclosure identifiers */
+{
+ ASSERT(scn && iprim < darray_prim_prop_size_get(&scn->prim_props));
+ ASSERT(encs);
+ encs[0] = darray_prim_prop_cdata_get(&scn->prim_props)[iprim].front_enclosure;
+ encs[1] = darray_prim_prop_cdata_get(&scn->prim_props)[iprim].back_enclosure;
+}
+
+static INLINE const struct enclosure*
+scene_get_enclosure(struct sdis_scene* scn, const unsigned ienc)
+{
+ const struct enclosure* enc = NULL;
+ ASSERT(scn);
+ enc = htable_enclosure_find(&scn->enclosures, &ienc);
+ return enc;
+}
+
static FINLINE int
scene_is_2d(const struct sdis_scene* scn)
{
diff --git a/src/sdis_solve_Xd.h b/src/sdis_solve_Xd.h
@@ -101,7 +101,7 @@ static const struct XD(rwalk) XD(RWALK_NULL) = {
struct XD(temperature) {
res_T (*func)/* Next function to invoke in order to compute the temperature */
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
@@ -114,7 +114,7 @@ static const struct XD(temperature) XD(TEMPERATURE_NULL) = { NULL, 0, 0 };
static res_T
XD(boundary_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
@@ -123,7 +123,7 @@ XD(boundary_temperature)
static res_T
XD(solid_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
@@ -132,7 +132,7 @@ XD(solid_temperature)
static res_T
XD(fluid_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
@@ -141,7 +141,7 @@ XD(fluid_temperature)
static res_T
XD(radiative_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
@@ -190,7 +190,7 @@ XD(check_rwalk_fragment_consistency)
static res_T
XD(trace_radiative_path)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const float ray_dir[3],
const double fp_to_meter,
const struct rwalk_context* ctx,
@@ -308,7 +308,7 @@ error:
res_T
XD(radiative_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
@@ -352,26 +352,142 @@ error:
res_T
XD(fluid_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
struct ssp_rng* rng,
struct XD(temperature)* T)
{
+ struct sXd(attrib) attr_P, attr_N;
+ struct sdis_interface_fragment frag;
+ const struct sdis_interface* interf;
+ const struct enclosure* enc;
+ unsigned enc_ids[2];
+ unsigned enc_id;
+ double rho; /* Volumic mass */
+ double hc; /* Convection coef */
+ double cp; /* Calorific capacity */
+ double mu;
+ double tau;
double tmp;
+#if DIM == 2
+ float st;
+#else
+ float st[2];
+#endif
(void)rng, (void)fp_to_meter, (void)ctx;
ASSERT(scn && fp_to_meter > 0 && ctx && rwalk && rng && T);
ASSERT(rwalk->mdm->type == SDIS_MEDIUM_FLUID);
tmp = fluid_get_temperature(rwalk->mdm, &rwalk->vtx);
- if(tmp < 0) {
- log_err(scn->dev, "%s: unknown fluid temperature is not supported.\n",
+ if(tmp >= 0) {
+ T->value += tmp;
+ T->done = 1;
+ return RES_OK;
+ }
+
+ if(SXD_HIT_NONE(&rwalk->hit)) { /* The path begins in the fluid */
+ const float range[2] = {0, FLT_MAX};
+ float dir[DIM] = {0};
+ float org[DIM];
+
+ dir[DIM-1] = 1;
+ fX_set_dX(org, rwalk->vtx.P);
+
+ /* Init the path hit field required to define the current enclosure and
+ * fetch the interface data */
+ SXD(scene_view_trace_ray(scn->sXd(view), org, dir, range, NULL, &rwalk->hit));
+
+ if(SXD_HIT_NONE(&rwalk->hit)) {
+ log_err(scn->dev,
+"%s: the position %g %g %g lise in the surrounding fluid whose temperature must \n"
+"be known.\n",
+ FUNC_NAME, SPLIT3(rwalk->vtx.P));
+ return RES_BAD_OP;
+ }
+ }
+
+ /* Setup the fragment of the interface */
+ XD(setup_interface_fragment)(&frag, &rwalk->vtx, &rwalk->hit);
+
+ /* Fetch the current interface and its associated enclosures */
+ interf = scene_get_interface(scn, rwalk->hit.prim.prim_id);
+ scene_get_enclosure_ids(scn, rwalk->hit.prim.prim_id, enc_ids);
+
+ /* Define 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];
+ } else {
+ ASSERT(rwalk->mdm == interf->medium_back);
+ enc_id = enc_ids[1];
+ }
+
+ /* Fetch the enclosure data */
+ enc = scene_get_enclosure(scn, enc_id);
+ if(!enc) {
+ log_err(scn->dev,
+"%s: invalid enclosure. The position %g %g %g may lie in the surrounding fluid.\n",
+ FUNC_NAME, SPLIT3(rwalk->vtx.P));
+ return RES_BAD_OP;
+ }
+
+ /* Fetch the physical properties */
+ hc = interface_get_convection_coef(interf, &frag);
+ cp = fluid_get_calorific_capacity(rwalk->mdm, &rwalk->vtx);
+ rho = fluid_get_volumic_mass(rwalk->mdm, &rwalk->vtx);
+
+ /* Sample the time.
+ * FIXME we assume that hc is constant for the whole enclosure */
+ mu = hc / (rho * cp) * enc->S_over_V;
+ tau = ssp_ran_exp(rng, mu);
+ rwalk->vtx.time = MMIN(rwalk->vtx.time - tau, 0);
+
+ /* Check the initial condition */
+ tmp = fluid_get_temperature(rwalk->mdm, &rwalk->vtx);
+ if(tmp >= 0) {
+ T->value += tmp;
+ T->done = 1;
+ return RES_OK;
+ }
+
+ /* The initial condition should be reached */
+ if(rwalk->vtx.time <=0) {
+ log_err(scn->dev,
+ "%s: undefined initial condition. "
+ "The time is null but the temperature remains unknown.\n",
FUNC_NAME);
return RES_BAD_OP;
}
- T->value += tmp;
- T->done = 1;
+
+ /* Uniformly sample the enclosure */
+#if DIM == 2
+ SXD(scene_view_sample
+ (enc->sXd(view),
+ ssp_rng_canonical_float(rng),
+ ssp_rng_canonical_float(rng),
+ &rwalk->hit.prim,
+ &rwalk->hit.u));
+ st = rwalk->hit.u;
+#else
+ SXD(scene_view_sample
+ (enc->sXd(view),
+ ssp_rng_canonical_float(rng),
+ ssp_rng_canonical_float(rng),
+ ssp_rng_canonical_float(rng),
+ &rwalk->hit.prim,
+ rwalk->hit.uv));
+ f2_set(st, rwalk->hit.uv);
+#endif
+ rwalk->hit.distance = 0;
+ SXD(primitive_get_attrib(&rwalk->hit.prim, SXD_POSITION, st, &attr_P));
+ SXD(primitive_get_attrib(&rwalk->hit.prim, SXD_GEOMETRY_NORMAL, st, &attr_N));
+ dX_set_fX(rwalk->vtx.P, attr_P.value);
+ fX(set)(rwalk->hit.normal, attr_N.value);
+
+ T->func = XD(boundary_temperature);
+ rwalk->mdm = NULL; /* The random walk is at an interface between 2 media */
return RES_OK;
}
@@ -530,7 +646,7 @@ XD(solid_fluid_boundary_temperature)
res_T
XD(boundary_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
@@ -574,7 +690,7 @@ XD(boundary_temperature)
res_T
XD(solid_temperature)
- (const struct sdis_scene* scn,
+ (struct sdis_scene* scn,
const double fp_to_meter,
const struct rwalk_context* ctx,
struct XD(rwalk)* rwalk,
