stardis-solver

Solve coupled heat transfers
git clone git://git.meso-star.fr/stardis-solver.git
Log | Files | Refs | README | LICENSE

commit ecad4c7363faff375b940ba2cd0f36b8871f50b3
parent 79af4c75fe5d943450d28a880d5cb2270ea33a04
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Mon,  9 Apr 2018 15:17:32 +0200

Test the volumic power property of a solid in 2D and 3D

Diffstat:
Mcmake/CMakeLists.txt | 1+
Msrc/test_sdis_solve_probe_boundary.c | 8+++++---
Asrc/test_sdis_volumic_power.c | 367+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
3 files changed, 373 insertions(+), 3 deletions(-)

diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -135,6 +135,7 @@ if(NOT NO_TEST) new_test(test_sdis_solve_probe2_2d) new_test(test_sdis_solve_probe3_2d) new_test(test_sdis_solve_probe_boundary) + new_test(test_sdis_volumic_power) target_link_libraries(test_sdis_solve_probe3 Star3DUT) target_link_libraries(test_sdis_solve_camera Star3DUT) diff --git a/src/test_sdis_solve_probe_boundary.c b/src/test_sdis_solve_probe_boundary.c @@ -32,14 +32,14 @@ * * 3D 2D * - * (1,1,1) (1,1) + * ///// (1,1,1) ///// (1,1) * +-------+ +-------+ * /' /| _\ | | _\ * +-------+ | / / Tf Tb | / / Tf * Tb +.....|.+ \__/ | | \__/ * |, |/ +-------+ - * +-------+ (0,0) - * (0,0,0) + * +-------+ (0,0) ///// + * (0,0,0) ///// */ #define UNKNOWN_TEMPERATURE -1 @@ -365,6 +365,7 @@ main(int argc, char** argv) printf("Boundary temperature of the box at (%g %g %g) = %g ~ %g +/- %g\n", SPLIT3(pos), ref, T.E, T.SE); + printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(eq_eps(T.E, ref, T.SE)); uv[0] = 0.5; @@ -381,6 +382,7 @@ main(int argc, char** argv) printf("Boundary temperature of the square at (%g %g) = %g ~ %g +/- %g\n", SPLIT2(pos), ref, T.E, T.SE); + printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); CHK(eq_eps(T.E, ref, T.SE)); #undef SOLVE diff --git a/src/test_sdis_volumic_power.c b/src/test_sdis_volumic_power.c @@ -0,0 +1,367 @@ +/* Copyright (C) 2016-2018 |Meso|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.h" +#include "test_sdis_utils.h" + +#include <rsys/double3.h> +#include <rsys/math.h> + +/* + * The scene is composed of a solid cube/square whose temperature is unknown. + * The convection coefficient with the surrounding fluid is null everywhere The + * Temperature of the +/- X faces are fixed to T0, and the solid has a volumic + * power of P0. This test computes the temperature of a probe position pos into + * the solid and check that it is is equal to: + * + * T(pos) = P0 / (2*LAMBDA) * (A^2/4 - (pos-0.5)^2) + T0 + * + * with LAMBDA the conductivity of the cube and A the size of the cube/square, + * i.e. 1. + * + * 3D 2D + * + * ///// (1,1,1) ///// (1,1) + * +-------+ +-------+ + * /' /| | | + * +-------+ T0 T0 T0 + * T0 +.....|.+ | | + * |, |/ +-------+ + * +-------+ (0,0) ///// + * (0,0,0) ///// + */ + +#define UNKNOWN_TEMPERATURE -1 +#define N 10000 /* #realisations */ + +#define T0 320 +#define LAMBDA 0.1 +#define P0 10 + +/******************************************************************************* + * Geometry 3D + ******************************************************************************/ +static void +box_get_indices(const size_t itri, size_t ids[3], void* context) +{ + (void)context; + CHK(ids); + ids[0] = box_indices[itri*3+0]; + ids[1] = box_indices[itri*3+1]; + ids[2] = box_indices[itri*3+2]; +} + +static void +box_get_position(const size_t ivert, double pos[3], void* context) +{ + (void)context; + CHK(pos); + pos[0] = box_vertices[ivert*3+0]; + pos[1] = box_vertices[ivert*3+1]; + pos[2] = box_vertices[ivert*3+2]; +} + +static void +box_get_interface(const size_t itri, struct sdis_interface** bound, void* context) +{ + struct sdis_interface** interfaces = context; + CHK(context && bound); + *bound = interfaces[itri]; +} + +/******************************************************************************* + * Geometry 2D + ******************************************************************************/ +static void +square_get_indices(const size_t iseg, size_t ids[2], void* context) +{ + (void)context; + CHK(ids); + ids[0] = square_indices[iseg*2+0]; + ids[1] = square_indices[iseg*2+1]; +} + +static void +square_get_position(const size_t ivert, double pos[2], void* context) +{ + (void)context; + CHK(pos); + pos[0] = square_vertices[ivert*2+0]; + pos[1] = square_vertices[ivert*2+1]; +} + +static void +square_get_interface + (const size_t iseg, struct sdis_interface** bound, void* context) +{ + struct sdis_interface** interfaces = context; + CHK(context && bound); + *bound = interfaces[iseg]; +} + +/******************************************************************************* + * Media + ******************************************************************************/ +static double +fluid_get_temperature + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return UNKNOWN_TEMPERATURE; +} + +static double +solid_get_calorific_capacity + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return 2.0; +} + +static double +solid_get_thermal_conductivity + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return LAMBDA; +} + +static double +solid_get_volumic_mass + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return 25.0; +} + +static double +solid_get_delta + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return 1.0/20.0; +} + +static double +solid_get_delta_boundary + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return 2.1/20.0; +} + +static double +solid_get_temperature + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return UNKNOWN_TEMPERATURE; +} + +static double +solid_get_volumic_power + (const struct sdis_rwalk_vertex* vtx, struct sdis_data* data) +{ + (void)data; + CHK(vtx != NULL); + return P0; +} + +/******************************************************************************* + * Interfaces + ******************************************************************************/ +struct interf { + double temperature; +}; + +static double +interface_get_temperature + (const struct sdis_interface_fragment* frag, struct sdis_data* data) +{ + const struct interf* interf = sdis_data_cget(data); + CHK(frag && data); + return interf->temperature; +} + +static double +interface_get_convection_coef + (const struct sdis_interface_fragment* frag, struct sdis_data* data) +{ + CHK(frag && data); + return 0; +} + +static double +interface_get_emissivity + (const struct sdis_interface_fragment* frag, struct sdis_data* data) +{ + CHK(frag && data); + return 0; +} + +static double +interface_get_specular_fraction + (const struct sdis_interface_fragment* frag, struct sdis_data* data) +{ + CHK(frag && data); + return 0; +} + +/******************************************************************************* + * Test + ******************************************************************************/ +int +main(int argc, char** argv) +{ + struct mem_allocator allocator; + struct sdis_mc T = SDIS_MC_NULL; + struct sdis_data* data = NULL; + struct sdis_device* dev = NULL; + struct sdis_medium* fluid = NULL; + struct sdis_medium* solid = NULL; + struct sdis_interface* interf_adiabatic = NULL; + struct sdis_interface* interf_T0 = NULL; + struct sdis_scene* box_scn = NULL; + struct sdis_scene* square_scn = NULL; + struct sdis_estimator* estimator = NULL; + struct sdis_fluid_shader fluid_shader = DUMMY_FLUID_SHADER; + struct sdis_solid_shader solid_shader = DUMMY_SOLID_SHADER; + struct sdis_interface_shader interf_shader = DUMMY_INTERFACE_SHADER; + struct sdis_interface* box_interfaces[12 /*#triangles*/]; + struct sdis_interface* square_interfaces[4/*#segments*/]; + struct interf* interf_props = NULL; + double pos[3]; + double x; + double ref; + size_t nreals; + size_t nfails; + (void)argc, (void)argv; + + CHK(mem_init_proxy_allocator(&allocator, &mem_default_allocator) == RES_OK); + CHK(sdis_device_create + (NULL, &allocator, SDIS_NTHREADS_DEFAULT, 0, &dev) == RES_OK); + + /* Create the fluid medium */ + fluid_shader.temperature = fluid_get_temperature; + CHK(sdis_fluid_create(dev, &fluid_shader, NULL, &fluid) == RES_OK); + + /* Create the solid_medium */ + solid_shader.calorific_capacity = solid_get_calorific_capacity; + solid_shader.thermal_conductivity = solid_get_thermal_conductivity; + solid_shader.volumic_mass = solid_get_volumic_mass; + solid_shader.delta_solid = solid_get_delta; + solid_shader.delta_boundary = solid_get_delta_boundary; + solid_shader.temperature = solid_get_temperature; + solid_shader.volumic_power = solid_get_volumic_power; + CHK(sdis_solid_create(dev, &solid_shader, NULL, &solid) == RES_OK); + + /* Setup the interface shader */ + interf_shader.temperature = interface_get_temperature; + interf_shader.convection_coef = interface_get_convection_coef; + interf_shader.emissivity = interface_get_emissivity; + interf_shader.specular_fraction = interface_get_specular_fraction; + + /* Create the adiabatic interface */ + CHK(sdis_data_create(dev, sizeof(struct interf), 16, NULL, &data) == RES_OK); + interf_props = sdis_data_get(data); + interf_props->temperature = UNKNOWN_TEMPERATURE; + CHK(sdis_interface_create + (dev, solid, fluid, &interf_shader, data, &interf_adiabatic) == RES_OK); + CHK(sdis_data_ref_put(data) == RES_OK); + + /* Create the T0 interface */ + CHK(sdis_data_create(dev, sizeof(struct interf), 16, NULL, &data) == RES_OK); + interf_props = sdis_data_get(data); + interf_props->temperature = T0; + CHK(sdis_interface_create + (dev, solid, fluid, &interf_shader, data, &interf_T0) == RES_OK); + CHK(sdis_data_ref_put(data) == RES_OK); + + /* Release the media */ + CHK(sdis_medium_ref_put(solid) == RES_OK); + CHK(sdis_medium_ref_put(fluid) == RES_OK); + + /* Map the interfaces to their box triangles */ + box_interfaces[0] = box_interfaces[1] = interf_adiabatic; /* Front */ + box_interfaces[2] = box_interfaces[3] = interf_T0; /* Left */ + box_interfaces[4] = box_interfaces[5] = interf_adiabatic; /* Back */ + box_interfaces[6] = box_interfaces[7] = interf_T0; /* Right */ + box_interfaces[8] = box_interfaces[9] = interf_adiabatic; /* Top */ + box_interfaces[10]= box_interfaces[11]= interf_adiabatic; /* Bottom */ + + /* Map the interfaces to their square segments */ + square_interfaces[0] = interf_adiabatic; /* Bottom */ + square_interfaces[1] = interf_T0; /* Lef */ + square_interfaces[2] = interf_adiabatic; /* Top */ + square_interfaces[3] = interf_T0; /* Right */ + + /* Create the box scene */ + CHK(sdis_scene_create(dev, box_ntriangles, box_get_indices, + box_get_interface, box_nvertices, box_get_position, box_interfaces, + &box_scn) == RES_OK); + + /* Create the square scene */ + CHK(sdis_scene_2d_create(dev, square_nsegments, square_get_indices, + square_get_interface, square_nvertices, square_get_position, + square_interfaces, &square_scn) == RES_OK); + + /* Release the interfaces */ + CHK(sdis_interface_ref_put(interf_adiabatic) == RES_OK); + CHK(sdis_interface_ref_put(interf_T0) == RES_OK); + + d3_splat(pos, 0.25); + x = pos[0] - 0.5; + ref = P0 / (2*LAMBDA) * (1.0/4.0 - x*x) + T0; + + /* Solve in 3D */ + CHK(sdis_solve_probe(box_scn, N, pos, INF, 1.0, 0, 0, &estimator) == RES_OK); + CHK(sdis_estimator_get_realisation_count(estimator, &nreals) == RES_OK); + CHK(sdis_estimator_get_failure_count(estimator, &nfails) == RES_OK); + CHK(nfails + nreals == N); + CHK(sdis_estimator_get_temperature(estimator, &T) == RES_OK); + CHK(sdis_estimator_ref_put(estimator) == RES_OK); + printf("Temperature of the box at (%g %g %g) = %g ~ %g +/- %g\n", + SPLIT3(pos), ref, T.E, T.SE); + printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); + CHK(eq_eps(T.E, ref, T.SE*2)); + + /* Solve in 2D */ + CHK(sdis_solve_probe(square_scn, N, pos, INF, 1.0, 0, 0, &estimator) == RES_OK); + CHK(sdis_estimator_get_realisation_count(estimator, &nreals) == RES_OK); + CHK(sdis_estimator_get_failure_count(estimator, &nfails) == RES_OK); + CHK(nfails + nreals == N); + CHK(sdis_estimator_get_temperature(estimator, &T) == RES_OK); + CHK(sdis_estimator_ref_put(estimator) == RES_OK); + printf("Temperature of the square at (%g %g) = %g ~ %g +/- %g\n", + SPLIT2(pos), ref, T.E, T.SE); + printf("#failures = %lu/%lu\n", (unsigned long)nfails, (unsigned long)N); + CHK(eq_eps(T.E, ref, T.SE*2.0)); + + CHK(sdis_scene_ref_put(box_scn) == RES_OK); + CHK(sdis_scene_ref_put(square_scn) == RES_OK); + CHK(sdis_device_ref_put(dev) == RES_OK); + + check_memory_allocator(&allocator); + mem_shutdown_proxy_allocator(&allocator); + CHK(mem_allocated_size() == 0); + return 0; +}