commit 2527eb556fd23d71ec627abe220904ed95f1a33c
parent f92eddb90dbe15af14c4b66237397f92470496f4
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Wed, 27 Feb 2019 10:46:11 +0100
Merge branch 'feature_paths' into develop
Diffstat:
35 files changed, 943 insertions(+), 125 deletions(-)
diff --git a/src/sdis.h b/src/sdis.h
@@ -51,12 +51,11 @@ struct senc2d_descriptor;
struct senc_descriptor;
/* Forward declaration of the Stardis opaque data types. These data types are
- * ref counted. Once created with the appropriated `sdis_<TYPE>_create'
- * function, the caller implicitly owns the created data, i.e. its reference
- * counter is set to 1. The sdis_<TYPE>_ref_<get|put> functions get or release
- * a reference on the data, i.e. they increment or decrement the reference
- * counter, respectively. When this counter reaches 0, the object is silently
- * destroyed and cannot be used anymore. */
+ * ref counted. Once created the caller implicitly owns the created data, i.e.
+ * its reference counter is set to 1. The sdis_<TYPE>_ref_<get|put> functions
+ * get or release a reference on the data, i.e. they increment or decrement the
+ * reference counter, respectively. When this counter reaches 0, the object is
+ * silently destroyed and cannot be used anymore. */
struct sdis_accum_buffer;
struct sdis_camera;
struct sdis_data;
@@ -67,6 +66,9 @@ struct sdis_interface;
struct sdis_medium;
struct sdis_scene;
+/* Forward declaration of non ref counted types */
+struct sdis_heat_path;
+
enum sdis_side {
SDIS_FRONT,
SDIS_BACK,
@@ -97,6 +99,19 @@ enum sdis_point_type {
SDIS_POINT_NONE = SDIS_POINT_TYPES_COUNT__
};
+enum sdis_heat_vertex_type {
+ SDIS_HEAT_VERTEX_CONDUCTION,
+ SDIS_HEAT_VERTEX_CONVECTION,
+ SDIS_HEAT_VERTEX_RADIATIVE
+};
+
+enum sdis_heat_path_flag {
+ SDIS_HEAT_PATH_SUCCEED = BIT(0),
+ SDIS_HEAT_PATH_FAILED = BIT(1),
+ SDIS_HEAT_PATH_ALL = SDIS_HEAT_PATH_SUCCEED | SDIS_HEAT_PATH_FAILED,
+ SDIS_HEAT_PATH_NONE = 0
+};
+
/* Random walk vertex, i.e. a spatiotemporal position at a given step of the
* random walk. */
struct sdis_rwalk_vertex {
@@ -245,6 +260,18 @@ typedef res_T
const size_t naccums[2], /* #accumulations in X and Y */
const struct sdis_accum* accums); /* List of row ordered accumulations */
+/* Vertex of heat path v*/
+struct sdis_heat_vertex {
+ double P[3];
+ double time;
+ double weight;
+ enum sdis_heat_vertex_type type;
+};
+#define SDIS_HEAT_VERTEX_NULL__ {{0,0,0}, 0, 0, SDIS_HEAT_VERTEX_CONDUCTION}
+static const struct sdis_heat_vertex SDIS_HEAT_VERTEX_NULL =
+ SDIS_HEAT_VERTEX_NULL__;
+
+/* Path used to estimate the green function */
struct sdis_green_path {
/* Internal data. Should not be accessed */
void* green__;
@@ -270,7 +297,7 @@ struct sdis_point {
#define SDIS_POINT_NULL__ { {{NULL, SDIS_RWALK_VERTEX_NULL__}}, SDIS_POINT_NONE}
static const struct sdis_point SDIS_POINT_NULL = SDIS_POINT_NULL__;
-/* Functor used to process the paths registered agains the green function */
+/* Functor used to process the paths registered against the green function */
typedef res_T
(*sdis_process_green_path_T)
(struct sdis_green_path* path,
@@ -293,6 +320,18 @@ typedef res_T
const double flux_term,
void* context);
+/* Functor used to process a heat path registered against the estimator */
+typedef res_T
+(*sdis_process_heat_path_T)
+ (const struct sdis_heat_path* path,
+ void* context);
+
+/* Functor used to process the vertices of a heat path */
+typedef res_T
+(*sdis_process_heat_vertex_T)
+ (const struct sdis_heat_vertex* vertex,
+ void* context);
+
BEGIN_DECLS
/*******************************************************************************
@@ -688,6 +727,23 @@ sdis_estimator_get_total_flux
(const struct sdis_estimator* estimator,
struct sdis_mc* flux);
+SDIS_API res_T
+sdis_estimator_get_paths_count
+ (const struct sdis_estimator* estimator,
+ size_t* npaths);
+
+SDIS_API res_T
+sdis_estimator_get_path
+ (const struct sdis_estimator* estimator,
+ const size_t ipath,
+ const struct sdis_heat_path** path);
+
+SDIS_API res_T
+sdis_estimator_for_each_path
+ (const struct sdis_estimator* estimator,
+ sdis_process_heat_path_T func,
+ void* context);
+
/*******************************************************************************
* The green function saves the estimation of the propagator
******************************************************************************/
@@ -754,6 +810,31 @@ sdis_green_path_for_each_flux_term
void* context);
/*******************************************************************************
+ * Heat path API
+ ******************************************************************************/
+SDIS_API res_T
+sdis_heat_path_get_vertices_count
+ (const struct sdis_heat_path* path,
+ size_t* nvertices);
+
+SDIS_API res_T
+sdis_heat_path_get_status
+ (const struct sdis_heat_path* path,
+ enum sdis_heat_path_flag* status);
+
+SDIS_API res_T
+sdis_heat_path_get_vertex
+ (const struct sdis_heat_path* path,
+ const size_t ivertex,
+ struct sdis_heat_vertex* vertex);
+
+SDIS_API res_T
+sdis_heat_path_for_each_vertex
+ (const struct sdis_heat_path* path,
+ sdis_process_heat_vertex_T func,
+ void* context);
+
+/*******************************************************************************
* Solvers
******************************************************************************/
SDIS_API res_T
@@ -765,6 +846,7 @@ sdis_solve_probe
const double fp_to_meter, /* Scale from floating point units to meters */
const double ambient_radiative_temperature, /* In Kelvin */
const double reference_temperature, /* In Kelvin */
+ const int register_path, /* Combination of enum sdis_heat_path_flag */
struct sdis_estimator** estimator);
SDIS_API res_T
diff --git a/src/sdis_Xd_begin.h b/src/sdis_Xd_begin.h
@@ -20,13 +20,15 @@
/* Forward declaration */
struct green_path_handle;
+struct sdis_heat_path;
struct rwalk_context {
struct green_path_handle* green_path;
+ struct sdis_heat_path* heat_path;
double Tarad; /* Ambient radiative temperature */
double Tref3; /* Reference temperature ^ 3 */
};
-#define RWALK_CONTEXT_NULL__ {NULL, 0, 0}
+#define RWALK_CONTEXT_NULL__ {NULL, NULL, 0, 0}
static const struct rwalk_context RWALK_CONTEXT_NULL = RWALK_CONTEXT_NULL__;
#endif /* SDIS_XD_BEGIN_H */
diff --git a/src/sdis_estimator.c b/src/sdis_estimator.c
@@ -17,6 +17,8 @@
#include "sdis_device_c.h"
#include "sdis_estimator_c.h"
+#include <rsys/mutex.h>
+
/*******************************************************************************
* Helper functions
******************************************************************************/
@@ -28,6 +30,8 @@ estimator_release(ref_T* ref)
ASSERT(ref);
estimator = CONTAINER_OF(ref, struct sdis_estimator, ref);
dev = estimator->dev;
+ darray_heat_path_release(&estimator->paths);
+ if(estimator->mutex) mutex_destroy(estimator->mutex);
MEM_RM(dev->allocator, estimator);
SDIS(device_ref_put(dev));
}
@@ -120,6 +124,56 @@ sdis_estimator_get_total_flux
return RES_OK;
}
+res_T
+sdis_estimator_get_paths_count
+ (const struct sdis_estimator* estimator, size_t* npaths)
+{
+ if(!estimator || !npaths) return RES_BAD_ARG;
+ *npaths = darray_heat_path_size_get(&estimator->paths);
+ return RES_OK;
+}
+
+SDIS_API res_T
+sdis_estimator_get_path
+ (const struct sdis_estimator* estimator,
+ const size_t ipath,
+ const struct sdis_heat_path** path)
+{
+ if(!estimator || !path
+ || ipath >= darray_heat_path_size_get(&estimator->paths))
+ return RES_BAD_ARG;
+ *path = darray_heat_path_cdata_get(&estimator->paths) + ipath;
+ return RES_OK;
+}
+
+res_T
+sdis_estimator_for_each_path
+ (const struct sdis_estimator* estimator,
+ sdis_process_heat_path_T func,
+ void* context)
+{
+ const struct sdis_heat_path* paths = NULL;
+ size_t i, n;
+ res_T res = RES_OK;
+
+ if(!estimator || !func) {
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ SDIS(estimator_get_paths_count(estimator, &n));
+ paths = darray_heat_path_cdata_get(&estimator->paths);
+ FOR_EACH(i, 0, n) {
+ res = func(paths+i, context);
+ if(res != RES_OK) goto error;
+ }
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
/*******************************************************************************
* Local functions
******************************************************************************/
@@ -127,8 +181,6 @@ res_T
estimator_create
(struct sdis_device* dev,
const enum sdis_estimator_type type,
- const size_t nrealisations,
- const size_t nsuccesses,
struct sdis_estimator** out_estimator)
{
struct sdis_estimator* estimator = NULL;
@@ -136,9 +188,6 @@ estimator_create
if(!dev
|| (unsigned)type >= SDIS_ESTIMATOR_TYPES_COUNT__
- || !nrealisations
- || !nsuccesses
- || nsuccesses > nrealisations
|| !out_estimator) {
res = RES_BAD_ARG;
goto error;
@@ -151,10 +200,18 @@ estimator_create
}
ref_init(&estimator->ref);
SDIS(device_ref_get(dev));
- estimator->nrealisations = nsuccesses;
- estimator->nfailures = nrealisations - nsuccesses;
+ estimator->nrealisations = 0;
+ estimator->nfailures = 0;
estimator->dev = dev;
estimator->type = type;
+ darray_heat_path_init(dev->allocator, &estimator->paths);
+
+ estimator->mutex = mutex_create();
+ if(!estimator->mutex) {
+ res = RES_MEM_ERR;
+ goto error;
+ }
+
exit:
if(out_estimator) *out_estimator = estimator;
return res;
@@ -166,3 +223,30 @@ error:
goto exit;
}
+res_T
+estimator_add_and_release_heat_path
+ (struct sdis_estimator* estimator, struct sdis_heat_path* path)
+{
+ struct sdis_heat_path* dst = NULL;
+ size_t i;
+ res_T res = RES_OK;
+ ASSERT(estimator && path);
+
+ mutex_lock(estimator->mutex);
+
+ i = darray_heat_path_size_get(&estimator->paths);
+
+ res = darray_heat_path_resize(&estimator->paths, i+1);
+ if(res != RES_OK) goto error;
+
+ dst = darray_heat_path_data_get(&estimator->paths) + i;
+ res = heat_path_copy_and_release(dst, path);
+ if(res != RES_OK) goto error;
+
+exit:
+ mutex_unlock(estimator->mutex);
+ return res;
+error:
+ goto exit;
+}
+
diff --git a/src/sdis_estimator_c.h b/src/sdis_estimator_c.h
@@ -16,6 +16,8 @@
#ifndef SDIS_ESTIMATOR_C_H
#define SDIS_ESTIMATOR_C_H
+#include "sdis_heat_path.h"
+
#include <rsys/math.h>
#include <rsys/ref_count.h>
@@ -37,11 +39,16 @@ struct sdis_estimator {
size_t nrealisations;
size_t nfailures;
+ struct mutex* mutex;
+ struct darray_heat_path paths; /* Tracked paths */
+
enum sdis_estimator_type type;
ref_T ref;
struct sdis_device* dev;
};
+struct sdis_estimator_handle;
+
/*******************************************************************************
* Estimator local API
******************************************************************************/
@@ -49,10 +56,26 @@ extern LOCAL_SYM res_T
estimator_create
(struct sdis_device* dev,
const enum sdis_estimator_type type,
- const size_t nrealisations,
- const size_t nsuccesses,
struct sdis_estimator** estimator);
+/* Thread safe */
+extern LOCAL_SYM res_T
+estimator_add_and_release_heat_path
+ (struct sdis_estimator* estimator,
+ struct sdis_heat_path* path);
+
+/* Must be invoked before any others "estimator_setup" functions */
+static INLINE void
+estimator_setup_realisations_count
+ (struct sdis_estimator* estimator,
+ const size_t nrealisations,
+ const size_t nsuccesses)
+{
+ ASSERT(estimator && nrealisations && nsuccesses && nsuccesses<=nrealisations);
+ estimator->nrealisations = nsuccesses;
+ estimator->nfailures = nrealisations - nsuccesses;
+}
+
static INLINE void
estimator_setup_temperature
(struct sdis_estimator* estim,
diff --git a/src/sdis_green.c b/src/sdis_green.c
@@ -191,6 +191,9 @@ struct sdis_green_function {
size_t npaths_valid;
size_t npaths_invalid;
+ struct ssp_rng_type rng_type;
+ FILE* rng_state;
+
ref_T ref;
struct sdis_device* dev;
};
@@ -408,6 +411,7 @@ green_function_release(ref_T* ref)
htable_medium_release(&green->media);
htable_interf_release(&green->interfaces);
darray_green_path_release(&green->paths);
+ if(green->rng_state) fclose(green->rng_state);
MEM_RM(dev->allocator, green);
SDIS(device_ref_put(dev));
}
@@ -452,15 +456,23 @@ sdis_green_function_solve
goto error;
}
- /* FIXME do not use a new RNG. Save the RNG state into the green function
- * after its estimation and initialize the following rng with this state */
- res = ssp_rng_create(green->dev->allocator, &ssp_rng_mt19937_64, &rng);
+ res = ssp_rng_create(green->dev->allocator, &green->rng_type, &rng);
+ if(res != RES_OK) goto error;
+
+ /* Avoid correlation by defining the RNG state from the final state of the
+ * RNG used to estimate the green function */
+ rewind(green->rng_state);
+ res = ssp_rng_read(rng, green->rng_state);
if(res != RES_OK) goto error;
npaths = darray_green_path_size_get(&green->paths);
+ /* Create the estimator */
+ res = estimator_create(green->dev, SDIS_ESTIMATOR_TEMPERATURE, &estimator);
+ if(res != RES_OK) goto error;
+
/* Solve the green function */
- FOR_EACH(ipath, 0, npaths) { /* TODO add multi-threading (?) */
+ FOR_EACH(ipath, 0, npaths) {
const double time = sample_time(rng, time_range);
double w;
@@ -473,12 +485,8 @@ sdis_green_function_solve
++N;
}
- /* Create the estimator */
- res = estimator_create
- (green->dev, SDIS_ESTIMATOR_TEMPERATURE, npaths, N, &estimator);
- if(res != RES_OK) goto error;
-
/* Setup the estimated temperature */
+ estimator_setup_realisations_count(estimator, npaths, N);
estimator_setup_temperature(estimator, accum, accum2);
exit:
@@ -685,6 +693,12 @@ green_function_create
green->npaths_valid = SIZE_MAX;
green->npaths_invalid = SIZE_MAX;
+ green->rng_state = tmpfile();
+ if(!green->rng_state) {
+ res = RES_IO_ERR;
+ goto error;
+ }
+
exit:
*out_green = green;
return res;
@@ -759,12 +773,24 @@ error:
}
res_T
-green_function_finalize(struct sdis_green_function* green)
+green_function_finalize
+ (struct sdis_green_function* green,
+ struct ssp_rng_proxy* proxy)
{
size_t i, n;
+ res_T res = RES_OK;
- if(!green) return RES_BAD_ARG;
+ if(!green || !proxy) {
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ /* Save the RNG state */
+ SSP(rng_proxy_get_type(proxy, &green->rng_type));
+ res = ssp_rng_proxy_write(proxy, green->rng_state);
+ if(res != RES_OK) goto error;
+ /* Compute the number of valid/invalid green paths */
green->npaths_valid = 0;
n = darray_green_path_size_get(&green->paths);
FOR_EACH(i, 0, n) {
@@ -772,7 +798,11 @@ green_function_finalize(struct sdis_green_function* green)
green->npaths_valid += path->limit_type != SDIS_POINT_NONE;
}
green->npaths_invalid = n - green->npaths_valid;
- return RES_OK;
+
+exit:
+ return res;
+error:
+ goto exit;
}
res_T
diff --git a/src/sdis_green.h b/src/sdis_green.h
@@ -20,6 +20,7 @@
/* Forward declaration */
struct sdis_green_function;
+struct ssp_rng_proxy;
struct green_path;
struct green_path_handle {
@@ -41,10 +42,12 @@ green_function_merge_and_clear
(struct sdis_green_function* dst,
struct sdis_green_function* src);
-/* Finalize the green function state (e.g.: computes the #paths & #failures) */
+/* Finalize the green function state (e.g.: computes the #paths & #failures,
+ * save the rng state, etc.) */
extern LOCAL_SYM res_T
green_function_finalize
- (struct sdis_green_function* green);
+ (struct sdis_green_function* green,
+ struct ssp_rng_proxy* rng_proxy); /* Proxy RNG used to estimate the function */
extern LOCAL_SYM res_T
green_function_create_path
diff --git a/src/sdis_heat_path.c b/src/sdis_heat_path.c
@@ -39,3 +39,67 @@
#define SDIS_XD_DIMENSION 3
#include "sdis_heat_path_boundary_Xd.h"
+/*******************************************************************************
+ * Exported functions
+ ******************************************************************************/
+res_T
+sdis_heat_path_get_vertices_count
+ (const struct sdis_heat_path* path, size_t* nvertices)
+{
+ if(!path || !nvertices) return RES_BAD_ARG;
+ *nvertices = darray_heat_vertex_size_get(&path->vertices);
+ return RES_OK;
+}
+
+res_T
+sdis_heat_path_get_status
+ (const struct sdis_heat_path* path, enum sdis_heat_path_flag* status)
+{
+ if(!path || !status) return RES_BAD_ARG;
+ *status = path->status;
+ return RES_OK;
+}
+
+res_T
+sdis_heat_path_get_vertex
+ (const struct sdis_heat_path* path,
+ const size_t ivertex,
+ struct sdis_heat_vertex* vertex)
+{
+ if(!path || !vertex
+ || ivertex >= darray_heat_vertex_size_get(&path->vertices)) {
+ return RES_BAD_ARG;
+ }
+
+ *vertex = darray_heat_vertex_cdata_get(&path->vertices)[ivertex];
+ return RES_OK;
+}
+
+res_T
+sdis_heat_path_for_each_vertex
+ (const struct sdis_heat_path* path,
+ sdis_process_heat_vertex_T func,
+ void* context)
+{
+ const struct sdis_heat_vertex* vertices;
+ size_t i, n;
+ res_T res = RES_OK;
+
+ if(!path || !func) {
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+ SDIS(heat_path_get_vertices_count(path, &n));
+ vertices = darray_heat_vertex_cdata_get(&path->vertices);
+ FOR_EACH(i, 0, n) {
+ res = func(vertices+i, context);
+ if(res != RES_OK) goto error;
+ }
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
diff --git a/src/sdis_heat_path.h b/src/sdis_heat_path.h
@@ -16,8 +16,12 @@
#ifndef SDIS_HEAT_PATH_H
#define SDIS_HEAT_PATH_H
+#include "sdis.h"
+
+#include <rsys/dynamic_array.h>
#include <rsys/rsys.h>
+/* Forward declarations */
struct rwalk_2d;
struct rwalk_3d;
struct rwalk_context;
@@ -26,6 +30,84 @@ struct ssp_rng;
struct temperature_2d;
struct temperature_3d;
+/* Generate the dynamic array of heat vertices */
+#define DARRAY_NAME heat_vertex
+#define DARRAY_DATA struct sdis_heat_vertex
+#include <rsys/dynamic_array.h>
+
+/*******************************************************************************
+ * Heat path data structure
+ ******************************************************************************/
+struct sdis_heat_path {
+ struct darray_heat_vertex vertices;
+ enum sdis_heat_path_flag status;
+};
+
+static INLINE void
+heat_path_init(struct mem_allocator* allocator, struct sdis_heat_path* path)
+{
+ ASSERT(path);
+ path->status = SDIS_HEAT_PATH_NONE;
+ darray_heat_vertex_init(allocator, &path->vertices);
+}
+
+static INLINE void
+heat_path_release(struct sdis_heat_path* path)
+{
+ ASSERT(path);
+ darray_heat_vertex_release(&path->vertices);
+}
+
+static INLINE res_T
+heat_path_copy(struct sdis_heat_path* dst, const struct sdis_heat_path* src)
+{
+ ASSERT(dst && src);
+ dst->status = src->status;
+ return darray_heat_vertex_copy(&dst->vertices, &src->vertices);
+}
+
+static INLINE res_T
+heat_path_copy_and_release(struct sdis_heat_path* dst, struct sdis_heat_path* src)
+{
+ ASSERT(dst && src);
+ dst->status = src->status;
+ return darray_heat_vertex_copy_and_release(&dst->vertices, &src->vertices);
+}
+
+static INLINE res_T
+heat_path_copy_and_clear(struct sdis_heat_path* dst, struct sdis_heat_path* src)
+{
+ ASSERT(dst && src);
+ dst->status = src->status;
+ return darray_heat_vertex_copy_and_clear(&dst->vertices, &src->vertices);
+}
+
+static INLINE res_T
+heat_path_add_vertex(struct sdis_heat_path* path, const struct sdis_heat_vertex* vtx)
+{
+ ASSERT(path && vtx);
+ return darray_heat_vertex_push_back(&path->vertices, vtx);
+}
+
+static INLINE struct sdis_heat_vertex*
+heat_path_get_last_vertex(struct sdis_heat_path* path)
+{
+ size_t sz;
+ ASSERT(path);
+ sz = darray_heat_vertex_size_get(&path->vertices);
+ ASSERT(sz);
+ return darray_heat_vertex_data_get(&path->vertices) + (sz-1);
+}
+
+/* Generate the dynamic array of heat paths */
+#define DARRAY_NAME heat_path
+#define DARRAY_DATA struct sdis_heat_path
+#define DARRAY_FUNCTOR_INIT heat_path_init
+#define DARRAY_FUNCTOR_RELEASE heat_path_release
+#define DARRAY_FUNCTOR_COPY heat_path_copy
+#define DARRAY_FUNCTOR_COPY_AND_RELEASE heat_path_copy_and_release
+#include <rsys/dynamic_array.h>
+
/*******************************************************************************
* Trace or pursue a radiative path
******************************************************************************/
diff --git a/src/sdis_heat_path_boundary_Xd.h b/src/sdis_heat_path_boundary_Xd.h
@@ -252,6 +252,11 @@ XD(solid_solid_boundary_path)
rwalk->hit_side = SDIS_SIDE_NULL__;
}
+ /* Register the new vertex against the heat path */
+ res = register_heat_vertex
+ (ctx->heat_path, &rwalk->vtx, T->value, SDIS_HEAT_VERTEX_CONDUCTION);
+ if(res != RES_OK) goto error;
+
exit:
return res;
error:
@@ -411,6 +416,11 @@ XD(solid_fluid_boundary_path)
rwalk->hit = SXD_HIT_NULL;
rwalk->hit_side = SDIS_SIDE_NULL__;
}
+
+ /* Register the new vertex against the heat path */
+ res = register_heat_vertex
+ (ctx->heat_path, &rwalk->vtx, T->value, SDIS_HEAT_VERTEX_CONDUCTION);
+ if(res != RES_OK) goto error;
}
exit:
@@ -541,8 +551,14 @@ XD(solid_boundary_with_flux_path)
rwalk->mdm = mdm;
rwalk->hit = SXD_HIT_NULL;
rwalk->hit_side = SDIS_SIDE_NULL__;
+
}
+ /* Register the new vertex against the heat path */
+ res = register_heat_vertex
+ (ctx->heat_path, &rwalk->vtx, T->value, SDIS_HEAT_VERTEX_CONDUCTION);
+ if(res != RES_OK) goto error;
+
exit:
return res;
error:
@@ -590,6 +606,9 @@ XD(boundary_path)
(ctx->green_path, interf, &frag);
if(res != RES_OK) goto error;
}
+ if(ctx->heat_path) {
+ heat_path_get_last_vertex(ctx->heat_path)->weight = T->value;
+ }
goto exit;
}
diff --git a/src/sdis_heat_path_conductive_Xd.h b/src/sdis_heat_path_conductive_Xd.h
@@ -80,9 +80,15 @@ XD(conductive_path)
T->done = 1;
if(ctx->green_path) {
- res = green_path_set_limit_vertex(ctx->green_path, rwalk->mdm, &rwalk->vtx);
+ res = green_path_set_limit_vertex
+ (ctx->green_path, rwalk->mdm, &rwalk->vtx);
if(res != RES_OK) goto error;
}
+
+ if(ctx->heat_path) {
+ heat_path_get_last_vertex(ctx->heat_path)->weight = T->value;
+ }
+
break;
}
@@ -176,7 +182,8 @@ XD(conductive_path)
}
}
- /* Register the power term for the green function */
+ /* Register the power term for the green function. Delay its registration
+ * until the end of the conductive path, i.e. the path is valid */
if(ctx->green_path && power != SDIS_VOLUMIC_POWER_NONE) {
green_power_factor += power_factor;
}
@@ -194,6 +201,13 @@ XD(conductive_path)
if(tmp >= 0) {
T->value += tmp;
T->done = 1;
+
+ if(ctx->heat_path) {
+ struct sdis_heat_vertex* vtx;
+ vtx = heat_path_get_last_vertex(ctx->heat_path);
+ vtx->time = rwalk->vtx.time;
+ vtx->weight = T->value;
+ }
break;
}
/* The initial condition should have been reached */
@@ -218,6 +232,11 @@ XD(conductive_path)
/* Update the random walk position */
XD(move_pos)(rwalk->vtx.P, dir0, delta);
+ /* Register the new vertex against the heat path */
+ res = register_heat_vertex
+ (ctx->heat_path, &rwalk->vtx, T->value, SDIS_HEAT_VERTEX_CONDUCTION);
+ if(res != RES_OK) goto error;
+
/* Fetch the current medium */
if(SXD_HIT_NONE(&rwalk->hit)) {
CHK(scene_get_medium(scn, rwalk->vtx.P, &info, &mdm) == RES_OK);
diff --git a/src/sdis_heat_path_convective_Xd.h b/src/sdis_heat_path_convective_Xd.h
@@ -23,6 +23,47 @@
#include "sdis_Xd_begin.h"
+/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+static res_T
+XD(register_heat_vertex_in_fluid)
+ (struct sdis_scene* scn,
+ const struct rwalk_context* ctx,
+ struct XD(rwalk)* rwalk,
+ const double weight)
+{
+ struct sdis_rwalk_vertex vtx = SDIS_RWALK_VERTEX_NULL;
+ const float empirical_dst = 0.1f;
+ const float range[2] = {0, FLT_MAX};
+ float org[DIM];
+ float dir[DIM];
+ float pos[DIM];
+ float dst;
+ struct sXd(hit) hit;
+
+ if(!ctx->heat_path) return RES_OK;
+
+ ASSERT(!SXD_HIT_NONE(&rwalk->hit));
+
+ fX_set_dX(org, rwalk->vtx.P);
+ fX(set)(dir, rwalk->hit.normal);
+ if(rwalk->hit_side == SDIS_BACK) fX(minus)(dir, dir);
+
+ SXD(scene_view_trace_ray(scn->sXd(view), org, dir, range, &rwalk->hit, &hit));
+ dst = SXD_HIT_NONE(&hit) ? empirical_dst : hit.distance * 0.5f;
+
+ vtx = rwalk->vtx;
+ fX(add)(pos, org, fX(mulf)(dir, dir, dst));
+ dX_set_fX(vtx.P, pos);
+
+ return register_heat_vertex
+ (ctx->heat_path, &vtx, weight, SDIS_HEAT_VERTEX_CONVECTION);
+}
+
+/*******************************************************************************
+ * Local functions
+ ******************************************************************************/
res_T
XD(convective_path)
(struct sdis_scene* scn,
@@ -33,7 +74,6 @@ XD(convective_path)
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];
@@ -62,6 +102,10 @@ XD(convective_path)
res = green_path_set_limit_vertex(ctx->green_path, rwalk->mdm, &rwalk->vtx);
if(res != RES_OK) goto error;
}
+
+ res = XD(register_heat_vertex_in_fluid)(scn, ctx, rwalk, T->value);
+ if(res != RES_OK) goto error;
+
goto exit;
}
@@ -147,30 +191,10 @@ XD(convective_path)
goto error;
}
- /* A trick to force first r test result.
- * TODO fix this workaround that seems useless */
- r = 1;
-
/* Sample time until init condition is reached or a true convection occurs. */
for(;;) {
+ struct sdis_interface_fragment frag;
struct sXd(primitive) prim;
- /* Setup the fragment of the interface. */
- XD(setup_interface_fragment)(&frag, &rwalk->vtx, &rwalk->hit, rwalk->hit_side);
-
- /* Fetch hc. */
- hc = interface_get_convection_coef(interf, &frag);
- if(hc > enc->hc_upper_bound) {
- log_err(scn->dev,
- "%s: hc (%g) exceeds its provided upper bound (%g) at %g %g %g.\n",
- FUNC_NAME, hc, enc->hc_upper_bound, SPLIT3(rwalk->vtx.P));
- res = RES_BAD_OP;
- goto error;
- }
-
- if(r < hc / enc->hc_upper_bound) {
- /* True convection. Always true if hc == bound. */
- break;
- }
/* Fetch other physical properties. */
cp = fluid_get_calorific_capacity(rwalk->mdm, &rwalk->vtx);
@@ -183,12 +207,20 @@ XD(convective_path)
tau = ssp_ran_exp(rng, mu);
t0 = ctx->green_path ? -INF : fluid_get_t0(rwalk->mdm);
rwalk->vtx.time = MMAX(rwalk->vtx.time - tau, t0);
+
+ /* Register the new vertex against the heat path */
+ res = XD(register_heat_vertex_in_fluid)(scn, ctx, rwalk, T->value);
+ if(res != RES_OK) goto error;
+
if(rwalk->vtx.time == t0) {
/* Check the initial condition. */
tmp = fluid_get_temperature(rwalk->mdm, &rwalk->vtx);
if(tmp >= 0) {
T->value += tmp;
T->done = 1;
+ if(ctx->heat_path) { /* Update the weight of the last heat vertex */
+ heat_path_get_last_vertex(ctx->heat_path)->weight = T->value;
+ }
goto exit;
}
/* The initial condition should have been reached. */
@@ -238,8 +270,29 @@ XD(convective_path)
FATAL("Unexpected fluid interface.\n");
}
- /* Renew r for next loop. */
+ /* Register the new vertex against the heat path */
+ res = register_heat_vertex
+ (ctx->heat_path, &rwalk->vtx, T->value, SDIS_HEAT_VERTEX_CONVECTION);
+ if(res != RES_OK) goto error;
+
+ /* Setup the fragment of the sampled position into the enclosure. */
+ XD(setup_interface_fragment)(&frag, &rwalk->vtx, &rwalk->hit, rwalk->hit_side);
+
+ /* Fetch the convection coefficient of the sampled position */
+ hc = interface_get_convection_coef(interf, &frag);
+ if(hc > enc->hc_upper_bound) {
+ log_err(scn->dev,
+ "%s: hc (%g) exceeds its provided upper bound (%g) at %g %g %g.\n",
+ FUNC_NAME, hc, enc->hc_upper_bound, SPLIT3(rwalk->vtx.P));
+ res = RES_BAD_OP;
+ goto error;
+ }
+
r = ssp_rng_canonical_float(rng);
+ if(r < hc / enc->hc_upper_bound) {
+ /* True convection. Always true if hc == bound. */
+ break;
+ }
}
rwalk->hit.distance = 0;
diff --git a/src/sdis_heat_path_radiative_Xd.h b/src/sdis_heat_path_radiative_Xd.h
@@ -39,7 +39,7 @@ XD(trace_radiative_path)
struct XD(temperature)* T)
{
/* The radiative random walk is always performed in 3D. In 2D, the geometry
- * are assumed to be extruded to the infinty along the Z dimension. */
+ * are assumed to be extruded to the infinity along the Z dimension. */
float N[3] = {0, 0, 0};
float dir[3] = {0, 0, 0};
res_T res = RES_OK;
@@ -83,6 +83,17 @@ XD(trace_radiative_path)
res = green_path_set_limit_vertex(ctx->green_path, rwalk->mdm, &vtx);
if(res != RES_OK) goto error;
}
+ if(ctx->heat_path) {
+ const float empirical_dst = 0.1f;
+ struct sdis_rwalk_vertex vtx;
+ vtx = rwalk->vtx;
+ vtx.P[0] += dir[0] * empirical_dst;
+ vtx.P[1] += dir[1] * empirical_dst;
+ vtx.P[2] += dir[2] * empirical_dst;
+ res = register_heat_vertex
+ (ctx->heat_path, &vtx, T->value, SDIS_HEAT_VERTEX_RADIATIVE);
+ if(res != RES_OK) goto error;
+ }
break;
} else {
log_err(scn->dev,
@@ -107,6 +118,11 @@ XD(trace_radiative_path)
/* Move the random walk to the hit position */
XD(move_pos)(rwalk->vtx.P, dir, rwalk->hit.distance);
+ /* Register the random walk vertex against the heat path */
+ res = register_heat_vertex
+ (ctx->heat_path, &rwalk->vtx, T->value, SDIS_HEAT_VERTEX_RADIATIVE);
+ if(res != RES_OK) goto error;
+
/* Fetch the new interface and setup the hit fragment */
interf = scene_get_interface(scn, rwalk->hit.prim.prim_id);
XD(setup_interface_fragment)(&frag, &rwalk->vtx, &rwalk->hit, rwalk->hit_side);
@@ -179,7 +195,7 @@ XD(radiative_path)
struct XD(temperature)* T)
{
/* The radiative random walk is always performed in 3D. In 2D, the geometry
- * are assumed to be extruded to the infinty along the Z dimension. */
+ * are assumed to be extruded to the infinity along the Z dimension. */
float N[3] = {0, 0, 0};
float dir[3] = {0, 0, 0};
res_T res = RES_OK;
diff --git a/src/sdis_misc.h b/src/sdis_misc.h
@@ -87,4 +87,25 @@ sample_time(struct ssp_rng* rng, const double time_range[2])
return ssp_rng_uniform_double(rng, time_range[0], time_range[1]);
}
+static INLINE res_T
+register_heat_vertex
+ (struct sdis_heat_path* path,
+ const struct sdis_rwalk_vertex* vtx,
+ const double weight,
+ const enum sdis_heat_vertex_type type)
+{
+ struct sdis_heat_vertex heat_vtx = SDIS_HEAT_VERTEX_NULL;
+ ASSERT(vtx);
+
+ if(!path) return RES_OK;
+
+ heat_vtx.P[0] = vtx->P[0];
+ heat_vtx.P[1] = vtx->P[1];
+ heat_vtx.P[2] = vtx->P[2];
+ heat_vtx.time = vtx->time;
+ heat_vtx.weight = weight;
+ heat_vtx.type = type;
+ return heat_path_add_vertex(path, &heat_vtx);
+}
+
#endif /* SDIS_MISC_H */
diff --git a/src/sdis_realisation.h b/src/sdis_realisation.h
@@ -46,6 +46,7 @@ probe_realisation_2d
const double ambient_radiative_temperature,
const double reference_temperature,
struct green_path_handle* green_path,
+ struct sdis_heat_path* heat_path,
double* weight);
extern LOCAL_SYM res_T
@@ -59,6 +60,7 @@ probe_realisation_3d
const double ambient_radiative_temperature,
const double reference_temperature,
struct green_path_handle* green_path,
+ struct sdis_heat_path* heat_path,
double* weight);
/*******************************************************************************
diff --git a/src/sdis_realisation_Xd.h b/src/sdis_realisation_Xd.h
@@ -97,11 +97,13 @@ XD(probe_realisation)
const double ambient_radiative_temperature,
const double reference_temperature,
struct green_path_handle* green_path,
+ struct sdis_heat_path* heat_path,
double* weight)
{
struct rwalk_context ctx = RWALK_CONTEXT_NULL;
struct XD(rwalk) rwalk = XD(RWALK_NULL);
struct XD(temperature) T = XD(TEMPERATURE_NULL);
+ enum sdis_heat_vertex_type type;
double t0;
double (*get_initial_temperature)
(const struct sdis_medium* mdm,
@@ -125,6 +127,14 @@ XD(probe_realisation)
dX(set)(rwalk.vtx.P, position);
rwalk.vtx.time = time;
+
+ /* Register the starting position against the heat path */
+ type = medium->type == SDIS_SOLID
+ ? SDIS_HEAT_VERTEX_CONDUCTION
+ : SDIS_HEAT_VERTEX_CONVECTION;
+ res = register_heat_vertex(heat_path, &rwalk.vtx, 0, type);
+ if(res != RES_OK) goto error;
+
/* No initial condition with green */
if(!green_path && t0 >= rwalk.vtx.time) {
double tmp;
@@ -133,20 +143,22 @@ XD(probe_realisation)
tmp = get_initial_temperature(medium, &rwalk.vtx);
if(tmp >= 0) {
*weight = tmp;
- return RES_OK;
+ goto exit;
}
/* The initial condition should have been reached */
log_err(scn->dev,
"%s: undefined initial condition. "
"The time is %f but the temperature remains unknown.\n",
FUNC_NAME, t0);
- return RES_BAD_OP;
+ res = RES_BAD_OP;
+ goto error;
}
rwalk.hit = SXD_HIT_NULL;
rwalk.mdm = medium;
ctx.green_path = green_path;
+ ctx.heat_path = heat_path;
ctx.Tarad = ambient_radiative_temperature;
ctx.Tref3 =
reference_temperature
@@ -154,10 +166,14 @@ XD(probe_realisation)
* reference_temperature;
res = XD(compute_temperature)(scn, fp_to_meter, &ctx, &rwalk, rng, &T);
- if(res != RES_OK) return res;
+ if(res != RES_OK) goto error;
*weight = T.value;
- return RES_OK;
+
+exit:
+ return res;
+error:
+ goto exit;
}
res_T
diff --git a/src/sdis_solve.c b/src/sdis_solve.c
@@ -167,15 +167,16 @@ sdis_solve_probe
const double fp_to_meter,/* Scale factor from floating point unit to meter */
const double Tarad, /* Ambient radiative temperature */
const double Tref, /* Reference temperature */
+ const int register_paths, /* Combination of enum sdis_heat_path_flag */
struct sdis_estimator** out_estimator)
{
if(!scn) return RES_BAD_ARG;
if(scene_is_2d(scn)) {
return solve_probe_2d(scn, nrealisations, position, time_range,
- fp_to_meter, Tarad, Tref, NULL, out_estimator);
+ fp_to_meter, Tarad, Tref, register_paths, NULL, out_estimator);
} else {
return solve_probe_3d(scn, nrealisations, position, time_range,
- fp_to_meter, Tarad, Tref, NULL, out_estimator);
+ fp_to_meter, Tarad, Tref, register_paths, NULL, out_estimator);
}
}
@@ -192,10 +193,10 @@ sdis_solve_probe_green_function
if(!scn) return RES_BAD_ARG;
if(scene_is_2d(scn)) {
return solve_probe_2d(scn, nrealisations, position, NULL,
- fp_to_meter, Tarad, Tref, out_green, NULL);
+ fp_to_meter, Tarad, Tref, SDIS_HEAT_PATH_NONE, out_green, NULL);
} else {
return solve_probe_3d(scn, nrealisations, position, NULL,
- fp_to_meter, Tarad, Tref, out_green, NULL);
+ fp_to_meter, Tarad, Tref, SDIS_HEAT_PATH_NONE, out_green, NULL);
}
}
diff --git a/src/sdis_solve_Xd.h b/src/sdis_solve_Xd.h
@@ -157,6 +157,7 @@ XD(solve_probe)
const double fp_to_meter,/* Scale factor from floating point unit to meter */
const double Tarad, /* Ambient radiative temperature */
const double Tref, /* Reference temperature */
+ const int register_paths, /* Combination of enum sdis_heat_path_flag */
struct sdis_green_function** out_green, /* May be NULL <=> No green func */
struct sdis_estimator** out_estimator)
{
@@ -225,6 +226,12 @@ XD(solve_probe)
}
+ /* Create the estimator */
+ if(out_estimator) {
+ res = estimator_create(scn->dev, SDIS_ESTIMATOR_TEMPERATURE, &estimator);
+ if(res != RES_OK) goto error;
+ }
+
/* Here we go! Launch the Monte Carlo estimation */
omp_set_num_threads((int)scn->dev->nthreads);
#pragma omp parallel for schedule(static) reduction(+:weight,sqr_weight,N)
@@ -236,11 +243,17 @@ XD(solve_probe)
struct ssp_rng* rng = rngs[ithread];
struct green_path_handle* pgreen_path = NULL;
struct green_path_handle green_path = GREEN_PATH_HANDLE_NULL;
+ struct sdis_heat_path* pheat_path = NULL;
+ struct sdis_heat_path heat_path;
- if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occured */
+ if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */
if(!out_green) {
time = sample_time(rng, time_range);
+ if(register_paths) {
+ heat_path_init(scn->dev->allocator, &heat_path);
+ pheat_path = &heat_path;
+ }
} else {
/* Do not take care of the submitted time when registering the green
* function. Simply takes 0 as relative time */
@@ -252,7 +265,7 @@ XD(solve_probe)
}
res_local = XD(probe_realisation)(scn, rng, medium, position, time,
- fp_to_meter, Tarad, Tref, pgreen_path, &w);
+ fp_to_meter, Tarad, Tref, pgreen_path, pheat_path, &w);
if(res_local != RES_OK) {
if(res_local != RES_BAD_OP) { ATOMIC_SET(&res, res_local); continue; }
} else {
@@ -260,16 +273,26 @@ XD(solve_probe)
sqr_weight += w*w;
++N;
}
+
+ if(pheat_path) {
+ pheat_path->status = res_local == RES_OK
+ ? SDIS_HEAT_PATH_SUCCEED
+ : SDIS_HEAT_PATH_FAILED;
+
+ /* Check if the path must be saved regarding the register_paths mask */
+ if(!(register_paths & (int)pheat_path->status)) {
+ heat_path_release(pheat_path);
+ } else { /* Register the sampled path */
+ res_local = estimator_add_and_release_heat_path(estimator, pheat_path);
+ if(res_local != RES_OK) { ATOMIC_SET(&res, res_local); continue; }
+ }
+ }
}
if(res != RES_OK) goto error;
+ /* Setup the estimated temperature */
if(out_estimator) {
- /* Create the estimator */
- res = estimator_create
- (scn->dev, SDIS_ESTIMATOR_TEMPERATURE, nrealisations, N, &estimator);
- if(res != RES_OK) goto error;
-
- /* Setup the estimated temperature */
+ estimator_setup_realisations_count(estimator, nrealisations, N);
estimator_setup_temperature(estimator, weight, sqr_weight);
}
@@ -282,7 +305,7 @@ XD(solve_probe)
}
/* Finalize the estimated green */
- res = green_function_finalize(green);
+ res = green_function_finalize(green, rng_proxy);
if(res != RES_OK) goto error;
}
@@ -407,6 +430,10 @@ XD(solve_probe_boundary)
if(res != RES_OK) goto error;
}
+ /* Create the estimator */
+ res = estimator_create(scn->dev, SDIS_ESTIMATOR_TEMPERATURE, &estimator);
+ if(res != RES_OK) goto error;
+
/* Here we go! Launch the Monte Carlo estimation */
omp_set_num_threads((int)scn->dev->nthreads);
#pragma omp parallel for schedule(static) reduction(+:weight,sqr_weight,N)
@@ -436,12 +463,8 @@ XD(solve_probe_boundary)
}
if(res != RES_OK) goto error;
- /* Create the estimator */
- res = estimator_create
- (scn->dev, SDIS_ESTIMATOR_TEMPERATURE, nrealisations, N, &estimator);
- if(res != RES_OK) goto error;
-
/* Setup the estimated temperature */
+ estimator_setup_realisations_count(estimator, nrealisations, N);
estimator_setup_temperature(estimator, weight, sqr_weight);
exit:
@@ -564,6 +587,10 @@ XD(solve_boundary)
if(res != RES_OK) goto error;
}
+ /* Create the estimator */
+ res = estimator_create(scn->dev, SDIS_ESTIMATOR_TEMPERATURE, &estimator);
+ if(res != RES_OK) goto error;
+
omp_set_num_threads((int)scn->dev->nthreads);
#pragma omp parallel for schedule(static) reduction(+:weight,sqr_weight,N)
for(irealisation=0; irealisation<(int64_t)nrealisations; ++irealisation) {
@@ -621,12 +648,8 @@ XD(solve_boundary)
}
}
- /* Create the estimator */
- res = estimator_create
- (scn->dev, SDIS_ESTIMATOR_TEMPERATURE, nrealisations, N, &estimator);
- if(res != RES_OK) goto error;
-
/* Setup the estimated temperature */
+ estimator_setup_realisations_count(estimator, nrealisations, N);
estimator_setup_temperature(estimator, weight, sqr_weight);
exit:
@@ -761,6 +784,10 @@ XD(solve_probe_boundary_flux)
res = XD(interface_prebuild_fragment)
(&frag, scn, (unsigned)iprim, uv, fluid_side);
+ /* Create the estimator */
+ res = estimator_create(scn->dev, SDIS_ESTIMATOR_FLUX, &estimator);
+ if(res != RES_OK) goto error;
+
/* Here we go! Launch the Monte Carlo estimation */
omp_set_num_threads((int)scn->dev->nthreads);
#pragma omp parallel for schedule(static) reduction(+:weight_t,sqr_weight_t,\
@@ -814,12 +841,8 @@ XD(solve_probe_boundary_flux)
}
if(res != RES_OK) goto error;
- /* Create the estimator */
- res = estimator_create
- (scn->dev, SDIS_ESTIMATOR_FLUX, nrealisations, N, &estimator);
- if(res != RES_OK) goto error;
-
/* Setup the estimated values */
+ estimator_setup_realisations_count(estimator, nrealisations, N);
estimator_setup_temperature(estimator, weight_t, sqr_weight_t);
estimator_setup_flux(estimator, FLUX_CONVECTIVE, weight_fc, sqr_weight_fc);
estimator_setup_flux(estimator, FLUX_RADIATIVE, weight_fr, sqr_weight_fr);
@@ -947,6 +970,10 @@ XD(solve_boundary_flux)
if(res != RES_OK) goto error;
}
+ /* Create the estimator */
+ res = estimator_create(scn->dev, SDIS_ESTIMATOR_FLUX, &estimator);
+ if(res != RES_OK) goto error;
+
omp_set_num_threads((int)scn->dev->nthreads);
#pragma omp parallel for schedule(static) reduction(+:weight_t,sqr_weight_t,\
weight_fc,sqr_weight_fc,weight_fr,sqr_weight_fr,weight_f,sqr_weight_f,N)
@@ -1045,12 +1072,8 @@ XD(solve_boundary_flux)
}
if(res != RES_OK) goto error;
- /* Create the estimator */
- res = estimator_create
- (scn->dev, SDIS_ESTIMATOR_FLUX, nrealisations, N, &estimator);
- if(res != RES_OK) goto error;
-
/* Setup the estimated values */
+ estimator_setup_realisations_count(estimator, nrealisations, N);
estimator_setup_temperature(estimator, weight_t, sqr_weight_t);
estimator_setup_flux(estimator, FLUX_CONVECTIVE, weight_fc, sqr_weight_fc);
estimator_setup_flux(estimator, FLUX_RADIATIVE, weight_fr, sqr_weight_fr);
diff --git a/src/test_sdis_conducto_radiative.c b/src/test_sdis_conducto_radiative.c
@@ -398,7 +398,7 @@ main(int argc, char** argv)
pos[1] = ssp_rng_uniform_double(rng, -0.9, 0.9);
pos[2] = ssp_rng_uniform_double(rng, -0.9, 0.9);
- OK(sdis_solve_probe(scn, N, pos, time_range, 1, -1, Tref, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1, -1, Tref, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
OK(sdis_estimator_get_temperature(estimator, &T));
@@ -422,6 +422,10 @@ main(int argc, char** argv)
OK(sdis_estimator_ref_put(estimator));
OK(sdis_estimator_ref_put(estimator2));
OK(sdis_green_function_ref_put(green));
+
+ OK(sdis_solve_probe(scn, 10, pos, time_range, 1, -1, Tref,
+ SDIS_HEAT_PATH_ALL, &estimator));
+ OK(sdis_estimator_ref_put(estimator));
}
/* Release memory */
diff --git a/src/test_sdis_conducto_radiative_2d.c b/src/test_sdis_conducto_radiative_2d.c
@@ -403,7 +403,7 @@ main(int argc, char** argv)
pos[0] = ssp_rng_uniform_double(rng, -0.9, 0.9);
pos[1] = ssp_rng_uniform_double(rng, -0.9, 0.9);
- OK(sdis_solve_probe(scn, 10000, pos, time_range, 1, -1, Tref, &estimator));
+ OK(sdis_solve_probe(scn, 10000, pos, time_range, 1, -1, Tref, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
OK(sdis_estimator_get_temperature(estimator, &T));
@@ -427,6 +427,10 @@ main(int argc, char** argv)
OK(sdis_estimator_ref_put(estimator));
OK(sdis_estimator_ref_put(estimator2));
OK(sdis_green_function_ref_put(green));
+
+ OK(sdis_solve_probe
+ (scn, 10, pos, time_range, 1, -1, Tref, SDIS_HEAT_PATH_ALL, &estimator));
+ OK(sdis_estimator_ref_put(estimator));
}
/* Release memory */
diff --git a/src/test_sdis_convection.c b/src/test_sdis_convection.c
@@ -277,7 +277,7 @@ main(int argc, char** argv)
*((int*)sdis_data_get(is_stationary)) = IS_INF(time);
/* Solve in 3D */
- OK(sdis_solve_probe(box_scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(box_scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
OK(sdis_estimator_get_temperature(estimator, &T));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
@@ -313,7 +313,7 @@ main(int argc, char** argv)
*((int*)sdis_data_get(is_stationary)) = IS_INF(time);
/* Solve in 2D */
- OK(sdis_solve_probe(square_scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(square_scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
CHK(nfails + nreals == N);
diff --git a/src/test_sdis_convection_non_uniform.c b/src/test_sdis_convection_non_uniform.c
@@ -293,7 +293,7 @@ main(int argc, char** argv)
*((int*)sdis_data_get(is_stationary)) = IS_INF(time);
/* Solve in 3D */
- OK(sdis_solve_probe(box_scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(box_scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
CHK(nfails + nreals == N);
@@ -327,7 +327,7 @@ main(int argc, char** argv)
*((int*)sdis_data_get(is_stationary)) = IS_INF(time);
- OK(sdis_solve_probe(square_scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(square_scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
CHK(nfails + nreals == N);
diff --git a/src/test_sdis_flux.c b/src/test_sdis_flux.c
@@ -145,7 +145,7 @@ solve(struct sdis_scene* scn, const double pos[])
ref = T0 + (1 - pos[0]) * PHI/LAMBDA;
time_current(&t0);
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
time_sub(&t0, time_current(&t1), &t0);
time_dump(&t0, TIME_ALL, NULL, dump, sizeof(dump));
diff --git a/src/test_sdis_solve_probe.c b/src/test_sdis_solve_probe.c
@@ -167,6 +167,80 @@ interface_get_specular_fraction
}
/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+struct dump_path_context {
+ FILE* stream;
+ size_t offset;
+ size_t nfailures;
+ size_t nsuccesses;
+};
+static const struct dump_path_context DUMP_PATH_CONTEXT_NULL = {NULL, 0, 0, 0};
+
+static res_T
+dump_vertex_pos(const struct sdis_heat_vertex* vert, void* context)
+{
+ struct dump_path_context* ctx = context;
+ CHK(vert && context);
+ fprintf(ctx->stream, "v %g %g %g\n", SPLIT3(vert->P));
+ return RES_OK;
+}
+
+static res_T
+process_heat_path(const struct sdis_heat_path* path, void* context)
+{
+ struct dump_path_context* ctx = context;
+ struct sdis_heat_vertex vert = SDIS_HEAT_VERTEX_NULL;
+ enum sdis_heat_path_flag status = SDIS_HEAT_PATH_NONE;
+ size_t i;
+ size_t n;
+ (void)context;
+
+ CHK(path && context);
+
+ BA(sdis_heat_path_get_vertices_count(NULL, &n));
+ BA(sdis_heat_path_get_vertices_count(path, NULL));
+ OK(sdis_heat_path_get_vertices_count(path, &n));
+ CHK(n != 0);
+
+ BA(sdis_heat_path_get_status(NULL, &status));
+ BA(sdis_heat_path_get_status(path, NULL));
+ OK(sdis_heat_path_get_status(path, &status));
+ CHK(status == SDIS_HEAT_PATH_SUCCEED || status == SDIS_HEAT_PATH_FAILED);
+
+ switch(status) {
+ case SDIS_HEAT_PATH_FAILED: ++ctx->nfailures; break;
+ case SDIS_HEAT_PATH_SUCCEED: ++ctx->nsuccesses; break;
+ default: FATAL("Unreachable code.\n"); break;
+ }
+
+ BA(sdis_heat_path_get_vertex(NULL, 0, &vert));
+ BA(sdis_heat_path_get_vertex(path, n, &vert));
+ BA(sdis_heat_path_get_vertex(path, 0, NULL));
+
+ FOR_EACH(i, 0, n) {
+ OK(sdis_heat_path_get_vertex(path, i, &vert));
+ CHK(vert.type == SDIS_HEAT_VERTEX_CONVECTION
+ || vert.type == SDIS_HEAT_VERTEX_CONDUCTION
+ || vert.type == SDIS_HEAT_VERTEX_RADIATIVE);
+ }
+
+ BA(sdis_heat_path_for_each_vertex(NULL, dump_vertex_pos, context));
+ BA(sdis_heat_path_for_each_vertex(path, NULL, context));
+ OK(sdis_heat_path_for_each_vertex(path, dump_vertex_pos, context));
+
+ FOR_EACH(i, 0, n-1) {
+ fprintf(ctx->stream, "l %lu %lu\n",
+ (unsigned long)(i+1 + ctx->offset),
+ (unsigned long)(i+2 + ctx->offset));
+ }
+
+ ctx->offset += n;
+
+ return RES_OK;
+}
+
+/*******************************************************************************
* Test
******************************************************************************/
int
@@ -184,9 +258,11 @@ main(int argc, char** argv)
struct sdis_estimator* estimator = NULL;
struct sdis_estimator* estimator2 = NULL;
struct sdis_green_function* green = NULL;
+ const struct sdis_heat_path* path = NULL;
struct sdis_fluid_shader fluid_shader = DUMMY_FLUID_SHADER;
struct sdis_solid_shader solid_shader = DUMMY_SOLID_SHADER;
struct sdis_interface_shader interface_shader = SDIS_INTERFACE_SHADER_NULL;
+ struct dump_path_context dump_ctx = DUMP_PATH_CONTEXT_NULL;
struct context ctx;
struct fluid* fluid_param;
struct solid* solid_param;
@@ -196,12 +272,14 @@ main(int argc, char** argv)
double time_range[2];
double ref;
const size_t N = 1000;
+ const size_t N_dump = 10;
size_t nreals;
size_t nfails;
+ size_t n;
(void)argc, (void)argv;
OK(mem_init_proxy_allocator(&allocator, &mem_default_allocator));
- OK(sdis_device_create(NULL, &allocator, SDIS_NTHREADS_DEFAULT, 1, &dev));
+ OK(sdis_device_create(NULL, &allocator, SDIS_NTHREADS_DEFAULT, 0, &dev));
/* Create the fluid medium */
OK(sdis_data_create
@@ -263,13 +341,13 @@ main(int argc, char** argv)
pos[1] = 0.5;
pos[2] = 0.5;
time_range[0] = time_range[1] = INF;
- BA(sdis_solve_probe(NULL, N, pos, time_range, 1.0, 0, 0, &estimator));
- BA(sdis_solve_probe(scn, 0, pos, time_range, 1.0, 0, 0, &estimator));
- BA(sdis_solve_probe(scn, N, NULL, time_range, 1.0, 0, 0, &estimator));
- BA(sdis_solve_probe(scn, N, pos, time_range, 0, 0, 0, &estimator));
- BA(sdis_solve_probe(scn, N, pos, time_range, 0, 0, -1, &estimator));
- BA(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, NULL));
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ BA(sdis_solve_probe(NULL, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
+ BA(sdis_solve_probe(scn, 0, pos, time_range, 1.0, 0, 0, 0, &estimator));
+ BA(sdis_solve_probe(scn, N, NULL, time_range, 1.0, 0, 0, 0, &estimator));
+ BA(sdis_solve_probe(scn, N, pos, time_range, 0, 0, 0, 0, &estimator));
+ BA(sdis_solve_probe(scn, N, pos, time_range, 0, 0, -1, 0, &estimator));
+ BA(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, NULL));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
BA(sdis_estimator_get_type(estimator, NULL));
BA(sdis_estimator_get_type(NULL, &type));
@@ -318,10 +396,10 @@ main(int argc, char** argv)
/* The external fluid cannot have an unknown temperature */
fluid_param->temperature = -1;
- BA(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ BA(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
fluid_param->temperature = 300;
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
BA(sdis_solve_probe_green_function(NULL, N, pos, 1.0, 0, 0, &green));
BA(sdis_solve_probe_green_function(scn, 0, pos, 1.0, 0, 0, &green));
@@ -348,6 +426,29 @@ main(int argc, char** argv)
OK(sdis_estimator_ref_put(estimator));
OK(sdis_estimator_ref_put(estimator2));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
+ BA(sdis_estimator_get_paths_count(NULL, &n));
+ BA(sdis_estimator_get_paths_count(estimator, NULL));
+ OK(sdis_estimator_get_paths_count(estimator, &n));
+ CHK(n == 0);
+ OK(sdis_estimator_ref_put(estimator));
+
+ OK(sdis_solve_probe(scn, N_dump, pos, time_range, 1.0, 0, 0,
+ SDIS_HEAT_PATH_ALL, &estimator));
+ OK(sdis_estimator_get_paths_count(estimator, &n));
+ CHK(n == N_dump);
+
+ BA(sdis_estimator_get_path(NULL, 0, &path));
+ BA(sdis_estimator_get_path(estimator, n, &path));
+ BA(sdis_estimator_get_path(estimator, 0, NULL));
+ OK(sdis_estimator_get_path(estimator, 0, &path));
+
+ dump_ctx.stream = stderr;
+ BA(sdis_estimator_for_each_path(NULL, process_heat_path, &dump_ctx));
+ BA(sdis_estimator_for_each_path(estimator, NULL, &dump_ctx));
+ OK(sdis_estimator_for_each_path(estimator, process_heat_path, &dump_ctx));
+
+ OK(sdis_estimator_ref_put(estimator));
OK(sdis_scene_ref_put(scn));
OK(sdis_device_ref_put(dev));
diff --git a/src/test_sdis_solve_probe2.c b/src/test_sdis_solve_probe2.c
@@ -242,7 +242,7 @@ main(int argc, char** argv)
pos[1] = 0.5;
pos[2] = 0.5;
time_range[0] = time_range[1] = INF;
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, -1, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, -1, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
OK(sdis_estimator_get_temperature(estimator, &T));
diff --git a/src/test_sdis_solve_probe2_2d.c b/src/test_sdis_solve_probe2_2d.c
@@ -240,7 +240,7 @@ main(int argc, char** argv)
pos[0] = 0.5;
pos[1] = 0.5;
time_range[0] = time_range[1] = INF;
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, -1, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, -1, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
OK(sdis_estimator_get_temperature(estimator, &T));
diff --git a/src/test_sdis_solve_probe3.c b/src/test_sdis_solve_probe3.c
@@ -298,7 +298,7 @@ main(int argc, char** argv)
pos[1] = 0.5;
pos[2] = 0.5;
time_range[0] = time_range[1] = INF;
- OK(sdis_solve_probe( scn, N, pos, time_range, 1.0, -1, 0, &estimator));
+ OK(sdis_solve_probe( scn, N, pos, time_range, 1.0, -1, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
OK(sdis_estimator_get_temperature(estimator, &T));
diff --git a/src/test_sdis_solve_probe3_2d.c b/src/test_sdis_solve_probe3_2d.c
@@ -291,7 +291,7 @@ main(int argc, char** argv)
pos[0] = 0.5;
pos[1] = 0.5;
time_range[0] = time_range[1] = INF;
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, -1, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, -1, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
OK(sdis_estimator_get_temperature(estimator, &T));
diff --git a/src/test_sdis_solve_probe_2d.c b/src/test_sdis_solve_probe_2d.c
@@ -201,7 +201,7 @@ main(int argc, char** argv)
pos[0] = 0.5;
pos[1] = 0.5;
time_range[0] = time_range[1] = INF;
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
@@ -228,7 +228,7 @@ main(int argc, char** argv)
/* The external fluid cannot have an unknown temperature */
fluid_param->temperature = -1;
- BA(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ BA(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
OK(sdis_scene_ref_put(scn));
OK(sdis_device_ref_put(dev));
diff --git a/src/test_sdis_utils.c b/src/test_sdis_utils.c
@@ -16,6 +16,12 @@
#include "test_sdis_utils.h"
#include <rsys/math.h>
+enum heat_vertex_attrib {
+ HEAT_VERTEX_WEIGHT,
+ HEAT_VERTEX_TIME,
+ HEAT_VERTEX_TYPE
+};
+
/*******************************************************************************
* Helper functions
******************************************************************************/
@@ -141,6 +147,71 @@ solve_green_path(struct sdis_green_path* path, void* ctx)
return RES_OK;
}
+static void
+dump_heat_path_position(FILE* stream, const struct sdis_heat_path* path)
+{
+ size_t nverts;
+ size_t ivert;
+
+ CHK(stream && path);
+
+ OK(sdis_heat_path_get_vertices_count(path, &nverts));
+ FOR_EACH(ivert, 0, nverts) {
+ struct sdis_heat_vertex vtx;
+ OK(sdis_heat_path_get_vertex(path, ivert, &vtx));
+ fprintf(stream, "%g %g %g\n", SPLIT3(vtx.P));
+ }
+}
+
+static void
+dump_heat_path_segments
+ (FILE* stream, const struct sdis_heat_path* path, const size_t offset)
+{
+ size_t nverts, ivert;
+
+ CHK(stream);
+
+ OK(sdis_heat_path_get_vertices_count(path, &nverts));
+ fprintf(stream, "%lu", (unsigned long)nverts);
+ FOR_EACH(ivert, 0, nverts) {
+ fprintf(stream, " %lu", (unsigned long)(ivert + offset));
+ }
+ fprintf(stream, "\n");
+}
+
+static void
+dump_heat_path_vertex_attribs
+ (FILE* stream,
+ const struct sdis_heat_path* path,
+ const enum heat_vertex_attrib attr)
+{
+ size_t nverts, ivert;
+ CHK(stream && path);
+
+ OK(sdis_heat_path_get_vertices_count(path, &nverts));
+ FOR_EACH(ivert, 0, nverts) {
+ struct sdis_heat_vertex vtx;
+ OK(sdis_heat_path_get_vertex(path, ivert, &vtx));
+ switch(attr) {
+ case HEAT_VERTEX_WEIGHT:
+ fprintf(stream, "%g\n", vtx.weight);
+ break;
+ case HEAT_VERTEX_TIME:
+ fprintf(stream, "%g\n", IS_INF(vtx.time) ? FLT_MAX : vtx.time);
+ break;
+ case HEAT_VERTEX_TYPE:
+ switch(vtx.type) {
+ case SDIS_HEAT_VERTEX_CONDUCTION: fprintf(stream, "0.0\n"); break;
+ case SDIS_HEAT_VERTEX_CONVECTION: fprintf(stream, "0.5\n"); break;
+ case SDIS_HEAT_VERTEX_RADIATIVE: fprintf(stream, "1.0\n"); break;
+ default: FATAL("Unreachable code.\n"); break;
+ }
+ break;
+ default: FATAL("Unreachable code.\n"); break;
+ }
+ }
+}
+
/*******************************************************************************
* Local function
******************************************************************************/
@@ -194,3 +265,97 @@ check_green_function(struct sdis_green_function* green)
OK(sdis_estimator_ref_put(estimator));
}
+void
+dump_heat_paths(FILE* stream, struct sdis_estimator* estimator)
+{
+ const struct sdis_heat_path* path;
+ size_t ipath;
+ size_t npaths;
+ size_t nvertices;
+ size_t offset;
+ size_t n;
+ CHK(stream && estimator);
+
+ OK(sdis_estimator_get_paths_count(estimator, &npaths));
+ CHK(npaths);
+
+ /* Header */
+ fprintf(stream, "# vtk DataFile Version 2.0\n");
+ fprintf(stream, "Heat paths\n");
+ fprintf(stream, "ASCII\n");
+ fprintf(stream, "DATASET POLYDATA\n");
+
+ /* Compute the overall number of vertices */
+ nvertices = 0;
+ FOR_EACH(ipath, 0, npaths) {
+ OK(sdis_estimator_get_path(estimator, ipath, &path));
+ OK(sdis_heat_path_get_vertices_count(path, &n));
+ nvertices += n;
+ }
+
+ /* Write path positions */
+ fprintf(stream, "POINTS %lu double\n", (unsigned long)nvertices);
+ FOR_EACH(ipath, 0, npaths) {
+ OK(sdis_estimator_get_path(estimator, ipath, &path));
+ dump_heat_path_position(stream, path);
+ }
+
+ /* Write the segment of the paths */
+ fprintf(stream, "LINES %lu %lu\n",
+ (unsigned long)npaths, (unsigned long)(npaths + nvertices));
+ offset = 0;
+ FOR_EACH(ipath, 0, npaths) {
+ OK(sdis_estimator_get_path(estimator, ipath, &path));
+ dump_heat_path_segments(stream, path, offset);
+ OK(sdis_heat_path_get_vertices_count(path, &n));
+ offset += n;
+ }
+
+ fprintf(stream, "POINT_DATA %lu\n", (unsigned long)nvertices);
+
+ /* Write the type of the random walk vertices */
+ fprintf(stream, "SCALARS Vertex_Type float 1\n");
+ fprintf(stream, "LOOKUP_TABLE vertex_type\n");
+ FOR_EACH(ipath, 0, npaths) {
+ OK(sdis_estimator_get_path(estimator, ipath, &path));
+ dump_heat_path_vertex_attribs(stream, path, HEAT_VERTEX_TYPE);
+ }
+ fprintf(stream, "LOOKUP_TABLE vertex_type 3\n");
+ fprintf(stream, "0.0 1.0 1.0 1.0\n"); /* 0.0 = Magenta: conduction */
+ fprintf(stream, "1.0 1.0 0.0 1.0\n"); /* 0.5 = Yellow: convection */
+ fprintf(stream, "1.0 0.0 1.0 1.0\n"); /* 1.0 = Purple: radiative */
+
+ /* Write the weights of the random walk vertices */
+ fprintf(stream, "SCALARS Weight double 1\n");
+ fprintf(stream, "LOOKUP_TABLE default\n");
+ FOR_EACH(ipath, 0, npaths) {
+ OK(sdis_estimator_get_path(estimator, ipath, &path));
+ dump_heat_path_vertex_attribs(stream, path, HEAT_VERTEX_WEIGHT);
+ }
+
+ /* Write the time of the random walk vertices */
+ fprintf(stream, "SCALARS Time double 1\n");
+ fprintf(stream, "LOOKUP_TABLE default\n");
+ FOR_EACH(ipath, 0, npaths) {
+ OK(sdis_estimator_get_path(estimator, ipath, &path));
+ dump_heat_path_vertex_attribs(stream, path, HEAT_VERTEX_TIME);
+ }
+
+ /* Write path type */
+ fprintf(stream, "CELL_DATA %lu\n", (unsigned long)npaths);
+ fprintf(stream, "SCALARS Path_Type float 1\n");
+ fprintf(stream, "LOOKUP_TABLE path_type\n");
+ FOR_EACH(ipath, 0, npaths) {
+ enum sdis_heat_path_flag status = SDIS_HEAT_PATH_NONE;
+ OK(sdis_estimator_get_path(estimator, ipath, &path));
+ OK(sdis_heat_path_get_status(path, &status));
+ switch(status) {
+ case SDIS_HEAT_PATH_SUCCEED: fprintf(stream, "0.0\n"); break;
+ case SDIS_HEAT_PATH_FAILED: fprintf(stream, "1.0\n"); break;
+ default: FATAL("Unreachable code.\n"); break;
+ }
+ }
+ fprintf(stream, "LOOKUP_TABLE path_type 2\n");
+ fprintf(stream, "0.0 0.0 1.0 1.0\n"); /* 0.0 = Bleu: success */
+ fprintf(stream, "1.0 0.0 0.0 1.0\n"); /* 1.0 = Red: failure */
+}
diff --git a/src/test_sdis_utils.h b/src/test_sdis_utils.h
@@ -278,7 +278,7 @@ static INLINE void
check_memory_allocator(struct mem_allocator* allocator)
{
if(MEM_ALLOCATED_SIZE(allocator)) {
- char dump[128];
+ char dump[1024];
MEM_DUMP(allocator, dump, sizeof(dump));
fprintf(stderr, "%s\n", dump);
FATAL("Memory leaks.\n");
@@ -289,5 +289,10 @@ extern LOCAL_SYM void
check_green_function
(struct sdis_green_function* green);
+extern LOCAL_SYM void
+dump_heat_paths
+ (FILE* stream,
+ struct sdis_estimator* estimator);
+
#endif /* TEST_SDIS_UTILS_H */
diff --git a/src/test_sdis_volumic_power.c b/src/test_sdis_volumic_power.c
@@ -169,7 +169,7 @@ solve(struct sdis_scene* scn, const double pos[])
ref = P0 / (2*LAMBDA) * (1.0/4.0 - x*x) + T0;
time_current(&t0);
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.0, 0, 0, 0, &estimator));
time_sub(&t0, time_current(&t1), &t0);
time_dump(&t0, TIME_ALL, NULL, dump, sizeof(dump));
diff --git a/src/test_sdis_volumic_power2.c b/src/test_sdis_volumic_power2.c
@@ -242,7 +242,7 @@ check(struct sdis_scene* scn, const struct reference refs[], const size_t nrefs)
pos[1] = refs[i].pos[1];
pos[2] = refs[i].pos[2];
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, 0, &estimator));
OK(sdis_estimator_get_temperature(estimator, &T));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
@@ -256,7 +256,6 @@ check(struct sdis_scene* scn, const struct reference refs[], const size_t nrefs)
}
}
-
int
main(int argc, char** argv)
{
diff --git a/src/test_sdis_volumic_power2_2d.c b/src/test_sdis_volumic_power2_2d.c
@@ -261,7 +261,7 @@ check(struct sdis_scene* scn, const struct reference refs[], const size_t nrefs)
pos[0] = refs[i].pos[0];
pos[1] = refs[i].pos[1];
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, 0, &estimator));
OK(sdis_estimator_get_temperature(estimator, &T));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
diff --git a/src/test_sdis_volumic_power3_2d.c b/src/test_sdis_volumic_power3_2d.c
@@ -441,7 +441,7 @@ main(int argc, char** argv)
FATAL("Unreachable code.\n");
}
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, 0, &estimator));
OK(sdis_estimator_get_temperature(estimator, &T));
OK(sdis_estimator_get_failure_count(estimator, &nfails));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
diff --git a/src/test_sdis_volumic_power4_2d.c b/src/test_sdis_volumic_power4_2d.c
@@ -348,7 +348,7 @@ main(int argc, char** argv)
Tref = T2 + (T1-T2)/L * (pos[1] + vertices[3]);
#endif
- OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, &estimator));
+ OK(sdis_solve_probe(scn, N, pos, time_range, 1.f, -1, 0, 0, &estimator));
OK(sdis_estimator_get_temperature(estimator, &T));
OK(sdis_estimator_get_realisation_count(estimator, &nreals));
OK(sdis_estimator_get_failure_count(estimator, &nfails));