stardis-solver

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

test_sdis_primkey_2d.c (10121B)


      1 /* Copyright (C) 2016-2025 |Méso|Star> (contact@meso-star.com)
      2  *
      3  * This program is free software: you can redistribute it and/or modify
      4  * it under the terms of the GNU General Public License as published by
      5  * the Free Software Foundation, either version 3 of the License, or
      6  * (at your option) any later version.
      7  *
      8  * This program is distributed in the hope that it will be useful,
      9  * but WITHOUT ANY WARRANTY; without even the implied warranty of
     10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
     11  * GNU General Public License for more details.
     12  *
     13  * You should have received a copy of the GNU General Public License
     14  * along with this program. If not, see <http://www.gnu.org/licenses/>. */
     15 
     16 #include "sdis.h"
     17 #include "test_sdis_utils.h"
     18 #include "test_sdis_mesh.h"
     19 
     20 #include <star/s2d.h>
     21 
     22 #include <rsys/float2.h>
     23 
     24 /*
     25  * The system is a solid supershape whose boundary temperature is set to a
     26  * constant. The temperature of the solid is therefore this same temperature.
     27  * This simplistic test case is not used to verify a Monte Carlo estimate, but
     28  * to ensure that the caller can recover the internal representation of the
     29  * geometric primitives from his own data.
     30  *
     31  *             /\
     32  *         ___/  \___
     33  *         \        /
     34  *         /_  __  _\
     35  *           \/  \/
     36  *
     37  */
     38 
     39 /*******************************************************************************
     40  * Super shape
     41  ******************************************************************************/
     42 static struct mesh
     43 create_super_shape(void)
     44 {
     45   struct mesh sshape = MESH_NULL;
     46 
     47   const unsigned nslices = 128;
     48   const double a = 1.0;
     49   const double b = 1.0;
     50   const double n1 = 1.0;
     51   const double n2 = 1.0;
     52   const double n3 = 1.0;
     53   const double m = 6.0;
     54   size_t i = 0;
     55 
     56   FOR_EACH(i, 0, nslices) {
     57     const double theta = (double)i * (2.0*PI / (double)nslices);
     58     const double tmp0 = pow(fabs(1.0/a * cos(m/4.0*theta)), n2);
     59     const double tmp1 = pow(fabs(1.0/b * sin(m/4.0*theta)), n3);
     60     const double tmp2 = pow(tmp0 + tmp1, 1.0/n1);
     61     const double r = 1.0 / tmp2;
     62     const double x = cos(theta) * r;
     63     const double y = sin(theta) * r;
     64     const size_t j = (i + 1) % nslices;
     65 
     66     sa_push(sshape.positions, x);
     67     sa_push(sshape.positions, y);
     68     sa_push(sshape.indices, i);
     69     sa_push(sshape.indices, j);
     70   }
     71 
     72   return sshape;
     73 }
     74 
     75 /*******************************************************************************
     76  * Solid, i.e. medium of the super shape
     77  ******************************************************************************/
     78 #define SOLID_PROP(Prop, Val)                                                  \
     79   static double                                                                \
     80   solid_get_##Prop                                                             \
     81     (const struct sdis_rwalk_vertex* vtx,                                      \
     82      struct sdis_data* data)                                                   \
     83   {                                                                            \
     84     (void)vtx, (void)data; /* Avoid the "unused variable" warning */           \
     85     return Val;                                                                \
     86   }
     87 SOLID_PROP(calorific_capacity, 1) /* [J/K/Kg] */
     88 SOLID_PROP(thermal_conductivity, 1) /* [W/m/K] */
     89 SOLID_PROP(volumic_mass, 1) /* [kg/m^3] */
     90 SOLID_PROP(temperature, SDIS_TEMPERATURE_NONE) /* [K] */
     91 SOLID_PROP(delta, 1.0/20.0) /* [m/fp_to_meter] */
     92 #undef SOLID_PROP
     93 
     94 static struct sdis_medium*
     95 create_solid(struct sdis_device* sdis)
     96 {
     97   struct sdis_solid_shader shader = SDIS_SOLID_SHADER_NULL;
     98   struct sdis_medium* solid = NULL;
     99 
    100   shader.calorific_capacity = solid_get_calorific_capacity;
    101   shader.thermal_conductivity = solid_get_thermal_conductivity;
    102   shader.volumic_mass = solid_get_volumic_mass;
    103   shader.delta = solid_get_delta;
    104   shader.temperature = solid_get_temperature;
    105   OK(sdis_solid_create(sdis, &shader, NULL, &solid));
    106   return solid;
    107 }
    108 
    109 /*******************************************************************************
    110  * Dummy environment, i.e. environment surrounding the super shape. It is
    111  * defined only for Stardis compliance: in Stardis, an interface must divide 2
    112  * media.
    113  ******************************************************************************/
    114 static struct sdis_medium*
    115 create_dummy(struct sdis_device* sdis)
    116 {
    117   struct sdis_fluid_shader shader = SDIS_FLUID_SHADER_NULL;
    118   struct sdis_medium* dummy = NULL;
    119 
    120   shader.calorific_capacity = dummy_medium_getter;
    121   shader.volumic_mass = dummy_medium_getter;
    122   shader.temperature = dummy_medium_getter;
    123   OK(sdis_fluid_create(sdis, &shader, NULL, &dummy));
    124   return dummy;
    125 }
    126 
    127 /*******************************************************************************
    128  * Interface
    129  ******************************************************************************/
    130 static double
    131 interface_get_temperature
    132   (const struct sdis_interface_fragment* frag,
    133    struct sdis_data* data)
    134 {
    135   (void)frag, (void)data; /* Avoid the "unused variable" warning */
    136   return 300; /* [K] */
    137 }
    138 
    139 static struct sdis_interface*
    140 create_interface
    141   (struct sdis_device* sdis,
    142    struct sdis_medium* front,
    143    struct sdis_medium* back)
    144 {
    145   struct sdis_interface* interf = NULL;
    146   struct sdis_interface_shader shader = SDIS_INTERFACE_SHADER_NULL;
    147 
    148   shader.front.temperature = interface_get_temperature;
    149   shader.back.temperature = interface_get_temperature;
    150   OK(sdis_interface_create(sdis, front, back, &shader, NULL, &interf));
    151   return interf;
    152 }
    153 
    154 /*******************************************************************************
    155  * Scene, i.e. the system to simulate
    156  ******************************************************************************/
    157 struct scene_context {
    158   const struct mesh* sshape;
    159   struct sdis_interface* interf;
    160 };
    161 static const struct scene_context SCENE_CONTEXT_NULL = {NULL, NULL};
    162 
    163 static void
    164 scene_get_indices(const size_t iseg, size_t ids[2], void* ctx)
    165 {
    166   struct scene_context* context = ctx;
    167   CHK(ids && context);
    168   /* Flip the indices to ensure that the normal points into the super shape */
    169   ids[0] = (unsigned)context->sshape->indices[iseg*2+1];
    170   ids[1] = (unsigned)context->sshape->indices[iseg*2+0];
    171 }
    172 
    173 static void
    174 scene_get_interface(const size_t iseg, struct sdis_interface** interf, void* ctx)
    175 {
    176   struct scene_context* context = ctx;
    177   (void)iseg;
    178   CHK(interf && context);
    179   *interf = context->interf;
    180 }
    181 
    182 static void
    183 scene_get_position(const size_t ivert, double pos[2], void* ctx)
    184 {
    185   struct scene_context* context = ctx;
    186   CHK(pos && context);
    187   pos[0] = context->sshape->positions[ivert*2+0];
    188   pos[1] = context->sshape->positions[ivert*2+1];
    189 }
    190 
    191 static struct sdis_scene*
    192 create_scene
    193   (struct sdis_device* sdis,
    194    const struct mesh* sshape,
    195    struct sdis_interface* interf)
    196 {
    197   struct sdis_scene* scn = NULL;
    198   struct sdis_scene_create_args scn_args = SDIS_SCENE_CREATE_ARGS_DEFAULT;
    199   struct scene_context context = SCENE_CONTEXT_NULL;
    200 
    201   context.interf = interf;
    202   context.sshape = sshape;
    203 
    204   scn_args.get_indices = scene_get_indices;
    205   scn_args.get_interface = scene_get_interface;
    206   scn_args.get_position = scene_get_position;
    207   scn_args.nprimitives = mesh_2d_nsegments(sshape);
    208   scn_args.nvertices = mesh_2d_nvertices(sshape);
    209   scn_args.context = &context;
    210   OK(sdis_scene_2d_create(sdis, &scn_args, &scn));
    211   return scn;
    212 }
    213 
    214 /*******************************************************************************
    215  * Validation
    216  ******************************************************************************/
    217 static void
    218 check(struct sdis_scene* scn, const struct mesh* mesh)
    219 {
    220   struct s2d_primitive prim = S2D_PRIMITIVE_NULL;
    221   struct sdis_primkey key = SDIS_PRIMKEY_NULL;
    222   size_t iprim = 0;
    223   size_t nprims = 0;
    224 
    225   BA(sdis_scene_get_s2d_primitive(NULL, &key, &prim));
    226   BA(sdis_scene_get_s2d_primitive(scn, NULL, &prim));
    227   BA(sdis_scene_get_s2d_primitive(scn, &key, NULL));
    228   BA(sdis_scene_get_s2d_primitive(scn, &key, &prim));
    229 
    230   nprims = mesh_2d_nsegments(mesh);
    231   FOR_EACH(iprim, 0, nprims) {
    232     const double *v0, *v1;
    233     float v0f[2], v1f[2];
    234     struct s2d_attrib attr0, attr1;
    235 
    236     /* Check that a primitive can be obtained from the key constructed on the
    237      * user side */
    238     v0 = mesh->positions + mesh->indices[iprim*2 + 0]*2;
    239     v1 = mesh->positions + mesh->indices[iprim*2 + 1]*2;
    240     sdis_primkey_2d_setup(&key, v0, v1);
    241     OK(sdis_scene_get_s2d_primitive(scn, &key, &prim));
    242 
    243     /* Check that the primitive on the solver side is the same as that on the
    244      * user side. On the solver side, vertices are stored in simple precision in
    245      * Star-3D view. We therefore need to take care of this conversion to check
    246      * that the vertices are the same */
    247     OK(s2d_segment_get_vertex_attrib(&prim, 0, S2D_POSITION, &attr0));
    248     OK(s2d_segment_get_vertex_attrib(&prim, 1, S2D_POSITION, &attr1));
    249     f2_set_d2(v0f, v0);
    250     f2_set_d2(v1f, v1);
    251 
    252     /* The vertices have been inverted on the user's side to reverse the normal
    253      * orientation. Below it is taken into account */
    254     CHK(f2_eq(v0f, attr1.value));
    255     CHK(f2_eq(v1f, attr0.value));
    256   }
    257 }
    258 
    259 /*******************************************************************************
    260  * The test
    261  ******************************************************************************/
    262 int
    263 main(int argc, char** argv)
    264 {
    265   /* Stardis */
    266   struct sdis_device* sdis = NULL;
    267   struct sdis_interface* interf = NULL;
    268   struct sdis_medium* solid = NULL;
    269   struct sdis_medium* dummy = NULL; /* Medium surrounding the solid */
    270   struct sdis_scene* scn = NULL;
    271 
    272   /* Miscellaneous */
    273   struct mesh sshape = MESH_NULL;
    274   (void)argc, (void)argv;
    275 
    276   OK(sdis_device_create(&SDIS_DEVICE_CREATE_ARGS_DEFAULT, &sdis));
    277 
    278   sshape = create_super_shape();
    279   solid = create_solid(sdis);
    280   dummy = create_dummy(sdis);
    281   interf = create_interface(sdis, solid, dummy);
    282   scn = create_scene(sdis, &sshape, interf);
    283 
    284   check(scn, &sshape);
    285 
    286   mesh_release(&sshape);
    287   OK(sdis_device_ref_put(sdis));
    288   OK(sdis_interface_ref_put(interf));
    289   OK(sdis_medium_ref_put(solid));
    290   OK(sdis_medium_ref_put(dummy));
    291   OK(sdis_scene_ref_put(scn));
    292 
    293   CHK(mem_allocated_size() == 0);
    294   return 0;
    295 }