commit 086bfbe04ddee1765bdfe59541b5377f56bb13d3
parent b07f9bf0e5b655f3458e2ab4fb0c81a2bec31df3
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Mon, 18 Dec 2023 12:27:30 +0100
Start implementing external flux calculation
Add a new functor on the side of an interface used to sample an external
source. It returns a sampled direction, its pdf and the radiance of the
sampled source. The caller can define as many sources as required, as
long as it is able to provide a sampling function. External sources
therefore boil down to this single functor. Or to put it another way,
the data from the external sources are external to the solver.
We then call a new function on the solid/fluid interface when linear
radiative exchange is involved (i.e. a picard order of 1) to update the
Monte Carlo weight of the sampled coupled trajectory. Non-linear
radiative transfer is not supported for the same reasons as other
sources: data accumulation along the trajectory is impossible with
Picard recursion. In any case, implementation of this function has only
just begun, with a first draft of the direct incident flux calculation.
In short, this commitment is only a first step.
Diffstat:
7 files changed, 237 insertions(+), 6 deletions(-)
diff --git a/src/sdis.h b/src/sdis.h
@@ -171,15 +171,35 @@ enum sdis_medium_type {
* medium. */
typedef double
(*sdis_medium_getter_T)
- (const struct sdis_rwalk_vertex* vert,
- struct sdis_data* data);
+ (const struct sdis_rwalk_vertex* vert, /* Medium position */
+ struct sdis_data* data); /* User data */
/* Functor type used to retrieve the spatio temporal physical properties of an
* interface. */
typedef double
(*sdis_interface_getter_T)
- (const struct sdis_interface_fragment* frag,
- struct sdis_data* data);
+ (const struct sdis_interface_fragment* frag, /* Interface position */
+ struct sdis_data* data); /* User data */
+
+/* A sample of an external source */
+struct sdis_external_sources_sample {
+ double dir[3]; /* Direction _to_ the source */
+ double pdf; /* Pdf of sampled direction */
+ double distance; /* Distance to the source [m] */
+ double radiance; /* Constant source radiance [W/m^2/sr] */
+};
+#define SDIS_EXTERNAL_SOURCES_SAMPLE_NULL__ {{0,0,0}, 0, 0, 0}
+static const struct sdis_external_sources_sample
+SDIS_EXTERNAL_SOURCES_SAMPLE_NULL = SDIS_EXTERNAL_SOURCES_SAMPLE_NULL__;
+
+/* Functor for sampling external sources. The returned sample is used to
+ * evaluate an external flux at the interface. */
+typedef res_T
+(*sdis_interface_sample_external_sources_T)
+ (const struct sdis_interface_fragment* frag, /* Interface position */
+ struct ssp_rng* rng, /* Random Number Generator to use */
+ struct sdis_external_sources_sample* sample, /* Returned sample */
+ struct sdis_data* data); /* User data */
/* Define the physical properties of a solid */
struct sdis_solid_shader {
@@ -198,6 +218,7 @@ struct sdis_solid_shader {
* unknown for the submitted random walk vertex.
* This getter is always called at time >= t0 (see below). */
sdis_medium_getter_T temperature;
+
/* The time until the initial condition is maintained for this solid;
* can neither be negative nor infinity, default is 0. */
double t0;
@@ -238,8 +259,12 @@ struct sdis_interface_side_shader {
/* Reference temperature used in Picard 1 */
sdis_interface_getter_T reference_temperature;
+
+ /* Manage external sources, i.e. sample them to evaluate the corresponding
+ * external flow. Can be NULL <=> no external source */
+ sdis_interface_sample_external_sources_T sample_external_sources;
};
-#define SDIS_INTERFACE_SIDE_SHADER_NULL__ { NULL, NULL, NULL, NULL, NULL }
+#define SDIS_INTERFACE_SIDE_SHADER_NULL__ { NULL, NULL, NULL, NULL, NULL, NULL }
static const struct sdis_interface_side_shader SDIS_INTERFACE_SIDE_SHADER_NULL =
SDIS_INTERFACE_SIDE_SHADER_NULL__;
diff --git a/src/sdis_heat_path_boundary.c b/src/sdis_heat_path_boundary.c
@@ -41,3 +41,7 @@
#include "sdis_heat_path_boundary_Xd.h"
#define SDIS_XD_DIMENSION 3
#include "sdis_heat_path_boundary_Xd.h"
+#define SDIS_XD_DIMENSION 2
+#include "sdis_heat_path_boundary_Xd_handle_external_net_flux.h"
+#define SDIS_XD_DIMENSION 3
+#include "sdis_heat_path_boundary_Xd_handle_external_net_flux.h"
diff --git a/src/sdis_heat_path_boundary_Xd_handle_external_net_flux.h b/src/sdis_heat_path_boundary_Xd_handle_external_net_flux.h
@@ -0,0 +1,128 @@
+/* Copyright (C) 2016-2023 |Méso|Star> (contact@meso-star.com)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "sdis_heat_path_boundary_c.h"
+#include "sdis_interface_c.h"
+#include "sdis_log.h"
+#include "sdis_scene_c.h"
+
+#include <rsys/cstr.h> /* res_to_cstr */
+
+#include "sdis_Xd_begin.h"
+
+/*******************************************************************************
+ * Generic helper functions
+ ******************************************************************************/
+static INLINE res_T
+XD(check_handle_external_net_flux_args)
+ (const struct sdis_device* dev,
+ const char* func_name,
+ const struct XD(handle_external_net_flux_args)* args)
+{
+ sdis_interface_sample_external_sources_T functor = NULL;
+ res_T res = RES_OK;
+
+ /* Handle bugs */
+ ASSERT(dev && func_name && args);
+ ASSERT(args->interf && args->frag);
+ ASSERT(!SXD_HIT_NONE(args->hit));
+
+ functor = interface_side_get_external_sources_sampling_functor
+ (args->interf, args->frag);
+
+ if(functor && args->picard_order != 0) {
+ res = RES_BAD_ARG;
+ log_err(dev,
+ "%s: Impossible to process external fluxes when Picard order is not "
+ "equal to 1; Picard order is currently set to %lu.\n",
+ func_name, (unsigned long)args->picard_order);
+ return res;
+ }
+
+ return RES_OK;
+}
+
+/*******************************************************************************
+ * Local functions
+ ******************************************************************************/
+res_T
+XD(handle_external_net_flux)
+ (const struct sdis_scene* scn,
+ struct ssp_rng* rng,
+ const struct XD(handle_external_net_flux_args)* args,
+ struct XD(temperature)* T)
+{
+ /* Sampling external sources */
+ struct sdis_external_sources_sample sample = SDIS_EXTERNAL_SOURCES_SAMPLE_NULL;
+ sdis_interface_sample_external_sources_T sample_sources = NULL;
+
+ /* Ray tracing */
+ struct hit_filter_data filter_data = HIT_FILTER_DATA_NULL;
+ struct sXd(hit) hit = SXD_HIT_NULL;
+ float ray_org[DIM] = {0};
+ float ray_dir[DIM] = {0};
+ float ray_range[2] = {0};
+
+ /* External flux */
+ double incident_direct_flux = 0; /* [W/m^2/sr] TODO ? */
+
+ /* Miscellaneous */
+ res_T res = RES_OK;
+ ASSERT(scn && args && T);
+
+ res = XD(check_handle_external_net_flux_args)(scn->dev, FUNC_NAME, args);
+ if(res != RES_OK) goto error;
+
+ /* Retrieve the functor to sample external sources */
+ sample_sources = interface_side_get_external_sources_sampling_functor
+ (args->interf, args->frag);
+
+ /* No external sources <=> no external fluxes. Nothing to do */
+ if(!sample_sources) goto exit;
+
+ /* Sample an external sources */
+ res = sample_sources(args->frag, rng, &sample, args->interf->data);
+ if(res != RES_OK) {
+ log_err(scn->dev,
+ "%s: error when sampling external sources -- %s.\n",
+ FUNC_NAME, res_to_cstr(res));
+ goto error;
+ }
+
+ /* Check whether the sampled external source is occluded or not */
+ filter_data.XD(hit) = *args->hit;
+ fX_set_dX(ray_org, args->frag->P);
+ fX_set_dX(ray_dir, sample.dir);
+ ray_range[0] = 0;
+ ray_range[1] = (float)sample.distance;
+ SXD(scene_view_trace_ray
+ (scn->sXd(view), ray_org, ray_dir, ray_range, &filter_data, &hit));
+
+ /* If the source is not occluded, compute its direct contribution */
+ if(!SXD_HIT_NONE(&hit)) {
+ const double cos_theta = dX(dot)(args->frag->Ng, sample.dir);
+ incident_direct_flux = fabs(cos_theta) * sample.radiance / sample.pdf;
+ }
+
+ /* TODO handle diffuse contribution */
+ /* TODO update the weight */
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
+#include "sdis_Xd_end.h"
diff --git a/src/sdis_heat_path_boundary_Xd_solid_fluid_picard1.h b/src/sdis_heat_path_boundary_Xd_solid_fluid_picard1.h
@@ -122,6 +122,10 @@ XD(solid_fluid_boundary_picard1_path)
/* Input argument used to handle the net flux */
struct handle_net_flux_args handle_net_flux_args = HANDLE_NET_FLUX_ARGS_NULL;
+ /* Input argument used to handle the external net flux */
+ struct XD(handle_external_net_flux_args) handle_external_net_flux_args =
+ XD(HANDLE_EXTERNAL_NET_FLUX_ARGS_NULL);
+
/* Input/output arguments of the function used to sample a reinjection */
struct XD(sample_reinjection_step_args) samp_reinject_step_args =
XD(SAMPLE_REINJECTION_STEP_ARGS_NULL);
@@ -238,6 +242,14 @@ XD(solid_fluid_boundary_picard1_path)
res = XD(handle_net_flux)(scn, &handle_net_flux_args, T);
if(res != RES_OK) goto error;
+ /* Handle the external net flux if any */
+ handle_external_net_flux_args.interf = interf;
+ handle_external_net_flux_args.frag = frag;
+ handle_external_net_flux_args.hit = &rwalk->hit;
+ handle_external_net_flux_args.picard_order = get_picard_order(ctx);
+ res = XD(handle_external_net_flux)(scn, rng, &handle_external_net_flux_args, T);
+ if(res != RES_OK) goto error;
+
/* Fetch the last registered heat path vertex */
if(ctx->heat_path) hvtx = *heat_path_get_last_vertex(ctx->heat_path);
diff --git a/src/sdis_heat_path_boundary_c.h b/src/sdis_heat_path_boundary_c.h
@@ -169,6 +169,52 @@ handle_net_flux_3d
struct temperature_3d* T);
/*******************************************************************************
+ * Handle external flux
+ ******************************************************************************/
+struct handle_external_net_flux_args_2d {
+ struct sdis_interface* interf;
+ const struct sdis_interface_fragment* frag;
+ const struct s2d_hit* hit;
+
+ struct green_path_handle* green_path; /* Store the propagator */
+ struct sdis_heat_path* heat_path; /* Save paths */
+
+ size_t picard_order;
+};
+
+struct handle_external_net_flux_args_3d {
+ struct sdis_interface* interf;
+ const struct sdis_interface_fragment* frag;
+ const struct s3d_hit* hit;
+
+ struct green_path_handle* green_path; /* Store the propagator */
+ struct sdis_heat_path* heat_path; /* Save paths */
+
+ size_t picard_order;
+};
+
+#define HANDLE_EXTERNAL_NET_FLUX_ARGS_NULL___2d {NULL,NULL,NULL,NULL,NULL,0}
+#define HANDLE_EXTERNAL_NET_FLUX_ARGS_NULL___3d {NULL,NULL,NULL,NULL,NULL,0}
+static const struct handle_external_net_flux_args_2d
+HANDLE_EXTERNAL_NET_FLUX_ARGS_NULL_2d = HANDLE_EXTERNAL_NET_FLUX_ARGS_NULL___2d;
+static const struct handle_external_net_flux_args_3d
+HANDLE_EXTERNAL_NET_FLUX_ARGS_NULL_3d = HANDLE_EXTERNAL_NET_FLUX_ARGS_NULL___3d;
+
+extern LOCAL_SYM res_T
+handle_external_net_flux_2d
+ (const struct sdis_scene* scn,
+ struct ssp_rng* rng,
+ const struct handle_external_net_flux_args_2d* args,
+ struct temperature_2d* T);
+
+extern LOCAL_SYM res_T
+handle_external_net_flux_3d
+ (const struct sdis_scene* scn,
+ struct ssp_rng* rng,
+ const struct handle_external_net_flux_args_3d* args,
+ struct temperature_3d* T);
+
+/*******************************************************************************
* Boundary sub-paths
******************************************************************************/
extern LOCAL_SYM res_T
diff --git a/src/sdis_interface_c.h b/src/sdis_interface_c.h
@@ -183,5 +183,20 @@ interface_side_get_reference_temperature
? shader->reference_temperature(frag, interf->data) : -1;
}
+static INLINE sdis_interface_sample_external_sources_T
+interface_side_get_external_sources_sampling_functor
+ (const struct sdis_interface* interf,
+ const struct sdis_interface_fragment* frag)
+{
+ const struct sdis_interface_side_shader* shader;
+ ASSERT(interf && frag);
+ switch(frag->side) {
+ case SDIS_FRONT: shader = &interf->shader.front; break;
+ case SDIS_BACK: shader = &interf->shader.back; break;
+ default: FATAL("Unreachable code\n"); break;
+ }
+ return shader->sample_external_sources;
+}
+
#endif /* SDIS_INTERFACE_C_H */
diff --git a/src/test_sdis_utils.h b/src/test_sdis_utils.h
@@ -192,7 +192,8 @@ static const struct sdis_fluid_shader DUMMY_FLUID_SHADER = {
dummy_interface_getter, /* Flux */ \
dummy_interface_getter, /* Emissivity */ \
dummy_interface_getter, /* Specular fraction */ \
- dummy_interface_getter /* Reference temperature */ \
+ dummy_interface_getter, /* Reference temperature */ \
+ NULL, /* External sources */ \
}
static const struct sdis_interface_shader DUMMY_INTERFACE_SHADER = {
dummy_interface_getter, /* Convection coef */