test_sck_load.c (16058B)
1 /* Copyright (C) 2022, 2023 |Méso|Star> (contact@meso-star.com) 2 * Copyright (C) 2020, 2021 CNRS 3 * 4 * This program is free software: you can redistribute it and/or modify 5 * it under the terms of the GNU General Public License as published by 6 * the Free Software Foundation, either version 3 of the License, or 7 * (at your option) any later version. 8 * 9 * This program is distributed in the hope that it will be useful, 10 * but WITHOUT ANY WARRANTY; without even the implied warranty of 11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 12 * GNU General Public License for more details. 13 * 14 * You should have received a copy of the GNU General Public License 15 * along with this program. If not, see <http://www.gnu.org/licenses/>. */ 16 17 #define _POSIX_C_SOURCE 200112L /* nextarfter */ 18 19 #include "sck.h" 20 21 #include <rsys/math.h> 22 #include <rsys/mem_allocator.h> 23 #include <math.h> 24 #include <string.h> 25 26 static INLINE double 27 rand_canonic(void) 28 { 29 int r; 30 while((r = rand()) == RAND_MAX); 31 return (double)r / (double)RAND_MAX; 32 } 33 34 /******************************************************************************* 35 * Helper functions 36 ******************************************************************************/ 37 static void 38 check_sck_load 39 (const struct sck* sck, 40 const size_t nbands, 41 const size_t nnodes) 42 { 43 #define NMOMENTS 4 44 struct sck_band band = SCK_BAND_NULL; 45 struct sck_quad_pt qpt = SCK_QUAD_PT_NULL; 46 size_t iband; 47 size_t iqpt; 48 const size_t nsamps = 1000000; 49 size_t isamp; 50 51 CHK(sck); 52 CHK(nbands); 53 CHK(nnodes); 54 55 CHK(sck_validate(NULL) == RES_BAD_ARG); 56 CHK(sck_validate(sck) == RES_OK); 57 58 CHK(sck_get_bands_count(sck) == nbands); 59 CHK(sck_get_nodes_count(sck) == nnodes); 60 61 CHK(sck_get_band(NULL, 0, &band) == RES_BAD_ARG); 62 CHK(sck_get_band(sck, nbands, &band) == RES_BAD_ARG); 63 CHK(sck_get_band(sck, nbands, NULL) == RES_BAD_ARG); 64 CHK(sck_get_band(sck, 0, &band) == RES_OK); 65 66 CHK(band.quad_pts_count == 1); 67 CHK(sck_band_get_quad_pt(NULL, 0, &qpt) == RES_BAD_ARG); 68 CHK(sck_band_get_quad_pt(&band, 1, &qpt) == RES_BAD_ARG); 69 CHK(sck_band_get_quad_pt(&band, 0, NULL) == RES_BAD_ARG); 70 71 CHK(sck_band_sample_quad_pt(NULL, 0, &iqpt) == RES_BAD_ARG); 72 CHK(sck_band_sample_quad_pt(&band, 1, &iqpt) == RES_BAD_ARG); 73 CHK(sck_band_sample_quad_pt(&band, 0, NULL) == RES_BAD_ARG); 74 75 FOR_EACH(iband, 0, nbands) { 76 const double low = (double)(iband+1); 77 const double upp = (double)(iband+2); 78 const size_t nqpts = iband + 1; 79 double sum[NMOMENTS] = {0}; 80 double sum2[NMOMENTS] = {0}; 81 double E[NMOMENTS] = {0}; 82 double V[NMOMENTS] = {0}; 83 double SE[NMOMENTS] = {0}; 84 double ref; 85 size_t inode; 86 int imoment; 87 88 CHK(sck_get_band(sck, iband, &band) == RES_OK); 89 CHK(band.lower == low); 90 CHK(band.upper == upp); 91 CHK(band.id == iband); 92 CHK(band.quad_pts_count == nqpts); 93 94 FOR_EACH(inode, 0, nnodes) { 95 const float ks = (float)(iband*2000 + inode); 96 CHK(band.ks_list[inode] == ks); 97 } 98 99 FOR_EACH(iqpt, 0, nqpts) { 100 const double weight = 1.0/(double)nqpts; 101 102 CHK(sck_band_get_quad_pt(&band, iqpt, &qpt) == RES_OK); 103 CHK(qpt.weight == weight); 104 CHK(qpt.id == iqpt); 105 106 FOR_EACH(inode, 0, nnodes) { 107 const float ka = (float)(iband*1000 + iqpt*100 + inode); 108 CHK(qpt.ka_list[inode] == ka); 109 } 110 } 111 112 /* Check the sampling of the quadrature points */ 113 FOR_EACH(isamp, 0, nsamps) { 114 const double r = rand_canonic(); 115 double w[NMOMENTS]; 116 117 CHK(sck_band_sample_quad_pt(&band, r, &iqpt) == RES_OK); 118 ASSERT(iqpt < band.quad_pts_count); 119 CHK(sck_band_get_quad_pt(&band, iqpt, &qpt) == RES_OK); 120 FOR_EACH(imoment, 0, NMOMENTS) { 121 if(imoment == 0) { 122 w[imoment] = (double)(qpt.id + 1); 123 } else { 124 w[imoment] = (double)(qpt.id + 1) * w[imoment-1];; 125 } 126 sum[imoment] += w[imoment]; 127 sum2[imoment] += w[imoment]*w[imoment]; 128 } 129 } 130 131 FOR_EACH(imoment, 0, NMOMENTS) { 132 E[imoment] = sum[imoment] / (double)nsamps; 133 V[imoment] = sum2[imoment] / (double)nsamps - E[imoment]*E[imoment]; 134 SE[imoment] = sqrt(V[imoment]/(double)nsamps); 135 } 136 137 ref = (double)(nqpts+1)/2.0; 138 CHK(eq_eps(ref, E[0], 3*SE[0])); 139 140 ref = (double)(2*nqpts*nqpts + 3*nqpts + 1) / 6.0; 141 CHK(eq_eps(ref, E[1], 3*SE[1])); 142 143 ref = ((double)nqpts * pow((double)nqpts + 1.0, 2.0)) / 4.0; 144 CHK(eq_eps(ref, E[2], 3*SE[2])); 145 146 ref = (double)1.0/30.0 * (double) 147 (6*nqpts*nqpts*nqpts*nqpts + 15*nqpts*nqpts*nqpts + 10*nqpts*nqpts - 1); 148 CHK(eq_eps(ref, E[3], 3*SE[3])); 149 } 150 } 151 152 static void 153 write_sck 154 (FILE* fp, 155 const uint64_t pagesize, 156 const uint64_t nbands, 157 const uint64_t nnodes) 158 { 159 uint64_t iband; 160 const char byte = 0; 161 162 /* Write the header */ 163 CHK(fwrite(&pagesize, sizeof(pagesize), 1, fp) == 1); 164 CHK(fwrite(&nbands, sizeof(nbands), 1, fp) == 1); 165 CHK(fwrite(&nnodes, sizeof(nnodes), 1, fp) == 1); 166 167 FOR_EACH(iband, 0, nbands) { 168 const double low = (double)(iband+1); 169 const double upp = (double)(iband+2); 170 const uint64_t nqpts = iband + 1; 171 uint64_t iqpt; 172 173 /* Write band description */ 174 CHK(fwrite(&low, sizeof(low), 1, fp) == 1); 175 CHK(fwrite(&upp, sizeof(upp), 1, fp) == 1); 176 CHK(fwrite(&nqpts, sizeof(nqpts), 1, fp) == 1); 177 178 /* Write per band quadrature points */ 179 FOR_EACH(iqpt, 0, nqpts) { 180 const double weight = 1.0/(double)nqpts; 181 CHK(fwrite(&weight, sizeof(weight), 1, fp) == 1); 182 } 183 } 184 185 /* Write per band ks and per quadrature point ka */ 186 FOR_EACH(iband, 0, nbands) { 187 const uint64_t nqpts = iband + 1; 188 uint64_t inode; 189 uint64_t iqpt; 190 /* Padding */ 191 CHK(fseek(fp, (long)ALIGN_SIZE((size_t)ftell(fp), pagesize), SEEK_SET)==0); 192 193 FOR_EACH(inode, 0, nnodes) { 194 const float ks = (float)(iband*2000 + inode); 195 CHK(fwrite(&ks, sizeof(ks), 1, fp) == 1); 196 } 197 198 FOR_EACH(iqpt, 0, nqpts) { 199 200 /* Padding */ 201 CHK(fseek(fp, (long)ALIGN_SIZE((size_t)ftell(fp), pagesize), SEEK_SET)==0); 202 203 FOR_EACH(inode, 0, nnodes) { 204 const float ka = (float)(iband*1000 + iqpt*100 + inode); 205 CHK(fwrite(&ka, sizeof(ka), 1, fp) == 1); 206 } 207 } 208 } 209 210 /* Padding. Write one char to position the EOF indicator */ 211 CHK(fseek(fp, (long)ALIGN_SIZE((size_t)ftell(fp), pagesize)-1, SEEK_SET) == 0); 212 CHK(fwrite(&byte, sizeof(byte), 1, fp) == 1); 213 } 214 215 static void 216 test_load(struct sck* sck) 217 { 218 hash256_T hash0; 219 hash256_T hash1; 220 hash256_T band_hash0; 221 hash256_T band_hash1; 222 hash256_T pt_hash0; 223 hash256_T pt_hash1; 224 struct sck_load_args args = SCK_LOAD_ARGS_NULL; 225 struct sck_load_stream_args stream_args = SCK_LOAD_STREAM_ARGS_NULL; 226 struct sck_band band; 227 228 FILE* fp = NULL; 229 const char* filename = "test_file.sck"; 230 const uint64_t pagesize = 16384; 231 const uint64_t nbands = 11; 232 const uint64_t nnodes = 1000; 233 234 CHK(fp = fopen(filename, "w+")); 235 write_sck(fp, pagesize, nbands, nnodes); 236 rewind(fp); 237 238 stream_args.stream = fp; 239 stream_args.name = filename; 240 CHK(sck_load_stream(NULL, &stream_args) == RES_BAD_ARG); 241 CHK(sck_load_stream(sck, NULL) == RES_BAD_ARG); 242 stream_args.stream = NULL; 243 CHK(sck_load_stream(sck, NULL) == RES_BAD_ARG); 244 stream_args.stream = fp; 245 CHK(sck_load_stream(sck, &stream_args) == RES_OK); 246 247 CHK(!strcmp(sck_get_name(sck), filename)); 248 249 CHK(sck_compute_hash(sck, NULL) == RES_BAD_ARG); 250 CHK(sck_compute_hash(NULL, hash0) == RES_BAD_ARG); 251 CHK(sck_compute_hash(sck, hash0) == RES_OK); 252 253 CHK(sck_band_compute_hash(NULL, 0, band_hash0) == RES_BAD_ARG); 254 CHK(sck_band_compute_hash(sck, nbands, band_hash0) == RES_BAD_ARG); 255 CHK(sck_band_compute_hash(sck, 0, NULL) == RES_BAD_ARG); 256 CHK(sck_band_compute_hash(sck, 0, band_hash0) == RES_OK); 257 CHK(!hash256_eq(band_hash0, hash0)); 258 259 CHK(sck_get_band(sck, 0, &band) == RES_OK); 260 CHK(sck_quad_pt_compute_hash(NULL, 0, pt_hash0) == RES_BAD_ARG); 261 CHK(sck_quad_pt_compute_hash(&band, band.quad_pts_count, pt_hash0) == RES_BAD_ARG); 262 CHK(sck_quad_pt_compute_hash(&band, 0, NULL) == RES_BAD_ARG); 263 CHK(sck_quad_pt_compute_hash(&band, 0, pt_hash0) == RES_OK); 264 CHK(!hash256_eq(pt_hash0, band_hash0)); 265 CHK(!hash256_eq(pt_hash0, hash0)); 266 267 check_sck_load(sck, nbands, nnodes); 268 rewind(fp); 269 270 stream_args.name = NULL; 271 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 272 stream_args.name = SCK_LOAD_STREAM_ARGS_NULL.name; 273 stream_args.memory_mapping = 1; 274 CHK(sck_load_stream(sck, &stream_args) == RES_OK); 275 CHK(!strcmp(sck_get_name(sck), SCK_LOAD_STREAM_ARGS_NULL.name)); 276 277 args.path = "nop"; 278 CHK(sck_load(NULL, &args) == RES_BAD_ARG); 279 CHK(sck_load(sck, NULL) == RES_BAD_ARG); 280 CHK(sck_load(sck, &args) == RES_IO_ERR); 281 args.path = filename; 282 CHK(sck_load(sck, &args) == RES_OK); 283 CHK(!strcmp(sck_get_name(sck), args.path)); 284 check_sck_load(sck, nbands, nnodes); 285 286 CHK(sck_compute_hash(sck, hash1) == RES_OK); 287 CHK(hash256_eq(hash0, hash1)); 288 289 rewind(fp); 290 write_sck(fp, pagesize, nbands+1, nnodes); 291 rewind(fp); 292 293 stream_args.stream = fp; 294 stream_args.name = filename; 295 CHK(sck_load_stream(sck, &stream_args) == RES_OK); 296 CHK(sck_compute_hash(sck, hash1) == RES_OK); 297 CHK(!hash256_eq(hash0, hash1)); 298 299 CHK(sck_band_compute_hash(sck, 0, band_hash1) == RES_OK); 300 CHK(hash256_eq(band_hash0, band_hash1)); 301 302 CHK(sck_get_band(sck, 0, &band) == RES_OK); 303 CHK(sck_quad_pt_compute_hash(&band, 0, pt_hash1) == RES_OK); 304 CHK(hash256_eq(pt_hash0, pt_hash1)); 305 306 CHK(fclose(fp) == 0); 307 } 308 309 static void 310 test_load_fail(struct sck* sck) 311 { 312 struct sck_load_stream_args stream_args = SCK_LOAD_STREAM_ARGS_NULL; 313 FILE* fp = NULL; 314 double low; 315 double upp; 316 317 318 /* The pagesize is less than the operating system page size*/ 319 CHK(fp = tmpfile()); 320 write_sck(fp, 2048, 1, 1); 321 rewind(fp); 322 stream_args.stream = fp; 323 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 324 CHK(fclose(fp) == 0); 325 326 /* The pagesize is not a power of two */ 327 CHK(fp = tmpfile()); 328 write_sck(fp, 4100, 1, 1); 329 rewind(fp); 330 stream_args.stream = fp; 331 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 332 CHK(fclose(fp) == 0); 333 334 /* Wrong #bands */ 335 CHK(fp = tmpfile()); 336 write_sck(fp, 4096, 0, 1); 337 rewind(fp); 338 stream_args.stream = fp; 339 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 340 CHK(fclose(fp) == 0); 341 342 /* Wrong #nodes */ 343 CHK(fp = tmpfile()); 344 write_sck(fp, 4096, 1, 0); 345 rewind(fp); 346 stream_args.stream = fp; 347 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 348 CHK(fclose(fp) == 0); 349 350 /* Wrong band boundaries */ 351 low = 1; 352 upp = 0; 353 CHK(fp = tmpfile()); 354 write_sck(fp, 4096, 1, 1); 355 CHK(fseek(fp, 24, SEEK_SET) == 0); 356 CHK(fwrite(&low, sizeof(low), 1, fp) == 1); 357 CHK(fwrite(&upp, sizeof(upp), 1, fp) == 1); 358 rewind(fp); 359 stream_args.stream = fp; 360 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 361 CHK(fclose(fp) == 0); 362 363 /* Unsorted bands */ 364 CHK(fp = tmpfile()); 365 write_sck(fp, 4096, 2, 1); 366 CHK(fseek(fp, 24, SEEK_SET) == 0); 367 low = 1; upp = 2; 368 CHK(fwrite(&low, sizeof(low), 1, fp) == 1); 369 CHK(fwrite(&upp, sizeof(upp), 1, fp) == 1); 370 CHK(fseek(fp, 56, SEEK_SET) == 0); 371 low = 0; upp = 1; 372 CHK(fwrite(&low, sizeof(low), 1, fp) == 1); 373 CHK(fwrite(&upp, sizeof(upp), 1, fp) == 1); 374 rewind(fp); 375 stream_args.stream = fp; 376 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 377 CHK(fclose(fp) == 0); 378 379 /* Bands overlap */ 380 CHK(fp = tmpfile()); 381 write_sck(fp, 4096, 2, 1); 382 CHK(fseek(fp, 24, SEEK_SET) == 0); 383 low = 0; upp = 1; 384 CHK(fwrite(&low, sizeof(low), 1, fp) == 1); 385 CHK(fwrite(&upp, sizeof(upp), 1, fp) == 1); 386 CHK(fseek(fp, 56, SEEK_SET) == 0); 387 low = 0.5; upp = 1.5; 388 CHK(fwrite(&low, sizeof(low), 1, fp) == 1); 389 CHK(fwrite(&upp, sizeof(upp), 1, fp) == 1); 390 rewind(fp); 391 stream_args.stream = fp; 392 CHK(sck_load_stream(sck, &stream_args) == RES_BAD_ARG); 393 CHK(fclose(fp) == 0); 394 } 395 396 static void 397 test_load_files(struct sck* sck, int argc, char** argv) 398 { 399 int i; 400 CHK(sck); 401 FOR_EACH(i, 1, argc) { 402 hash256_T hash; 403 size_t nnodes; 404 size_t nbands; 405 size_t iband; 406 407 if(!strcmp(argv[i], "-")) { 408 struct sck_load_stream_args args = SCK_LOAD_STREAM_ARGS_NULL; 409 printf("Load from stdin\n"); 410 args.stream = stdin; 411 args.name = "stdin"; 412 args.memory_mapping = 1; 413 CHK(sck_load_stream(sck, &args) == RES_BAD_ARG); 414 args.memory_mapping = 0; 415 CHK(sck_load_stream(sck, &args) == RES_OK); 416 } else { 417 struct sck_load_args args = SCK_LOAD_ARGS_NULL; 418 printf("Load %s\n", argv[1]); 419 args.path = argv[i]; 420 args.memory_mapping = 1; 421 CHK(sck_load(sck, &args) == RES_OK); 422 } 423 nbands = sck_get_bands_count(sck); 424 nnodes = sck_get_nodes_count(sck); 425 CHK(nbands); 426 CHK(nnodes); 427 428 FOR_EACH(iband, 0, nbands) { 429 struct sck_band band = SCK_BAND_NULL; 430 431 CHK(sck_get_band(sck, iband, &band) == RES_OK); 432 printf("band %lu in [%g, %g[ nm with %lu quadrature points\n", 433 (unsigned long)band.id, 434 band.lower, band.upper, 435 (unsigned long)band.quad_pts_count); 436 } 437 438 CHK(sck_validate(sck) == RES_OK); 439 CHK(sck_compute_hash(sck, hash) == RES_OK); 440 } 441 } 442 443 static void 444 test_find(struct sck* sck) 445 { 446 struct sck_load_stream_args stream_args = SCK_LOAD_STREAM_ARGS_NULL; 447 size_t ibands[2]; 448 double range[2]; 449 FILE* fp; 450 451 CHK(fp = tmpfile()); 452 write_sck(fp, 4096, 10, 1); 453 rewind(fp); 454 stream_args.stream = fp; 455 stream_args.memory_mapping = 1; 456 CHK(sck_load_stream(sck, &stream_args) == RES_OK); 457 458 range[0] = 0; 459 range[1] = 10; 460 CHK(sck_find_bands(NULL, range, ibands) == RES_BAD_ARG); 461 CHK(sck_find_bands(sck, NULL, ibands) == RES_BAD_ARG); 462 CHK(sck_find_bands(sck, range, NULL) == RES_BAD_ARG); 463 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 464 CHK(ibands[0] == 0 && ibands[1] == 9); 465 466 range[0] = 10; 467 range[1] = 0; 468 CHK(sck_find_bands(sck, range, ibands) == RES_BAD_ARG); 469 470 range[0] = 11; 471 range[1] = 12; 472 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 473 CHK(ibands[0] > ibands[1]); 474 475 range[0] = 11; 476 range[1] = 11; 477 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 478 CHK(ibands[0] > ibands[1]); 479 480 range[0] = 0; 481 range[1] = nextafter(1, 0); 482 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 483 CHK(ibands[0] > ibands[1]); 484 485 range[0] = 0; 486 range[1] = 0; 487 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 488 CHK(ibands[0] > ibands[1]); 489 490 range[0] = nextafter(1, 0); 491 range[1] = nextafter(1, 0); 492 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 493 CHK(ibands[0] > ibands[1]); 494 495 range[0] = 2; 496 range[1] = 2; 497 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 498 CHK(ibands[0] == 1 && ibands[1] == 1); 499 500 range[0] = 0; 501 range[1] = 1; 502 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 503 CHK(ibands[0] == 0 && ibands[1] == 0); 504 505 range[0] = 1; 506 range[1] = nextafterf(2, 1); 507 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 508 CHK(ibands[0] == 0 && ibands[1] == 0); 509 510 range[0] = 2; 511 range[1] = 20; 512 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 513 CHK(ibands[0] == 1 && ibands[1] == 9); 514 515 range[0] = 1.5; 516 range[1] = 2; 517 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 518 CHK(ibands[0] == 0 && ibands[1] == 1); 519 520 range[0] = 3.1; 521 range[1] = nextafter(6, 0); 522 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 523 CHK(ibands[0] == 2 && ibands[1] == 4); 524 525 range[0] = 3.1; 526 range[1] = 7; 527 CHK(sck_find_bands(sck, range, ibands) == RES_OK); 528 CHK(ibands[0] == 2 && ibands[1] == 6); 529 530 CHK(fclose(fp) == 0); 531 } 532 533 /******************************************************************************* 534 * Main function 535 ******************************************************************************/ 536 int 537 main(int argc, char** argv) 538 { 539 struct sck_create_args args = SCK_CREATE_ARGS_DEFAULT; 540 struct sck* sck = NULL; 541 (void)argc, (void)argv; 542 543 args.verbose = 1; 544 CHK(sck_create(&args, &sck) == RES_OK); 545 if(argc > 1) { 546 test_load_files(sck, argc, argv); 547 } else { 548 test_load(sck); 549 test_load_fail(sck); 550 test_find(sck); 551 } 552 553 CHK(sck_ref_put(sck) == RES_OK); 554 CHK(mem_allocated_size() == 0); 555 return 0; 556 }