commit 221ea417664665c9d3411d6c5cdb168c4f8bc687
parent 3ee7be850308d343ea9d8028e63f12e40bd51e6f
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Fri, 5 Nov 2021 16:47:57 +0100
Change the sdis API that controls the picard order
Diffstat:
6 files changed, 82 insertions(+), 57 deletions(-)
diff --git a/src/sdis.h b/src/sdis.h
@@ -386,9 +386,8 @@ struct sdis_scene_create_args {
/* Min/max temperature used to linearise the radiative temperature */
double t_range[2];
- /* Maximum number of heat path branchings. Actually the Picard order
- * used to estimate the radiative temperature is max_branchings+1 */
- size_t max_branchings;
+ /* Picard order used to estimate the radiative temperature */
+ size_t picard_order;
};
#define SDIS_SCENE_CREATE_ARGS_DEFAULT__ { \
@@ -401,7 +400,7 @@ struct sdis_scene_create_args {
1.0, /* #Floating point to meter scale factor */ \
SDIS_AMBIENT_RADIATIVE_TEMPERATURE_NULL__,/* Ambient radiative temperature */\
{0.0, -1.0}, /* Temperature range */ \
- 0 /* Maximum branchings */ \
+ 1 /* Picard order */ \
}
static const struct sdis_scene_create_args SDIS_SCENE_CREATE_ARGS_DEFAULT =
SDIS_SCENE_CREATE_ARGS_DEFAULT__;
@@ -867,21 +866,19 @@ sdis_scene_set_temperature_range
(struct sdis_scene* scn,
const double t_range[2]);
-/* Get the maximum number of branchings for the sampled paths. */
+/* Get the picard recursion order. */
SDIS_API res_T
-sdis_scene_get_max_branchings
+sdis_scene_get_picard_order
(const struct sdis_scene* scn,
- size_t* max_branchings);
+ size_t* picard_order);
-/* Set the maximum number of branchings for the sampled paths. A value greater
- * than zero enables the estimation of T4 radiative transfer by using the
- * Picard's algorithm. In others words, 'max_branchings+1' is the order of the
- * Picard's recursion: at order one (i.e. max_branchings==0), the radiative
- * transfer is linearised */
-SDIS_API res_T
-sdis_scene_set_max_branchings
+/* Set the Picard recursion order to estimate the radiative temperature. An
+ * order of one means that the radiative temperature is linearized, while
+ * higher orders allow the estimation of the T4 radiative transfer. */
+ SDIS_API res_T
+sdis_scene_set_picard_order
(struct sdis_scene* scn,
- const size_t max_branchings);
+ const size_t picard_order);
/* Search the point onto the scene geometry that is the closest of `pos'. The
* `radius' parameter controls the maximum search distance around `pos'. The
diff --git a/src/sdis_scene.c b/src/sdis_scene.c
@@ -243,22 +243,22 @@ sdis_scene_set_temperature_range
}
res_T
-sdis_scene_get_max_branchings
+sdis_scene_get_picard_order
(const struct sdis_scene* scn,
- size_t* max_branchings)
+ size_t* picard_order)
{
- if(!scn || !max_branchings) return RES_BAD_ARG;
- *max_branchings = scn->max_branchings;
+ if(!scn || !picard_order) return RES_BAD_ARG;
+ *picard_order = scn->max_branchings+1;
return RES_OK;
}
res_T
-sdis_scene_set_max_branchings
+sdis_scene_set_picard_order
(struct sdis_scene* scn,
- const size_t max_branchings)
+ const size_t picard_order)
{
- if(!scn) return RES_BAD_ARG;
- scn->max_branchings = max_branchings;
+ if(!scn || picard_order < 1) return RES_BAD_ARG;
+ scn->max_branchings = scene_get_picard_order(scn);
return RES_OK;
}
diff --git a/src/sdis_scene_Xd.h b/src/sdis_scene_Xd.h
@@ -120,7 +120,8 @@ check_sdis_scene_create_args(const struct sdis_scene_create_args* args)
&& args->nprimitives < UINT_MAX
&& args->nvertices
&& args->nvertices < UINT_MAX
- && args->fp_to_meter > 0;
+ && args->fp_to_meter > 0
+ && args->picard_order > 0;
}
#endif /* SDIS_SCENE_XD_H */
@@ -916,7 +917,7 @@ XD(scene_create)
scn->tmin = args->t_range[0];
scn->tmax = args->t_range[1];
scn->outer_enclosure_id = UINT_MAX;
- scn->max_branchings = args->max_branchings;
+ scn->max_branchings = args->picard_order - 1;
darray_interf_init(dev->allocator, &scn->interfaces);
darray_medium_init(dev->allocator, &scn->media);
darray_prim_prop_init(dev->allocator, &scn->prim_props);
diff --git a/src/sdis_scene_c.h b/src/sdis_scene_c.h
@@ -303,5 +303,13 @@ scene_is_2d(const struct sdis_scene* scn)
return scn->s2d_view != NULL;
}
+static FINLINE size_t
+scene_get_picard_order(const struct sdis_scene* scn)
+{
+ ASSERT(scn);
+ return scn->max_branchings+1;
+}
+
+
#endif /* SDIS_SCENE_C_H */
diff --git a/src/test_sdis_scene.c b/src/test_sdis_scene.c
@@ -139,6 +139,9 @@ test_scene_3d(struct sdis_device* dev, struct sdis_interface* interf)
scn_args.fp_to_meter = 0;
BA(sdis_scene_create(dev, &scn_args, &scn));
scn_args.fp_to_meter = 1;
+ scn_args.picard_order = 0;
+ BA(sdis_scene_create(dev, &scn_args, &scn));
+ scn_args.picard_order = 1;
/* Duplicated vertex */
ctx.positions = duplicated_vertices;
ctx.indices = dup_vrtx_indices;
@@ -271,6 +274,7 @@ test_scene_2d(struct sdis_device* dev, struct sdis_interface* interf)
size_t i;
size_t iprim;
size_t dup_vrtx_indices[] = { 0, 1 };
+ size_t picard_order;
enum sdis_scene_dimension dim;
ctx.positions = square_vertices;
@@ -305,6 +309,9 @@ test_scene_2d(struct sdis_device* dev, struct sdis_interface* interf)
scn_args.fp_to_meter = 0;
BA(sdis_scene_2d_create(dev, &scn_args, &scn));
scn_args.fp_to_meter = 1;
+ scn_args.picard_order = 0;
+ BA(sdis_scene_create(dev, &scn_args, &scn));
+ scn_args.picard_order = 1;
/* Duplicated vertex */
ctx.positions = duplicated_vertices;
ctx.indices = dup_vrtx_indices;
@@ -389,6 +396,17 @@ test_scene_2d(struct sdis_device* dev, struct sdis_interface* interf)
CHK(t_range[0] == 1);
CHK(t_range[1] == 100);
+ BA(sdis_scene_get_picard_order(NULL, &picard_order));
+ BA(sdis_scene_get_picard_order(scn, NULL));
+ OK(sdis_scene_get_picard_order(scn, &picard_order));
+ CHK(picard_order == SDIS_SCENE_CREATE_ARGS_DEFAULT.picard_order);
+ CHK(picard_order == 1);
+ OK(sdis_scene_set_picard_order(scn, 3));
+ OK(sdis_scene_get_picard_order(scn, &picard_order));
+ CHK(picard_order == 3);
+ BA(sdis_scene_set_picard_order(scn, 0));
+ OK(sdis_scene_set_picard_order(scn, 1));
+
BA(sdis_scene_get_boundary_position(NULL, 1, &u0, pos));
BA(sdis_scene_get_boundary_position(scn, 4, &u0, pos));
BA(sdis_scene_get_boundary_position(scn, 1, NULL, pos));
diff --git a/src/test_sdis_unstationary_atm.c b/src/test_sdis_unstationary_atm.c
@@ -89,7 +89,7 @@
#define X_PROBE (XH + 0.2 * XE)
-#define DELTA (XE/30.0)
+#define DELTA (XE/40.0)
/*******************************************************************************
* Box geometry
@@ -114,11 +114,11 @@ static const size_t model3d_nvertices = sizeof(model3d_vertices)/(sizeof(double)
* triangle.
* ,3---,4---,5 ,3----4----5 ,4
* ,' | ,' | ,'/| ,'/| \ | \ | ,'/|
- * 9----10---11 / | 9' / | \ | \ | 10 / | Y
- * |', |', | / ,2 | / ,0---,1---,2 | / ,1 |
+ * 9----10---11 / | 9' / | \ | \ | 10 / | Y
+ * |', |', | / ,2 | / ,0---,1---,2 | / ,1 |
* | ',| ',|/,' |/,' | ,' | ,' |/,' o--X
- * 6----7----8' 6----7'---8' 7 /
- * Front, right Back, left and Internal Z
+ * 6----7----8' 6----7'---8' 7 /
+ * Front, right Back, left and Internal Z
* and Top faces bottom faces face */
static const size_t model3d_indices[22/*#triangles*/*3/*#indices per triangle*/] = {
0, 3, 1, 1, 3, 4, 1, 4, 2, 2, 4, 5, /* -Z */
@@ -421,10 +421,11 @@ solve_tbound1
size_t nreals;
size_t nfails;
enum sdis_scene_dimension dim;
- double t[] = { 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
- double ref[sizeof(t) / sizeof(*t)]
- = { 290.046375, 289.903935, 289.840490, 289.802690, 289.777215,
- 289.759034, 289.745710, 289.735826, 289.728448, 289.722921 };
+ const double t[] = { 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
+ const double ref[sizeof(t) / sizeof(*t)] = {
+ 290.046375, 289.903935, 289.840490, 289.802690, 289.777215, 289.759034,
+ 289.745710, 289.735826, 289.728448, 289.722921
+ };
const int nsimuls = sizeof(t) / sizeof(*t);
int isimul;
ASSERT(scn && rng);
@@ -475,9 +476,8 @@ solve_tbound1
printf("Elapsed time = %s\n", dump);
printf("Time per realisation (in usec) = %g +/- %g\n\n", time.E, time.SE);
- CHK(nfails + nreals == N);
- CHK(nfails <= N/1000);
CHK(eq_eps(T.E, ref[isimul], EPS));
+ /*CHK(eq_eps(T.E, ref[isimul], T.SE*3));*/
OK(sdis_estimator_ref_put(estimator));
}
@@ -498,10 +498,11 @@ solve_tbound2
size_t nreals;
size_t nfails;
enum sdis_scene_dimension dim;
- double t[] = { 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
- double ref[sizeof(t) / sizeof(*t)]
- = { 309.08032, 309.34626, 309.46525, 309.53625, 309.58408,
- 309.618121, 309.642928, 309.661167, 309.674614, 309.684524 };
+ const double t[] = { 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
+ const double ref[sizeof(t) / sizeof(*t)] = {
+ 309.08032, 309.34626, 309.46525, 309.53625, 309.58408, 309.618121,
+ 309.642928, 309.661167, 309.674614, 309.684524
+ };
const int nsimuls = sizeof(t) / sizeof(*t);
int isimul;
ASSERT(scn && rng);
@@ -555,6 +556,7 @@ solve_tbound2
CHK(nfails + nreals == N);
CHK(nfails <= N/1000);
CHK(eq_eps(T.E, ref[isimul], EPS));
+ /*CHK(eq_eps(T.E, ref[isimul], T.SE*3));*/
OK(sdis_estimator_ref_put(estimator));
}
@@ -574,10 +576,11 @@ solve_tsolid
size_t nreals;
size_t nfails;
enum sdis_scene_dimension dim;
- double t[] = { 0, 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
- double ref[sizeof(t) / sizeof(*t)]
- = { 300, 300.87408, 302.25832, 303.22164, 303.89954, 304.39030,
- 304.75041, 305.01595, 305.21193, 305.35641, 305.46271 };
+ const double t[] = { 0, 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
+ const double ref[sizeof(t) / sizeof(*t)] = {
+ 300, 300.87408, 302.25832, 303.22164, 303.89954, 304.39030, 304.75041,
+ 305.01595, 305.21193, 305.35641, 305.46271
+ };
const int nsimuls = sizeof(t) / sizeof(*t);
int isimul;
ASSERT(scn && rng);
@@ -621,6 +624,7 @@ solve_tsolid
CHK(nfails + nreals == N);
CHK(nfails <= N / 1000);
CHK(eq_eps(T.E, ref[isimul], EPS));
+ /*CHK(eq_eps(T.E, ref[isimul], T.SE*3));*/
OK(sdis_estimator_ref_put(estimator));
}
@@ -640,10 +644,11 @@ solve_tfluid
size_t nfails;
enum sdis_scene_dimension dim;
double eps;
- double t[] = { 0, 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
- double ref[sizeof(t) / sizeof(*t)]
- = { 300, 309.53905, 309.67273, 309.73241, 309.76798, 309.79194, 309.80899,
- 309.82141, 309.83055, 309.83728, 309.84224 };
+ const double t[] = { 0, 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 };
+ const double ref[sizeof(t) / sizeof(*t)] = {
+ 300, 309.53905, 309.67273, 309.73241, 309.76798, 309.79194, 309.80899,
+ 309.82141, 309.83055, 309.83728, 309.84224
+ };
const int nsimuls = sizeof(t) / sizeof(*t);
int isimul;
ASSERT(scn);
@@ -834,7 +839,7 @@ main(int argc, char** argv)
model3d_interfaces[10] = interf_TA;
model3d_interfaces[11] = interf_TA;
/* Top */
- model3d_interfaces[12] = interf_adiabatic_1;
+ model3d_interfaces[12] = interf_adiabatic_1;
model3d_interfaces[13] = interf_adiabatic_1;
model3d_interfaces[14] = interf_adiabatic_2;
model3d_interfaces[15] = interf_adiabatic_2;
@@ -869,10 +874,8 @@ main(int argc, char** argv)
scn_args.context = model3d_interfaces;
scn_args.trad.temperature = TR;
scn_args.trad.reference = TR;
- scn_args.tmax = MMAX(T0_FLUID, T0_SOLID);
- scn_args.tmax = MMAX(scn_args.tmax, TA);
- scn_args.tmax = MMAX(scn_args.tmax, TG);
- scn_args.tmax = MMAX(scn_args.tmax, TR);
+ scn_args.t_range[0] = MMIN(MMIN(MMIN(MMIN(T0_FLUID, T0_SOLID), TA), TG), TR);
+ scn_args.t_range[1] = MMAX(MMAX(MMAX(MMAX(T0_FLUID, T0_SOLID), TA), TG), TR);
OK(sdis_scene_create(dev, &scn_args, &box_scn));
/* Create the square scene */
@@ -884,10 +887,8 @@ main(int argc, char** argv)
scn_args.context = model2d_interfaces;
scn_args.trad.temperature = TR;
scn_args.trad.reference = TR;
- scn_args.tmax = MMAX(T0_FLUID, T0_SOLID);
- scn_args.tmax = MMAX(scn_args.tmax, TA);
- scn_args.tmax = MMAX(scn_args.tmax, TG);
- scn_args.tmax = MMAX(scn_args.tmax, TR);
+ scn_args.t_range[0] = MMIN(MMIN(MMIN(MMIN(T0_FLUID, T0_SOLID), TA), TG), TR);
+ scn_args.t_range[1] = MMAX(MMAX(MMAX(MMAX(T0_FLUID, T0_SOLID), TA), TG), TR);
OK(sdis_scene_2d_create(dev, &scn_args, &square_scn));
/* Release the interfaces */
