commit 684fa27d7043378380aab0f1a4ed5f05d3e8eaca
parent bbd08b5fb3c7b9998c30c38dce1de44efb312e70
Author: Benjamin Piaud <benjamin.piaud@meso-star.com>
Date: Thu, 21 Jan 2021 15:37:40 +0100
Add support of thermal contact resistance
Diffstat:
5 files changed, 66 insertions(+), 14 deletions(-)
diff --git a/src/sdis.h b/src/sdis.h
@@ -221,11 +221,16 @@ struct sdis_interface_shader {
* [0 convection_coef_upper_bound] */
double convection_coef_upper_bound;
+ /* May be NULL for solid/fluid or if the thermal contact resistance is 0 onto
+ * the whole interface. */
+ sdis_interface_getter_T thermal_contact_resistance; /* In K.m^2.W^-1 */
+
struct sdis_interface_side_shader front;
struct sdis_interface_side_shader back;
};
#define SDIS_INTERFACE_SHADER_NULL__ \
- {NULL, 0, SDIS_INTERFACE_SIDE_SHADER_NULL__, SDIS_INTERFACE_SIDE_SHADER_NULL__}
+ {NULL, 0, NULL, SDIS_INTERFACE_SIDE_SHADER_NULL__, \
+ SDIS_INTERFACE_SIDE_SHADER_NULL__}
static const struct sdis_interface_shader SDIS_INTERFACE_SHADER_NULL =
SDIS_INTERFACE_SHADER_NULL__;
diff --git a/src/sdis_heat_path_boundary_Xd.h b/src/sdis_heat_path_boundary_Xd.h
@@ -478,6 +478,7 @@ XD(solid_solid_boundary_path)
double tmp;
double r;
double power;
+ double tcr;
float dir0[DIM], dir1[DIM], dir2[DIM], dir3[DIM];
float dir_front[DIM], dir_back[DIM];
float* dir;
@@ -501,6 +502,9 @@ XD(solid_solid_boundary_path)
ASSERT(solid_front->type == SDIS_SOLID);
ASSERT(solid_back->type == SDIS_SOLID);
+ /* Retrieve the thermal contact resistance */
+ tcr = interface_get_thermal_contact_resistance(interf, frag);
+
/* Fetch the properties of the media */
lambda_front = solid_get_thermal_conductivity(solid_front, &rwalk->vtx);
lambda_back = solid_get_thermal_conductivity(solid_back, &rwalk->vtx);
@@ -560,19 +564,44 @@ XD(solid_solid_boundary_path)
goto error;
}
- /* Define the reinjection side. Note that the proba should be :
- * Lf/Df' / (Lf/Df' + Lb/Db')
- *
- * with L<f|b> the lambda of the <front|back> side and D<f|b>' the adjusted
- * delta of the <front|back> side, i.e. :
- * D<f|b>' = reinject_dst_<front|back> / sqrt(DIM)
- *
- * Anyway, one can avoid to compute the adjusted delta by directly using the
- * adjusted reinjection distance since the resulting proba is strictly the
- * same; sqrt(DIM) can be simplified. */
r = ssp_rng_canonical(rng);
- proba = (lambda_front/reinject_dst_front)
- / (lambda_front/reinject_dst_front + lambda_back/reinject_dst_back);
+ if(tcr == 0) { /* No thermal contact resistance */
+ /* Define the reinjection side. Note that the proba should be : Lf/Df' /
+ * (Lf/Df' + Lb/Db')
+ *
+ * with L<f|b> the lambda of the <front|back> side and D<f|b>' the adjusted
+ * delta of the <front|back> side, i.e. : D<f|b>' =
+ * reinject_dst_<front|back> / sqrt(DIM)
+ *
+ * Anyway, one can avoid to compute the adjusted delta by directly using the
+ * adjusted reinjection distance since the resulting proba is strictly the
+ * same; sqrt(DIM) can be simplified. */
+ proba = (lambda_front/reinject_dst_front)
+ / (lambda_front/reinject_dst_front + lambda_back/reinject_dst_back);
+ } else {
+ const double df = reinject_dst_front/sqrt(DIM);
+ const double db = reinject_dst_back/sqrt(DIM);
+ const double tmp_front = lambda_front/df;
+ const double tmp_back = lambda_back/db;
+ const double tmp_r = tcr*tmp_front*tmp_back;
+ switch(rwalk->hit_side) {
+ case SDIS_BACK:
+ /* When coming from the BACK side, the probability to be reinjected on
+ * the FRONT side depends on the thermal contact resistance: it decrases
+ * when the TCR increases (and tends to 0 when TCR -> +\infty) */
+ proba = (tmp_front) / (tmp_front + tmp_back + tmp_r);
+ break;
+ case SDIS_FRONT:
+ /* Same thing when coming from the FRONT side: the probability of
+ * reinjection on the FRONT side depends on the thermal contact
+ * resistance: it increases when the TCR increases (and tends to 1 when
+ * the TCR -> +\infty) */
+ proba = (tmp_front + tmp_r) / (tmp_front + tmp_back + tmp_r);
+ break;
+ default: FATAL("Unreachable code.\n"); break;
+ }
+ }
+
if(r < proba) { /* Reinject in front */
dir = dir_front;
hit = &hit0;
diff --git a/src/sdis_interface.c b/src/sdis_interface.c
@@ -56,7 +56,6 @@ check_interface_shader
" shader's pointer function for this attribute should be NULL.\n",
caller_name);
}
-
if(shader->convection_coef_upper_bound < 0) {
log_warn(dev,
"%s: Invalid upper bound for convection coefficient (%g).\n",
@@ -64,6 +63,14 @@ check_interface_shader
if(type[0] == SDIS_FLUID || type[1] == SDIS_FLUID) return 0;
}
+ if((type[0] != SDIS_SOLID || type[1] != SDIS_SOLID)
+ && shader->thermal_contact_resistance) {
+ log_warn(dev,
+ "%s: only solid/solid interface can have a thermal contact resistance. The "
+ " shader's pointer function for this attribute should be NULL.\n",
+ caller_name);
+ }
+
FOR_EACH(i, 0, 2) {
switch(type[i]) {
case SDIS_SOLID:
diff --git a/src/sdis_interface_c.h b/src/sdis_interface_c.h
@@ -89,6 +89,16 @@ interface_get_convection_coef
}
static INLINE double
+interface_get_thermal_contact_resistance
+ (const struct sdis_interface* interf,
+ const struct sdis_interface_fragment* frag)
+{
+ ASSERT(interf && frag);
+ return interf->shader.thermal_contact_resistance
+ ? interf->shader.thermal_contact_resistance(frag, interf->data) : 0;
+}
+
+static INLINE double
interface_get_convection_coef_upper_bound
(const struct sdis_interface* interf)
{
diff --git a/src/test_sdis_utils.h b/src/test_sdis_utils.h
@@ -191,6 +191,7 @@ static const struct sdis_fluid_shader DUMMY_FLUID_SHADER = {
static const struct sdis_interface_shader DUMMY_INTERFACE_SHADER = {
dummy_interface_getter,
0,
+ dummy_interface_getter,
DUMMY_INTERFACE_SIDE_SHADER__,
DUMMY_INTERFACE_SIDE_SHADER__
};
