test_sdis_utils.h (13096B)
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 #ifndef TEST_SDIS_UTILS_H 17 #define TEST_SDIS_UTILS_H 18 19 #include "sdis.h" 20 21 #include <rsys/double33.h> 22 #include <rsys/mem_allocator.h> 23 #include <stdio.h> 24 #include <string.h> 25 26 #ifdef SDIS_ENABLE_MPI 27 #include <mpi.h> 28 #endif 29 30 #define BOLTZMANN_CONSTANT 5.6696e-8 /* W/m^2/K^4 */ 31 32 #define OK(Cond) CHK((Cond) == RES_OK) 33 #define BA(Cond) CHK((Cond) == RES_BAD_ARG) 34 #define BO(Cond) CHK((Cond) == RES_BAD_OP) 35 36 /******************************************************************************* 37 * Box geometry 38 ******************************************************************************/ 39 static const double box_vertices[8/*#vertices*/*3/*#coords per vertex*/] = { 40 0.0, 0.0, 0.0, 41 1.0, 0.0, 0.0, 42 0.0, 1.0, 0.0, 43 1.0, 1.0, 0.0, 44 0.0, 0.0, 1.0, 45 1.0, 0.0, 1.0, 46 0.0, 1.0, 1.0, 47 1.0, 1.0, 1.0 48 }; 49 static const size_t box_nvertices = sizeof(box_vertices) / (sizeof(double)*3); 50 51 /* The following array lists the indices toward the 3D vertices of each 52 * triangle. 53 * ,2---,3 ,2----3 54 * ,' | ,'/| ,'/| \ | 55 * 6----7' / | 6' / | \ | Y 56 * |', | / ,1 | / ,0---,1 | 57 * | ',|/,' |/,' | ,' o--X 58 * 4----5' 4----5' / 59 * Front, right Back, left and Z 60 * and Top faces bottom faces */ 61 static const size_t box_indices[12/*#triangles*/*3/*#indices per triangle*/] = { 62 0, 2, 1, 1, 2, 3, /* -Z */ 63 0, 4, 2, 2, 4, 6, /* -X */ 64 4, 5, 6, 6, 5, 7, /* +Z */ 65 3, 7, 5, 5, 1, 3, /* +X */ 66 2, 6, 7, 7, 3, 2, /* +Y */ 67 0, 1, 5, 5, 4, 0 /* -Y */ 68 }; 69 static const size_t box_ntriangles = sizeof(box_indices) / (sizeof(size_t)*3); 70 71 static INLINE void 72 box_get_indices(const size_t itri, size_t ids[3], void* context) 73 { 74 (void)context; 75 CHK(ids); 76 CHK(itri < box_ntriangles); 77 ids[0] = box_indices[itri*3+0]; 78 ids[1] = box_indices[itri*3+1]; 79 ids[2] = box_indices[itri*3+2]; 80 } 81 82 static INLINE void 83 box_get_position(const size_t ivert, double pos[3], void* context) 84 { 85 (void)context; 86 CHK(pos); 87 CHK(ivert < box_nvertices); 88 pos[0] = box_vertices[ivert*3+0]; 89 pos[1] = box_vertices[ivert*3+1]; 90 pos[2] = box_vertices[ivert*3+2]; 91 } 92 93 static INLINE void 94 box_get_interface(const size_t itri, struct sdis_interface** bound, void* context) 95 { 96 struct sdis_interface** interfaces = context; 97 CHK(context && bound); 98 CHK(itri < box_ntriangles); 99 *bound = interfaces[itri]; 100 } 101 102 /******************************************************************************* 103 * Square geometry 104 ******************************************************************************/ 105 static const double square_vertices[4/*#vertices*/*2/*#coords per vertex*/] = { 106 1.0, 0.0, 107 0.0, 0.0, 108 0.0, 1.0, 109 1.0, 1.0 110 }; 111 static const size_t square_nvertices = sizeof(square_vertices)/(sizeof(double)*2); 112 113 static const size_t square_indices[4/*#segments*/*2/*#indices per segment*/]= { 114 0, 1, /* Bottom */ 115 1, 2, /* Left */ 116 2, 3, /* Top */ 117 3, 0 /* Right */ 118 }; 119 static const size_t square_nsegments = sizeof(square_indices)/(sizeof(size_t)*2); 120 121 static INLINE void 122 square_get_indices(const size_t iseg, size_t ids[2], void* context) 123 { 124 (void)context; 125 CHK(ids); 126 CHK(iseg < square_nsegments); 127 ids[0] = square_indices[iseg*2+0]; 128 ids[1] = square_indices[iseg*2+1]; 129 } 130 131 static INLINE void 132 square_get_position(const size_t ivert, double pos[2], void* context) 133 { 134 (void)context; 135 CHK(pos); 136 CHK(ivert < square_nvertices); 137 pos[0] = square_vertices[ivert*2+0]; 138 pos[1] = square_vertices[ivert*2+1]; 139 } 140 141 static INLINE void 142 square_get_interface 143 (const size_t iseg, struct sdis_interface** bound, void* context) 144 { 145 struct sdis_interface** interfaces = context; 146 CHK(context && bound); 147 CHK(iseg < square_nsegments); 148 *bound = interfaces[iseg]; 149 } 150 151 /******************************************************************************* 152 * Medium, interface and ray 153 ******************************************************************************/ 154 static INLINE double 155 dummy_medium_getter 156 (const struct sdis_rwalk_vertex* vert, struct sdis_data* data) 157 { 158 (void)data; 159 CHK(vert != NULL); 160 return 0; 161 } 162 163 static INLINE double 164 dummy_interface_getter 165 (const struct sdis_interface_fragment* frag, struct sdis_data* data) 166 { 167 (void)data; 168 CHK(frag != NULL); 169 return 0; 170 } 171 172 static INLINE double 173 dummy_radiative_interface_getter 174 (const struct sdis_interface_fragment* frag, 175 const unsigned source_id, 176 struct sdis_data* data) 177 { 178 (void)data, (void)source_id; 179 CHK(frag != NULL); 180 return 0; 181 } 182 183 static INLINE double 184 dummy_ray_getter(const struct sdis_radiative_ray* ray, struct sdis_data* data) 185 { 186 (void)data; 187 CHK(ray != NULL); 188 return 0; 189 } 190 191 static const struct sdis_solid_shader DUMMY_SOLID_SHADER = { 192 dummy_medium_getter, /* Calorific capacity */ 193 dummy_medium_getter, /* Thermal conductivity */ 194 dummy_medium_getter, /* Volumic mass */ 195 dummy_medium_getter, /* Delta */ 196 dummy_medium_getter, /* Volumic power */ 197 dummy_medium_getter, /* Temperature */ 198 NULL, /* sample path */ 199 0 /* Initial time */ 200 }; 201 202 static const struct sdis_fluid_shader DUMMY_FLUID_SHADER = { 203 dummy_medium_getter, /* Calorific capacity */ 204 dummy_medium_getter, /* Volumic mass */ 205 dummy_medium_getter, /* Temperature */ 206 0 /* Initial time */ 207 }; 208 209 #define DUMMY_INTERFACE_SIDE_SHADER__ { \ 210 dummy_interface_getter, /* Temperature */ \ 211 dummy_interface_getter, /* Flux */ \ 212 dummy_radiative_interface_getter, /* Emissivity */ \ 213 dummy_radiative_interface_getter, /* Specular fraction */ \ 214 dummy_interface_getter, /* Reference temperature */ \ 215 1 /* Handle external flux */ \ 216 } 217 static const struct sdis_interface_shader DUMMY_INTERFACE_SHADER = { 218 dummy_interface_getter, /* Convection coef */ 219 0, /* Upper bound of the convection coef */ 220 dummy_interface_getter, /* Thermal contact resistance */ 221 DUMMY_INTERFACE_SIDE_SHADER__, /* Front side */ 222 DUMMY_INTERFACE_SIDE_SHADER__ /* Back side */ 223 }; 224 225 static const struct sdis_radiative_env_shader DUMMY_RAY_SHADER = { 226 dummy_ray_getter, 227 dummy_ray_getter 228 }; 229 230 /******************************************************************************* 231 * Device creation 232 ******************************************************************************/ 233 #ifndef SDIS_ENABLE_MPI 234 235 static INLINE void 236 create_default_device 237 (int* argc, 238 char*** argv, 239 int* is_master_process, 240 struct sdis_device** dev) 241 { 242 (void)argc, (void)argv; 243 CHK(dev && is_master_process); 244 OK(sdis_device_create(&SDIS_DEVICE_CREATE_ARGS_DEFAULT, dev)); 245 *is_master_process = 1; 246 } 247 248 #else 249 250 static INLINE void 251 create_default_device 252 (int* pargc, 253 char*** pargv, 254 int* is_master_process, 255 struct sdis_device** out_dev) 256 { 257 struct sdis_device_create_args dev_args = SDIS_DEVICE_CREATE_ARGS_DEFAULT; 258 struct sdis_device* dev = NULL; 259 int mpi_thread_support; 260 int mpi_rank; 261 CHK(pargc && pargv && is_master_process && out_dev); 262 263 CHK(MPI_Init_thread 264 (pargc, pargv, MPI_THREAD_SERIALIZED, &mpi_thread_support) == MPI_SUCCESS); 265 CHK(mpi_thread_support >= MPI_THREAD_SERIALIZED); 266 267 dev_args.use_mpi = *pargc >= 2 && !strcmp((*pargv)[1], "mpi"); 268 OK(sdis_device_create(&dev_args, &dev)); 269 270 if(dev_args.use_mpi) { 271 OK(sdis_device_get_mpi_rank(dev, &mpi_rank)); 272 *is_master_process = mpi_rank == 0; 273 } else { 274 CHK(sdis_device_get_mpi_rank(dev, &mpi_rank) == RES_BAD_OP); 275 *is_master_process = 1; 276 } 277 *out_dev = dev; 278 } 279 #endif 280 281 static INLINE void 282 free_default_device(struct sdis_device* dev) 283 { 284 OK(sdis_device_ref_put(dev)); 285 #ifdef SDIS_ENABLE_MPI 286 CHK(MPI_Finalize() == MPI_SUCCESS); 287 #endif 288 } 289 290 /******************************************************************************* 291 * Miscellaneous 292 ******************************************************************************/ 293 static INLINE void 294 dump_mesh 295 (FILE* stream, 296 const double* pos, 297 const size_t npos, 298 const size_t* ids, 299 const size_t nids) 300 { 301 size_t i; 302 CHK(pos != NULL && npos != 0); 303 CHK(ids != NULL && nids != 0); 304 FOR_EACH(i, 0, npos) { 305 fprintf(stream, "v %g %g %g\n", SPLIT3(pos+i*3)); 306 } 307 FOR_EACH(i, 0, nids) { 308 fprintf(stream, "f %lu %lu %lu\n", 309 (unsigned long)(ids[i*3+0] + 1), 310 (unsigned long)(ids[i*3+1] + 1), 311 (unsigned long)(ids[i*3+2] + 1)); 312 } 313 fflush(stream); 314 } 315 316 static INLINE void 317 dump_segments 318 (FILE* stream, 319 const double* pos, 320 const size_t npos, 321 const size_t* ids, 322 const size_t nids) 323 { 324 size_t i; 325 CHK(pos != NULL && npos != 0); 326 CHK(ids != NULL && nids != 0); 327 FOR_EACH(i, 0, npos) { 328 fprintf(stream, "v %g %g 0\n", SPLIT2(pos+i*2)); 329 } 330 FOR_EACH(i, 0, nids) { 331 fprintf(stream, "l %lu %lu\n", 332 (unsigned long)(ids[i*2+0] + 1), 333 (unsigned long)(ids[i*2+1] + 1)); 334 } 335 fflush(stream); 336 } 337 338 static INLINE void 339 check_estimator_eq 340 (const struct sdis_estimator* e1, const struct sdis_estimator* e2) 341 { 342 struct sdis_mc mc1, mc2; 343 enum sdis_estimator_type type1, type2; 344 ASSERT(e1 && e2); 345 346 OK(sdis_estimator_get_type(e1, &type1)); 347 OK(sdis_estimator_get_type(e2, &type2)); 348 CHK(type1 == type2); 349 350 switch(type1) { 351 case SDIS_ESTIMATOR_TEMPERATURE: 352 OK(sdis_estimator_get_temperature(e1, &mc1)); 353 OK(sdis_estimator_get_temperature(e2, &mc2)); 354 CHK(mc1.E + 3*mc1.SE >= mc2.E - 3*mc2.SE); 355 CHK(mc1.E - 3*mc1.SE <= mc2.E + 3*mc2.SE); 356 break; 357 358 case SDIS_ESTIMATOR_FLUX: 359 OK(sdis_estimator_get_convective_flux(e1, &mc1)); 360 OK(sdis_estimator_get_convective_flux(e2, &mc2)); 361 CHK(mc1.E + 3*mc1.SE >= mc2.E - 3*mc2.SE); 362 CHK(mc1.E - 3*mc1.SE <= mc2.E + 3*mc2.SE); 363 364 OK(sdis_estimator_get_radiative_flux(e1, &mc1)); 365 OK(sdis_estimator_get_radiative_flux(e2, &mc2)); 366 CHK(mc1.E + 3*mc1.SE >= mc2.E - 3*mc2.SE); 367 CHK(mc1.E - 3*mc1.SE <= mc2.E + 3*mc2.SE); 368 369 OK(sdis_estimator_get_total_flux(e1, &mc1)); 370 OK(sdis_estimator_get_total_flux(e2, &mc2)); 371 CHK(mc1.E + 3*mc1.SE >= mc2.E - 3*mc2.SE); 372 CHK(mc1.E - 3*mc1.SE <= mc2.E + 3*mc2.SE); 373 break; 374 375 case SDIS_ESTIMATOR_POWER: 376 OK(sdis_estimator_get_power(e1, &mc1)); 377 OK(sdis_estimator_get_power(e2, &mc2)); 378 CHK(mc1.E + 3*mc1.SE >= mc2.E - 3*mc2.SE); 379 CHK(mc1.E - 3*mc1.SE <= mc2.E + 3*mc2.SE); 380 break; 381 382 default: FATAL("Unreachable code.\n"); break; 383 } 384 } 385 386 static INLINE void 387 check_estimator_eq_strict 388 (const struct sdis_estimator* e1, const struct sdis_estimator* e2) 389 { 390 struct sdis_mc mc1, mc2; 391 enum sdis_estimator_type type1, type2; 392 ASSERT(e1 && e2); 393 394 OK(sdis_estimator_get_type(e1, &type1)); 395 OK(sdis_estimator_get_type(e2, &type2)); 396 CHK(type1 == type2); 397 398 switch(type1) { 399 case SDIS_ESTIMATOR_TEMPERATURE: 400 OK(sdis_estimator_get_temperature(e1, &mc1)); 401 OK(sdis_estimator_get_temperature(e2, &mc2)); 402 CHK(mc1.E == mc2.E && mc1.V == mc2.V && mc1.SE == mc2.SE); 403 break; 404 405 case SDIS_ESTIMATOR_FLUX: 406 OK(sdis_estimator_get_convective_flux(e1, &mc1)); 407 OK(sdis_estimator_get_convective_flux(e2, &mc2)); 408 CHK(mc1.E == mc2.E && mc1.V == mc2.V && mc1.SE == mc2.SE); 409 410 OK(sdis_estimator_get_radiative_flux(e1, &mc1)); 411 OK(sdis_estimator_get_radiative_flux(e2, &mc2)); 412 CHK(mc1.E == mc2.E && mc1.V == mc2.V && mc1.SE == mc2.SE); 413 414 OK(sdis_estimator_get_total_flux(e1, &mc1)); 415 OK(sdis_estimator_get_total_flux(e2, &mc2)); 416 CHK(mc1.E == mc2.E && mc1.V == mc2.V && mc1.SE == mc2.SE); 417 break; 418 419 case SDIS_ESTIMATOR_POWER: 420 OK(sdis_estimator_get_power(e1, &mc1)); 421 OK(sdis_estimator_get_power(e2, &mc2)); 422 CHK(mc1.E == mc2.E && mc1.V == mc2.V && mc1.SE == mc2.SE); 423 break; 424 425 default: FATAL("Unreachable code.\n"); break; 426 } 427 } 428 429 static INLINE void 430 check_memory_allocator(struct mem_allocator* allocator) 431 { 432 if(MEM_ALLOCATED_SIZE(allocator)) { 433 char dump[1024]; 434 MEM_DUMP(allocator, dump, sizeof(dump)); 435 fprintf(stderr, "%s\n", dump); 436 FATAL("Memory leaks.\n"); 437 } 438 } 439 440 extern LOCAL_SYM void 441 check_green_function 442 (struct sdis_green_function* green); 443 444 extern LOCAL_SYM void 445 dump_heat_paths 446 (FILE* stream, 447 const struct sdis_estimator* estimator); 448 449 extern LOCAL_SYM void 450 check_green_serialization 451 (struct sdis_green_function* green, 452 struct sdis_scene* scn); 453 454 #endif /* TEST_SDIS_UTILS_H */