commit a5e41da8bd3e495a7f3d1f204bc12942900c1607
parent 8c6fcd5e3a3b97f9dafc019587d0de2189eb989a
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Tue, 25 Jun 2024 14:46:37 +0200
Add a test to the custom sampling of a conductive path
Use a programmable solid on which the conductive path is sampled on the
caller's side.
Diffstat:
4 files changed, 719 insertions(+), 3 deletions(-)
diff --git a/Makefile b/Makefile
@@ -20,6 +20,7 @@ include config.mk
TESTS =\
sadist_conducto_radiative\
+ sadist_custom_conductive_path\
sadist_external_flux\
sadist_probe_boundary\
sadist_unsteady
@@ -27,6 +28,7 @@ TESTS =\
LIBS =\
libsadist_radenv_1d.so\
libsadist_solid_fluid.so\
+ libsadist_solid_sphere.so\
libsadist_spherical_source.so\
libsadist_trilinear_profile.so\
libsadist_unsteady_profile.so\
@@ -37,8 +39,14 @@ default: .config $(TESTS) $(LIBS)
.config:
@if ! $(PKG_CONFIG) --atleast-version $(RSYS_VERSION) rsys; then \
echo "rsys $(RSYS_VERSION) not found" >&2; exit 1; fi
+ @if ! $(PKG_CONFIG) --atleast-version $(S3D_VERSION) s3d; then \
+ echo "s3d $(S3D_VERSION) not found" >&2; exit 1; fi
@if ! $(PKG_CONFIG) --atleast-version $(S3DUT_VERSION) s3dut; then \
echo "s3dut $(S3DUT_VERSION) not found" >&2; exit 1; fi
+ @if ! $(PKG_CONFIG) --atleast-version $(SSP_VERSION) star-sp; then \
+ echo "star-sp $(SSP_VERSION) not found" >&2; exit 1; fi
+ @if ! $(PKG_CONFIG) --atleast-version $(SSTL_VERSION) sstl; then \
+ echo "sstl $(SSTL_VERSION) not found" >&2; exit 1; fi
@echo "config done" > $@
clean:
@@ -72,6 +80,15 @@ libsadist_solid_fluid.so: config.mk src/sadist_lib_solid_fluid.o
$(CC) $(CFLAGS_SO) $(RSYS_CFLAGS) -o $@ src/sadist_lib_solid_fluid.o \
$(LDFLAGS_SO) $(RSYS_LIBS)
+src/sadist_lib_solid_sphere.o: config.mk src/sadist_lib_solid_sphere.c
+ $(CC) $(CFLAGS_SO) $(S3D_CFLAGS) $(SSP_CFLAGS) $(SSTL_CFLAGS) $(RSYS_CFLAGS) \
+ -c $(@:.o=.c) -o $@
+
+libsadist_solid_sphere.so: config.mk src/sadist_lib_solid_sphere.o
+ $(CC) $(CFLAGS_SO) $(S3D_CFLAGS) $(SSP_CFLAGS) $(SSTL_CFLAGS) $(RSYS_CFLAGS) \
+ -o $@ src/sadist_lib_solid_sphere.o $(LDFLAGS_SO) $(S3D_LIBS) $(SSP_LIBS) \
+ $(SSTL_LIBS) $(RSYS_LIBS) -lm
+
src/sadist_lib_spherical_source.o: config.mk src/sadist_lib_spherical_source.c
$(CC) $(CFLAGS_SO) $(RSYS_CFLAGS) -c $(@:.o=.c) -o $@
@@ -105,6 +122,14 @@ src/sadist_conducto_radiative.o:\
sadist_conducto_radiative: config.mk src/sadist_conducto_radiative.o
$(CC) $(CFLAGS_EXE) -o $@ src/$@.o $(LDFLAGS_EXE) $(RSYS_LIBS)
+src/sadist_custom_conductive_path.o:\
+ config.mk src/sadist.h src/sadist_custom_conductive_path.c
+ $(CC) $(CFLAGS_EXE) $(S3DUT_CFLAGS) $(RSYS_CFLAGS) -c $(@:.o=.c) -o $@
+
+sadist_custom_conductive_path: config.mk src/sadist_custom_conductive_path.o
+ $(CC) $(CFLAGS_EXE) $(S3DUT_CFLAGS) $(RSYS_CFLAGS) -o $@ src/$@.o \
+ $(LDFLAGS_EXE) $(S3DUT_LIBS) $(RSYS_LIBS)
+
src/sadist_external_flux.o:\
config.mk src/sadist.h src/sadist_external_flux.c
$(CC) $(CFLAGS_EXE) -c $(@:.o=.c) -o $@
diff --git a/config.mk b/config.mk
@@ -4,8 +4,8 @@ PREFIX = /usr/local
LIB_TYPE = SHARED
#LIB_TYPE = STATIC
-BUILD_TYPE = RELEASE
-#BUILD_TYPE = DEBUG
+#BUILD_TYPE = RELEASE
+BUILD_TYPE = DEBUG
################################################################################
# Tools
@@ -25,11 +25,22 @@ RSYS_VERSION = 0.14
RSYS_CFLAGS = $$($(PKG_CONFIG) $(PCFLAGS) --cflags rsys)
RSYS_LIBS = $$($(PKG_CONFIG) $(PCFLAGS) --libs rsys)
-# For tests
+S3D_VERSION = 0.10
+S3D_CFLAGS = $$($(PKG_CONFIG) $(PCFLAGS) --cflags s3d)
+S3D_LIBS = $$($(PKG_CONFIG) $(PCFLAGS) --libs s3d)
+
S3DUT_VERSION = 0.4
S3DUT_CFLAGS = $$($(PKG_CONFIG) $(PCFLAGS) --cflags s3dut)
S3DUT_LIBS = $$($(PKG_CONFIG) $(PCFLAGS) --libs s3dut)
+SSP_VERSION = 0.14
+SSP_CFLAGS = $$($(PKG_CONFIG) $(PCFLAGS) --cflags star-sp)
+SSP_LIBS = $$($(PKG_CONFIG) $(PCFLAGS) --libs star-sp)
+
+SSTL_VERSION = 0.5
+SSTL_CFLAGS = $$($(PKG_CONFIG) $(PCFLAGS) --cflags sstl)
+SSTL_LIBS = $$($(PKG_CONFIG) $(PCFLAGS) --libs sstl)
+
################################################################################
# Compilation options
################################################################################
diff --git a/src/sadist_custom_conductive_path.c b/src/sadist_custom_conductive_path.c
@@ -0,0 +1,283 @@
+/* Copyright (C) 2024 |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/>. */
+
+#define _POSIX_C_SOURCE 200112L /* popen */
+
+#include "sadist.h"
+
+#include <stardis/stardis-prog-properties.h>
+
+#include <star/s3dut.h>
+
+#include <rsys/double2.h>
+#include <rsys/double3.h>
+#include <rsys/math.h>
+
+#include <errno.h>
+#include <float.h>
+#include <getopt.h>
+#include <string.h> /* strerror */
+#include <wait.h> /* WIFEXITED, WEXITSTATUS */
+
+/*
+ * The system is a trilinear profile of the temperature at steady state, i.e. at
+ * each point of the system we can calculate the temperature analytically. Two
+ * forms are immersed in this temperature field: a super shape and a sphere
+ * included in the super shape. On the Monte Carlo side, the temperature is
+ * unknown everywhere except on the surface of the super shape whose
+ * temperature is defined from the aformentionned trilinear profile.
+ *
+ * We will estimate the temperature at the position of a probe in solids by
+ * providing a user-side function to sample the conductive path in the sphere.
+ * We should find the temperature of the trilinear profile at the probe position
+ * by Monte Carlo, independently of this coupling with an external path sampling
+ * routine.
+ *
+ *
+ * /\ <-- T(x,y,z)
+ * ___/ \___
+ * T(z) \ __ /
+ * | T(y) T=?/. \ /
+ * |/ / \__/ \
+ * o--- T(x) /_ __ _\
+ * \/ \/
+ */
+
+#define FILENAME_SSHAPE "sshape.stl"
+#define FILENAME_SPHERE "sphere.stl"
+#define FILENAME_SCENE "scene.txt"
+
+/* Temperature at the upper bound of the X, Y and Z axis. The temperature at the
+ * lower bound is implicitly 0 K */
+#define TX 333.0 /* [K] */
+#define TY 432.0 /* [K] */
+#define TZ 579.0 /* [K] */
+
+/* Probe position */
+#define PX 0.125
+#define PY 0.250
+#define PZ 0.375
+
+/* Commands to calculate probe temperature */
+#define COMMAND "stardis -V3 -p "STR(PX)","STR(PY)","STR(PZ)" -M "FILENAME_SCENE
+
+/* Axis Aligned Bounding Box */
+struct aabb {
+ double lower[3];
+ double upper[3];
+};
+#define AABB_NULL__ {{DBL_MAX,DBL_MAX,DBL_MAX}, {-DBL_MAX,-DBL_MAX,-DBL_MAX}}
+static const struct aabb AABB_NULL = AABB_NULL__;
+
+
+/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+static void
+aabb_update(struct aabb* aabb, const struct s3dut_mesh_data* mesh)
+{
+ size_t ivertex = 0;
+ ASSERT(aabb);
+
+ FOR_EACH(ivertex, 0, mesh->nvertices) {
+ const double* vertex = mesh->positions + ivertex*3;
+ aabb->lower[0] = MMIN(aabb->lower[0], vertex[0]);
+ aabb->lower[1] = MMIN(aabb->lower[1], vertex[1]);
+ aabb->lower[2] = MMIN(aabb->lower[2], vertex[2]);
+ aabb->upper[0] = MMAX(aabb->upper[0], vertex[0]);
+ aabb->upper[1] = MMAX(aabb->upper[1], vertex[1]);
+ aabb->upper[2] = MMAX(aabb->upper[2], vertex[2]);
+ }
+}
+
+static void
+setup_sshape(FILE* stream, struct aabb* aabb)
+{
+ struct s3dut_mesh* sshape = NULL;
+ struct s3dut_mesh_data sshape_data;
+ struct s3dut_super_formula f0 = S3DUT_SUPER_FORMULA_NULL;
+ struct s3dut_super_formula f1 = S3DUT_SUPER_FORMULA_NULL;
+ const double radius = 2;
+ const unsigned nslices = 256;
+
+ f0.A = 1.5; f0.B = 1; f0.M = 11.0; f0.N0 = 1; f0.N1 = 1; f0.N2 = 2.0;
+ f1.A = 1.0; f1.B = 2; f1.M = 3.6; f1.N0 = 1; f1.N1 = 2; f1.N2 = 0.7;
+ S3DUT(create_super_shape(NULL, &f0, &f1, radius, nslices, nslices/2, &sshape));
+ S3DUT(mesh_get_data(sshape, &sshape_data));
+
+ aabb_update(aabb, &sshape_data);
+ sadist_write_stl(stream, sshape_data.positions, sshape_data.nvertices,
+ sshape_data.indices, sshape_data.nprimitives);
+
+ S3DUT(mesh_ref_put(sshape));
+}
+
+static void
+setup_sphere(FILE* stream, struct aabb* aabb)
+{
+ struct s3dut_mesh* sphere = NULL;
+ struct s3dut_mesh_data sphere_data;
+ const double radius = 1;
+ const unsigned nslices = 128;
+
+ S3DUT(create_sphere(NULL, radius, nslices, nslices/2, &sphere));
+ S3DUT(mesh_get_data(sphere, &sphere_data));
+
+ aabb_update(aabb, &sphere_data);
+ sadist_write_stl(stream, sphere_data.positions, sphere_data.nvertices,
+ sphere_data.indices, sphere_data.nprimitives);
+
+ S3DUT(mesh_ref_put(sphere));
+}
+
+static void
+setup_scene
+ (FILE* fp,
+ const char* sshape,
+ const char* sphere,
+ const struct aabb* aabb,
+ struct sadist_trilinear_profile* profile)
+{
+ double low, upp;
+ ASSERT(sshape && sphere && aabb && profile);
+
+ fprintf(fp, "PROGRAM solid_sphere libsadist_solid_sphere.so\n");
+ fprintf(fp, "PROGRAM trilinear_profile libsadist_trilinear_profile.so\n");
+ fprintf(fp, "SOLID SuperShape 25 7500 500 0.05 0 UNKNOWN 0 BACK %s FRONT %s\n",
+ sshape, sphere);
+ fprintf(fp, "SOLID_PROG Sphere solid_sphere BACK %s ", sphere);
+
+ /* StL lambda rho cp radius */
+ fprintf(fp, "PROG_PARAMS -m %s -l25 -r7500 -c500 -R1\n", FILENAME_SPHERE);
+
+ low = MMIN(MMIN(aabb->lower[0], aabb->lower[1]), aabb->lower[2]);
+ upp = MMAX(MMAX(aabb->upper[0], aabb->upper[1]), aabb->upper[2]);
+ fprintf(fp, "T_BOUNDARY_FOR_SOLID_PROG Dirichlet trilinear_profile %s", sshape);
+ fprintf(fp, " PROG_PARAMS -b %g,%g -t %g,%g,%g\n", low, upp, TX, TY, TZ);
+
+ d3_splat(profile->lower, low);
+ d3_splat(profile->upper, upp);
+ d2(profile->a, 0, TX);
+ d2(profile->b, 0, TY);
+ d2(profile->c, 0, TZ);
+}
+
+static int
+init(struct sadist_trilinear_profile* profile)
+{
+ struct aabb aabb = AABB_NULL;
+ FILE* fp_sshape = NULL;
+ FILE* fp_sphere = NULL;
+ FILE* fp_scene = NULL;
+ int err = 0;
+
+ #define FOPEN(Fp, Filename) \
+ if(((Fp) = fopen((Filename), "w")) == NULL) { \
+ fprintf(stderr, "Error opening `%s' -- file %s\n", \
+ (Filename), strerror(errno)); \
+ err = errno; \
+ goto error; \
+ } (void)0
+ FOPEN(fp_sshape, FILENAME_SSHAPE);
+ FOPEN(fp_sphere, FILENAME_SPHERE);
+ FOPEN(fp_scene, FILENAME_SCENE);
+ #undef FOPEN
+
+ setup_sshape(fp_sshape, &aabb);
+ setup_sphere(fp_sphere, &aabb);
+ setup_scene(fp_scene, FILENAME_SSHAPE, FILENAME_SPHERE, &aabb, profile);
+
+exit:
+ #define FCLOSE(Fp) \
+ if((Fp) && fclose(Fp)) { perror("fclose"); if(!err) err = errno; } (void)0
+ FCLOSE(fp_sshape);
+ FCLOSE(fp_sphere);
+ FCLOSE(fp_scene);
+ #undef FCLOSE
+ return err;
+error:
+ goto exit;
+}
+
+static int
+run(const struct sadist_trilinear_profile* profile)
+{
+ const double P[3] = {PX, PY, PZ};
+
+ FILE* output = NULL;
+ double E = 0;
+ double SE = 0;
+ double ref = 0;
+
+ int n = 0;
+ int err = 0;
+ int status = 0;
+
+ printf("%s\n", COMMAND);
+
+ if(!(output = popen(COMMAND, "r"))) {
+ fprintf(stderr, "Error executing stardis -- %s\n", strerror(errno));
+ err = errno;
+ goto error;
+ }
+
+ if((n = fscanf(output, "%lf %lf %*d %*d", &E, &SE)), n != 2 && n != EOF) {
+ fprintf(stderr, "Error reading the output stream -- %s\n", strerror(errno));
+ err = errno;
+ goto error;
+ }
+
+ /* Check command exit status */
+ if((status=pclose(output)), output=NULL, status) {
+ if(status == -1) err = errno;
+ else if(WIFEXITED(status)) err = WEXITSTATUS(status);
+ else if(WIFSIGNALED(status)) err = WTERMSIG(status);
+ else if(WIFSTOPPED(status)) err = WSTOPSIG(status);
+ else FATAL("Unreachable code.\n");
+ goto error;
+ }
+
+ ref = sadist_trilinear_profile_temperature(profile, P);
+ printf("T = %g ~ %g +/- %g\n", ref, E, SE);
+ if(!eq_eps(ref, E, SE*3)) {
+ err = 1;
+ goto error;
+ }
+
+exit:
+ if(output) pclose(output);
+ return err;
+error:
+ goto exit;
+}
+
+/*******************************************************************************
+ * The test
+ ******************************************************************************/
+int
+main(int argc, char** argv)
+{
+ struct sadist_trilinear_profile profile = SADIST_TRILINEAR_PROFILE_NULL;
+ int err;
+ (void)argc, (void)argv;
+
+ if((err = init(&profile))) goto error;
+ if((err = run(&profile))) goto error;
+
+exit:
+ return err;
+error:
+ goto exit;
+}
diff --git a/src/sadist_lib_solid_sphere.c b/src/sadist_lib_solid_sphere.c
@@ -0,0 +1,397 @@
+/* Copyright (C) 2024 |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 <stardis/stardis-prog-properties.h>
+
+#include <star/s3d.h>
+#include <star/ssp.h>
+#include <star/sstl.h>
+
+#include <rsys/cstr.h>
+#include <rsys/double3.h>
+#include <rsys/mem_allocator.h>
+
+#include <getopt.h>
+
+struct args {
+ const char* stl;
+ double lambda; /* [W/m/K] */
+ double rho; /* [kg/m^3] */
+ double cp; /* [J/K/kg] */
+ double radius; /* [m/fp_to_meter] */
+};
+#define ARGS_DEFAULT__ {NULL,1,1,1,1}
+static const struct args ARGS_DEFAULT = ARGS_DEFAULT__;
+
+/* Define a solid sphere */
+struct solid_sphere {
+ struct s3d_scene_view* view;
+ double lambda; /* [W/m/K] */
+ double rho; /* [kg/m^3] */
+ double cp; /* [J/K/kg] */
+ double radius; /* [m/fp_to_meter] */
+};
+#define SOLID_SPHERE_NULL__ {NULL,0,0,0,0}
+static const struct solid_sphere SOLID_SPHERE_NULL = SOLID_SPHERE_NULL__;
+
+/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+static void
+get_indices(const unsigned itri, unsigned ids[3], void* ctx)
+{
+ const struct sstl_desc* desc = ctx;
+ ASSERT(desc && itri < desc->triangles_count);
+ ids[0] = desc->indices[itri*3+0];
+ ids[1] = desc->indices[itri*3+1];
+ ids[2] = desc->indices[itri*3+2];
+}
+
+static void
+get_position(const unsigned ivert, float pos[3], void* ctx)
+{
+ const struct sstl_desc* desc = ctx;
+ ASSERT(desc && ivert < desc->vertices_count);
+ pos[0] = desc->vertices[ivert*3+0];
+ pos[1] = desc->vertices[ivert*3+1];
+ pos[2] = desc->vertices[ivert*3+2];
+}
+
+static struct s3d_scene_view*
+create_scene_view(const char* filename)
+{
+ /* Star-STL */
+ struct sstl_desc desc;
+ struct sstl* sstl = NULL;
+
+ /* Star3D */
+ struct s3d_vertex_data vdata = S3D_VERTEX_DATA_NULL;
+ struct s3d_device* dev = NULL;
+ struct s3d_scene* scn = NULL;
+ struct s3d_shape* shape = NULL;
+ struct s3d_scene_view* view = NULL;
+
+ res_T res = RES_OK;
+ ASSERT(filename);
+
+ if((res = sstl_create(NULL, NULL, 1, &sstl)) != RES_OK) goto error;
+ if((res = sstl_load(sstl, filename)) != RES_OK) goto error;
+ if((res = sstl_get_desc(sstl, &desc)) != RES_OK) goto error;
+
+ if((res = s3d_device_create(NULL, NULL, 0, &dev)) != RES_OK) goto error;
+ if((res = s3d_scene_create(dev, &scn)) != RES_OK) goto error;
+ if((res = s3d_shape_create_mesh(dev, &shape)) != RES_OK) goto error;
+ if((res = s3d_scene_attach_shape(scn, shape)) != RES_OK) goto error;
+
+ vdata.usage = S3D_POSITION;
+ vdata.type = S3D_FLOAT3;
+ vdata.get = get_position;
+ res = s3d_mesh_setup_indexed_vertices(shape,
+ (unsigned int)desc.triangles_count, get_indices,
+ (unsigned int)desc.vertices_count, &vdata, 1, &desc);
+ if(res != RES_OK) goto error;
+
+ res = s3d_scene_view_create(scn, S3D_TRACE|S3D_GET_PRIMITIVE, &view);
+ if(res != RES_OK) goto error;
+
+exit:
+ if(sstl) SSTL(ref_put(sstl));
+ if(dev) S3D(device_ref_put(dev));
+ if(scn) S3D(scene_ref_put(scn));
+ if(shape) S3D(shape_ref_put(shape));
+ return view;
+error:
+ if(view) { S3D(scene_view_ref_put(view)); view = NULL; }
+ goto exit;
+}
+
+static void
+print_usage(FILE* stream, const char* name)
+{
+ ASSERT(name);
+ fprintf(stream,
+ "usage: %s -m mesh [-c capa] [-h] [-l lambda] [-R radius] [-r rho]\n",
+ name);
+}
+
+static res_T
+parse_args
+ (const struct stardis_description_create_context* ctx,
+ struct args* args,
+ const int argc,
+ char* argv[])
+{
+ int opt = 0;
+ res_T res = RES_OK;
+
+ /* Check pre-conditions */
+ ASSERT(ctx && args);
+ *args = ARGS_DEFAULT;
+
+ optind = 1;
+ while((opt = getopt(argc, argv, "c:hl:m:R:r:")) != -1) {
+ switch(opt) {
+ case 'c':
+ res = cstr_to_double(optarg, &args->cp);
+ if(res == RES_OK && args->cp <= 0) res = RES_BAD_ARG;
+ break;
+ case 'h':
+ print_usage(stdout, argv[0]);
+ break;
+ case 'l':
+ res = cstr_to_double(optarg, &args->lambda);
+ if(res == RES_OK && args->lambda <= 0) res = RES_BAD_ARG;
+ break;
+ case 'm':
+ args->stl = optarg;
+ break;
+ case 'R':
+ res = cstr_to_double(optarg, &args->radius);
+ if(res == RES_OK && args->radius <= 0) res = RES_BAD_ARG;
+ break;
+ case 'r':
+ res = cstr_to_double(optarg, &args->rho);
+ if(res == RES_OK && args->rho <= 0) res = RES_BAD_ARG;
+ break;
+ default: res = RES_BAD_ARG; break;
+ }
+ if(res != RES_OK) {
+ if(optarg) {
+ fprintf(stderr, "%s: invalid option argument '%s' -- '%c'\n",
+ ctx->name, optarg, opt);
+ }
+ goto error;
+ }
+ }
+
+ if(!args->stl) {
+ fprintf(stderr, "%s: missing mesh -- option '-m'\n", argv[0]);
+ res = RES_BAD_ARG;
+ goto error;
+ }
+
+exit:
+ return res;
+error:
+ goto exit;
+}
+
+/*******************************************************************************
+ * Create data
+ ******************************************************************************/
+void*
+stardis_create_data
+ (const struct stardis_description_create_context *ctx,
+ void* libdata,
+ size_t argc,
+ char* argv[])
+{
+ struct args args = ARGS_DEFAULT__;
+ struct solid_sphere* solid = NULL;
+ res_T res = RES_OK;
+ (void)libdata;
+
+ solid = mem_alloc(sizeof(*solid));
+ if(!solid) {
+ fprintf(stderr, "%s:%d: error allocating the solid sphere.\n",
+ __FILE__, __LINE__);
+ goto error;
+ }
+ *solid = SOLID_SPHERE_NULL;
+
+ res = parse_args(ctx, &args, (int)argc, argv);
+ if(res != RES_OK) goto error;
+
+ solid->lambda = args.lambda;
+ solid->rho = args.rho;
+ solid->cp = args.cp;
+ solid->radius = args.radius;
+
+ if(!(solid->view = create_scene_view(args.stl))) goto error;
+
+exit:
+ return solid;
+error:
+ if(solid) {
+ stardis_release_data(solid);
+ solid = NULL;
+ }
+ goto exit;
+}
+
+void
+stardis_release_data(void* data)
+{
+ struct solid_sphere* solid = data;
+ ASSERT(solid);
+ if(solid->view) S3D(scene_view_ref_put(solid->view));
+ mem_rm(solid);
+}
+
+/*******************************************************************************
+ * Solid properties
+ ******************************************************************************/
+double
+stardis_calorific_capacity(const struct stardis_vertex* vtx, void* data)
+{
+ const struct solid_sphere* sphere = data;
+ (void)vtx;
+ return sphere->cp;
+}
+
+double
+stardis_volumic_mass(const struct stardis_vertex* vtx, void* data)
+{
+ const struct solid_sphere* sphere = data;
+ (void)vtx;
+ return sphere->rho;
+}
+
+double
+stardis_conductivity(const struct stardis_vertex* vtx, void* data)
+{
+ const struct solid_sphere* sphere = data;
+ (void)vtx;
+ return sphere->lambda;
+}
+
+double
+stardis_delta_solid(const struct stardis_vertex* vtx, void* data)
+{
+ const struct solid_sphere* sphere = data;
+ (void)vtx;
+ return sphere->radius/20.0;
+}
+
+double
+stardis_volumic_power(const struct stardis_vertex* vtx, void* data)
+{
+ (void)vtx, (void)data;
+ return STARDIS_VOLUMIC_POWER_NONE;
+}
+
+double
+stardis_medium_temperature(const struct stardis_vertex* vtx, void* data)
+{
+ (void)vtx, (void)data;
+ return STARDIS_TEMPERATURE_NONE;
+}
+
+int
+stardis_sample_conductive_path
+ (struct sdis_scene* scn,
+ struct ssp_rng* rng,
+ struct stardis_path* path,
+ void* data)
+{
+ const struct solid_sphere* sphere = data;
+ struct s3d_hit hit = S3D_HIT_NULL;
+ struct s3d_attrib attr0, attr1, attr2;
+
+ double pos[3];
+ double epsilon = 0;
+ ASSERT(scn && rng && path && data);
+ (void)scn;
+
+ epsilon = sphere->radius * 1.e-6; /* Epsilon shell */
+
+ d3_set(pos, path->vtx.P);
+
+ do {
+ /* Distance from the geometry center to the current position */
+ const double dst = d3_len(pos);
+
+ double dir[3] = {0,0,0};
+ double r = 0; /* Radius */
+
+ r = sphere->radius - dst;
+ CHK(dst > 0);
+
+ if(r > epsilon) {
+ /* Uniformly sample a new position on the surrounding sphere */
+ ssp_ran_sphere_uniform(rng, dir, NULL);
+
+ /* Move to the new position */
+ d3_muld(dir, dir, r);
+ d3_add(pos, pos, dir);
+
+ /* The current position is in the epsilon shell:
+ * move it to the nearest interface position */
+ } else {
+ float posf[3];
+
+ d3_set(dir, pos);
+ d3_normalize(dir, dir);
+ d3_muld(pos, dir, sphere->radius);
+
+ /* Map the position to the sphere geometry */
+ f3_set_d3(posf, pos);
+ S3D(scene_view_closest_point(sphere->view, posf, (float)INF, NULL, &hit));
+ }
+
+ /* The calculation is performed in steady state, so the path necessarily stops
+ * at a boundary */
+ } while(S3D_HIT_NONE(&hit));
+
+ /* Setup the path state */
+ d3_set(path->vtx.P, pos);
+ path->weight = 0;
+ path->at_limit = 0;
+
+ /* Configure the intersecting primitive. The vertices must be exactly the same
+ * (i.e. bit-compatible) as those used internally by Stardis, as they are used
+ * as the key to retrieve the corresponding Stardis triangle. Note that
+ * Stardis also uses Star-STL internally to load its geometry. As a result,
+ * the in-memory representation of the loaded data is the same between Stardis
+ * and the current code, i.e. single-precision floating point. And Star-3D
+ * also uses the same representation. We can therefore use Star-3D data
+ * directly to define the coordinates of the vertices of the intersecting
+ * primitive */
+ S3D(triangle_get_vertex_attrib(&hit.prim, 0, S3D_POSITION, &attr0));
+ S3D(triangle_get_vertex_attrib(&hit.prim, 1, S3D_POSITION, &attr1));
+ S3D(triangle_get_vertex_attrib(&hit.prim, 2, S3D_POSITION, &attr2));
+ d3_set_f3(path->tri.vtx0, attr0.value);
+ d3_set_f3(path->tri.vtx1, attr1.value);
+ d3_set_f3(path->tri.vtx2, attr2.value);
+
+ return RES_OK;
+}
+
+/*******************************************************************************
+ * Legal notices
+ ******************************************************************************/
+const char*
+get_copyright_notice(void* data)
+{
+ (void)data; /* Avoid "unused variable" warnings */
+ return "Copyright (C) 2024 |Méso|Star> (contact@meso-star.com)";
+}
+
+const char*
+get_license_short(void* data)
+{
+ (void)data; /* Avoid "unused variable" warnings */
+ return "GNU GPL version 3 or later <http://www.gnu.org/licenses/>";
+}
+
+const char*
+get_license_text(void* data)
+{
+ (void)data; /* Avoid "unused variable" warnings */
+ return
+ "This is free software released under the GPL v3+ license: GNU GPL\n"
+ "version 3 or later. You are welcome to redistribute it under certain\n"
+ "conditions; refer to <http://www.gnu.org/licenses/> for details.";
+}
