green-compute.c (19600B)
1 /* Copyright (C) 2020-2022, 2024 |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 #define _POSIX_C_SOURCE 200112L /* strtok_r */ 17 18 #include "green-compute.h" 19 #include "green-types.h" 20 21 #include <rsys/str.h> 22 #include <rsys/cstr.h> 23 #include <rsys/logger.h> 24 #include <rsys/text_reader.h> 25 #include <rsys/dynamic_array.h> 26 #include <rsys/mem_allocator.h> 27 28 #include <math.h> 29 #include <string.h> 30 #include <ctype.h> 31 #include <omp.h> 32 33 enum log_record_type { 34 RECORD_INVALID = BIT(0), 35 RECORD_UNKNOWN = BIT(1), 36 RECORD_UNUSED = BIT(2), 37 RECORD_DEFAULT = BIT(3), 38 RECORD_TYPES_UNDEF__ = BIT(4) 39 }; 40 41 struct log_record { 42 enum log_record_type type; 43 int line_num, var_num; 44 struct str name; 45 const char* keep_line; 46 double default_value; 47 struct mem_allocator* allocator; 48 }; 49 50 static FINLINE res_T 51 log_record_create 52 (struct mem_allocator* allocator, 53 const enum log_record_type type, 54 struct log_record** out_record) 55 { 56 struct log_record* record; 57 ASSERT(allocator && out_record); 58 record = MEM_ALLOC(allocator, sizeof(*record)); 59 if(!record) return RES_MEM_ERR; 60 record->type = type; 61 record->line_num = -1; 62 record->var_num = -1; 63 record->keep_line = NULL; 64 record->default_value = INF; 65 str_init(allocator, &record->name); 66 record->allocator = allocator; 67 *out_record = record; 68 return RES_OK; 69 } 70 71 static FINLINE void 72 log_record_release 73 (struct log_record* data) 74 { 75 str_release(&data->name); 76 MEM_RM(data->allocator, data); 77 } 78 79 #define DARRAY_NAME logs 80 #define DARRAY_DATA struct log_record* 81 #include <rsys/dynamic_array.h> 82 83 void 84 check_green_table_variables_use 85 (struct green* green) 86 { 87 size_t i; 88 unsigned j; 89 90 ASSERT(green); 91 ASSERT(!green->references_checked); /* Fix double call in caller */ 92 93 FOR_EACH(i, 0, green->header.ok_count) { 94 const struct sample* sample = green->samples + i; 95 unsigned id = sample->header.sample_end_description_id; 96 /* Ambient ID is description_count */ 97 ASSERT(id <= green->header.description_count); 98 ASSERT(id == green->header.description_count 99 || green->descriptions[id].type == GREEN_MAT_SOLID 100 || green->descriptions[id].type == GREEN_MAT_FLUID 101 || green->descriptions[id].type == GREEN_BOUND_H 102 || green->descriptions[id].type == GREEN_BOUND_T); 103 if(sample->header.at_initial) 104 green->table[id].initial_T_weight++; 105 else 106 green->table[id].imposed_T_weight++; 107 108 FOR_EACH(j, 0, sample->header.pw_count) { 109 const double w = sample->pw_weights[j]; 110 id = sample->pw_ids[j]; 111 ASSERT(id < green->header.description_count); 112 ASSERT(green->descriptions[id].type == GREEN_MAT_SOLID 113 || green->descriptions[id].type == GREEN_MAT_FLUID); 114 green->table[id].other_weight += w; 115 } 116 FOR_EACH(j, 0, sample->header.fx_count) { 117 const double w = sample->fx_weights[j]; 118 id = sample->fx_ids[j]; 119 ASSERT(id < green->header.description_count); 120 ASSERT(green->descriptions[id].type == GREEN_BOUND_F); 121 green->table[id].other_weight += w; 122 } 123 } 124 green->references_checked = 1; 125 } 126 127 #define INSERT(Elt, Field) { \ 128 struct variable_data vd___; \ 129 vd___.idx = ((Elt) * sizeof(struct applied_settings_elt) \ 130 + offsetof(struct applied_settings_elt, Field ## _value)) / sizeof(double); \ 131 vd___.used = (green->table[Elt].Field ## _weight != 0); \ 132 vd___.unknown = green->table[Elt].Field ## _unknown; \ 133 if(vd___.unknown) green->unknown_variables = 1; \ 134 else if(!vd___.used) green->unused_variables = 1; \ 135 if(htable_variable_ptr_find(&green->variable_ptrs, \ 136 &(green->table[Elt].Field ## _name))) {\ 137 /* Should have been detected by stardis: corrupted file??? */ \ 138 logger_print(green->logger, LOG_ERROR, \ 139 "Duplicate description name found.\n"); \ 140 res = RES_BAD_ARG; \ 141 goto error; \ 142 } \ 143 ERR(htable_variable_ptr_set(&green->variable_ptrs, \ 144 &(green->table[Elt].Field ## _name), &vd___)); \ 145 } (void)0 146 147 #define MK_VAR(Str, Desc) \ 148 for(;;) { \ 149 size_t n = (size_t)snprintf(buf, bufsz, (Str), (Desc).name);\ 150 if(n < bufsz) break; \ 151 buf = MEM_REALLOC(green->allocator, buf, n + 1); \ 152 bufsz = n + 1; \ 153 } (void)0 154 155 res_T 156 build_green_table 157 (struct green* green) 158 { 159 res_T res = RES_OK; 160 size_t i, ambient_id; 161 char* buf = NULL; 162 size_t bufsz = 32; 163 164 ASSERT(green); 165 166 green->table = MEM_CALLOC(green->allocator, 1+green->header.description_count, 167 sizeof(*green->table)); 168 buf = MEM_ALLOC(green->allocator, bufsz); 169 if(!green->table || ! buf) { 170 res = RES_MEM_ERR; 171 goto error; 172 } 173 FOR_EACH(i, 0, 1 + green->header.description_count) 174 init_table_elt(green->allocator, green->table + i); 175 176 check_green_table_variables_use(green); 177 178 FOR_EACH(i, 0, green->header.description_count) { 179 struct green_description* desc = green->descriptions + i; 180 switch(desc->type) { 181 case GREEN_BOUND_T: 182 MK_VAR("%s.T", desc->d.t_boundary); 183 ERR(str_set(&green->table[i].imposed_T_name, buf)); 184 ASSERT(desc->d.t_boundary.imposed_temperature >= 0); 185 green->table[i].imposed_T_value = desc->d.t_boundary.imposed_temperature; 186 green->table[i].imposed_T_defined = 1; 187 INSERT(i, imposed_T); 188 break; 189 case GREEN_BOUND_H: 190 MK_VAR("%s.T", desc->d.h_boundary); 191 ERR(str_set(&green->table[i].imposed_T_name, buf)); 192 ASSERT(desc->d.h_boundary.imposed_temperature >= 0); 193 green->table[i].imposed_T_value = desc->d.h_boundary.imposed_temperature; 194 green->table[i].imposed_T_defined = 1; 195 INSERT(i, imposed_T); 196 break; 197 case GREEN_BOUND_F: 198 MK_VAR("%s.F", desc->d.f_boundary); 199 ERR(str_set(&green->table[i].other_name, buf)); 200 green->table[i].other_value = desc->d.f_boundary.imposed_flux; 201 green->table[i].other_defined = 1; 202 INSERT(i, other); 203 break; 204 case GREEN_MAT_SOLID: 205 MK_VAR("%s.T", desc->d.solid); 206 ERR(str_set(&green->table[i].imposed_T_name, buf)); 207 green->table[i].imposed_T_value = desc->d.solid.imposed_temperature; 208 green->table[i].imposed_T_defined = 1; 209 if(green->table[i].imposed_T_value < 0) 210 green->table[i].imposed_T_unknown = 1; 211 INSERT(i, imposed_T); 212 MK_VAR("%s.Tinit", desc->d.solid); 213 ERR(str_set(&green->table[i].initial_T_name, buf)); 214 green->table[i].initial_T_value = desc->d.solid.initial_temperature; 215 green->table[i].initial_T_defined = 1; 216 if(green->table[i].initial_T_value < 0) 217 green->table[i].initial_T_unknown = 1; 218 INSERT(i, initial_T); 219 MK_VAR("%s.VP", desc->d.solid); 220 ERR(str_set(&green->table[i].other_name, buf)); 221 green->table[i].other_value = desc->d.solid.volumic_power; 222 green->table[i].other_defined = 1; 223 INSERT(i, other); 224 break; 225 case GREEN_MAT_FLUID: 226 MK_VAR("%s.T", desc->d.fluid); 227 ERR(str_set(&green->table[i].imposed_T_name, buf)); 228 green->table[i].imposed_T_value = desc->d.fluid.imposed_temperature; 229 green->table[i].imposed_T_defined = 1; 230 if(green->table[i].imposed_T_value < 0) 231 green->table[i].imposed_T_unknown = 1; 232 INSERT(i, imposed_T); 233 MK_VAR("%s.Tinit", desc->d.fluid); 234 ERR(str_set(&green->table[i].initial_T_name, buf)); 235 green->table[i].initial_T_value = desc->d.fluid.initial_temperature; 236 green->table[i].initial_T_defined = 1; 237 if(green->table[i].initial_T_value < 0) 238 green->table[i].initial_T_unknown = 1; 239 INSERT(i, initial_T); 240 break; 241 case GREEN_SOLID_SOLID_CONNECT: 242 case GREEN_SOLID_FLUID_CONNECT: 243 /* No variables here */ 244 break; 245 default: 246 logger_print(green->logger, LOG_ERROR, "Invalid description type.\n"); 247 res = RES_BAD_ARG; 248 goto error; 249 } 250 } 251 /* Ambient ID is description_count */ 252 ambient_id = green->header.description_count; 253 ERR(str_set(&green->table[ambient_id].imposed_T_name, "AMBIENT")); 254 green->table[ambient_id].imposed_T_value 255 = green->header.ambient_radiative_temperature; 256 green->table[ambient_id].imposed_T_defined = 1; 257 INSERT(ambient_id, imposed_T); 258 259 end: 260 MEM_RM(green->allocator, buf); 261 return res; 262 error: 263 goto end; 264 } 265 266 #undef INSERT 267 268 static void 269 applied_settings_set_defaults 270 (const struct green* green, 271 struct applied_settings_elt* settings) 272 { 273 unsigned i; 274 FOR_EACH(i, 0, 1 + green->header.description_count) { 275 settings[i].imposed_T_value = green->table[i].imposed_T_value; 276 settings[i].initial_T_value = green->table[i].initial_T_value; 277 settings[i].other_value = green->table[i].other_value; 278 } 279 } 280 281 static void 282 green_compute_1 283 (const struct green* green, 284 const struct applied_settings_elt* settings, 285 double* E, 286 double* STD) 287 { 288 size_t i; 289 unsigned j; 290 double V, Ti, T = 0, T2 = 0; 291 292 ASSERT(green && E && STD); 293 294 FOR_EACH(i, 0, green->header.ok_count) { 295 const struct sample* sample = green->samples + i; 296 unsigned id = sample->header.sample_end_description_id; 297 ASSERT(id <= green->header.description_count); /* Ambient ID is description_count */ 298 ASSERT(id == green->header.description_count 299 || (green->descriptions[id].type == GREEN_BOUND_T) 300 || (green->descriptions[id].type == GREEN_BOUND_H) 301 || (green->descriptions[id].type == GREEN_MAT_SOLID) 302 || (green->descriptions[id].type == GREEN_MAT_FLUID)); 303 if(sample->header.at_initial) { 304 ASSERT(green->table[id].initial_T_defined 305 && !green->table[id].initial_T_unknown); 306 Ti = settings[id].initial_T_value; 307 } else { 308 ASSERT(green->table[id].imposed_T_defined 309 && !green->table[id].imposed_T_unknown); 310 Ti = settings[id].imposed_T_value; 311 } 312 313 FOR_EACH(j, 0, sample->header.pw_count) { 314 const double w = sample->pw_weights[j]; 315 id = sample->pw_ids[j]; 316 ASSERT(id < green->header.description_count); 317 ASSERT(green->descriptions[id].type == GREEN_MAT_SOLID 318 || green->descriptions[id].type == GREEN_MAT_FLUID); 319 ASSERT(green->table[id].other_defined && !green->table[id].other_unknown); 320 Ti += w * settings[id].other_value; 321 } 322 FOR_EACH(j, 0, sample->header.fx_count) { 323 const double w = sample->fx_weights[j]; 324 id = sample->fx_ids[j]; 325 ASSERT(id < green->header.description_count); 326 ASSERT(green->descriptions[id].type == GREEN_BOUND_F); 327 ASSERT(green->table[id].other_defined && !green->table[id].other_unknown); 328 Ti += w * settings[id].other_value; 329 } 330 T += Ti; 331 T2 += Ti * Ti; 332 } 333 *E = T / (double)green->header.ok_count; 334 V = T2 / (double)green->header.ok_count - *E * *E; 335 *STD = (V > 0) ? sqrt(V / (double)green->header.ok_count) : 0; 336 } 337 338 static int 339 cmp_records 340 (void const* r1, void const* r2) 341 { 342 struct log_record* record1; 343 struct log_record* record2; 344 ASSERT(r1 && r2); 345 record1 = *(struct log_record**)r1; 346 record2 = *(struct log_record**)r2; 347 if(record1->line_num == record2->line_num) 348 return record1->var_num - record2->var_num; 349 return record1->line_num - record2->line_num; 350 } 351 352 static void 353 do_log 354 (const char* file_name, 355 struct green* green, 356 struct darray_logs* logs) 357 { 358 size_t i; 359 int prev_line_num = INT_MAX; 360 enum log_type log_type = LOG_WARNING; 361 ASSERT(file_name && green && logs); 362 363 if(darray_logs_size_get(logs) == 0) return; 364 365 /* Sort records by #line */ 366 qsort(darray_logs_data_get(logs), 367 darray_logs_size_get(logs), 368 sizeof(*darray_logs_data_get(logs)), 369 cmp_records); 370 371 FOR_EACH(i, 0, darray_logs_size_get(logs)) { 372 struct log_record* record = darray_logs_data_get(logs)[i]; 373 if(record->type & (RECORD_INVALID | RECORD_UNKNOWN)) { 374 log_type = LOG_ERROR; 375 break; 376 } 377 } 378 logger_print(green->logger, log_type, "In file '%s':\n", file_name); 379 FOR_EACH(i, 0, darray_logs_size_get(logs)) { 380 struct log_record* record = darray_logs_data_get(logs)[i]; 381 switch (record->type) { 382 case RECORD_INVALID: 383 if(record->line_num != prev_line_num) { 384 logger_print(green->logger, LOG_ERROR, 385 "In line #%d: '%s':\n", 386 record->line_num, record->keep_line); 387 logger_print(green->logger, LOG_ERROR, 388 (str_is_empty(&record->name) 389 ? "Invalid data (missing name)\n" : "Invalid data (missing value)\n")); 390 } else 391 logger_print(green->logger, LOG_ERROR, 392 (str_is_empty(&record->name) 393 ? "Invalid data (missing name)\n" : "Invalid data (missing value)\n")); 394 break; 395 case RECORD_UNKNOWN: 396 if(record->line_num != prev_line_num) { 397 logger_print(green->logger, LOG_ERROR, 398 "In line #%d: '%s':\n", 399 record->line_num, record->keep_line); 400 logger_print(green->logger, LOG_ERROR, 401 "Unknown variable name: '%s'\n", 402 str_cget(&record->name)); 403 } else 404 logger_print(green->logger, LOG_ERROR, 405 "\tUnknown variable name: '%s'\n", 406 str_cget(&record->name)); 407 break; 408 case RECORD_UNUSED: 409 if(record->line_num != prev_line_num) { 410 logger_print(green->logger, LOG_WARNING, 411 "In line #%d: '%s':\n", 412 record->line_num, record->keep_line); 413 logger_print(green->logger, LOG_WARNING, 414 "Attempt to change unused variable '%s'.\n", 415 str_cget(&record->name)); 416 } else 417 logger_print(green->logger, LOG_WARNING, 418 "\tAttempt to change unused variable '%s'.\n", 419 str_cget(&record->name)); 420 break; 421 case RECORD_DEFAULT: 422 if(record->line_num != prev_line_num) { 423 logger_print(green->logger, LOG_WARNING, 424 "In line #%d: '%s':\n", 425 record->line_num, record->keep_line); 426 logger_print(green->logger, LOG_WARNING, 427 "Change variable '%s' to its default value %g.\n", 428 str_cget(&record->name), record->default_value); 429 } else 430 logger_print(green->logger, LOG_WARNING, 431 "Change variable '%s' to its default value %g.\n", 432 str_cget(&record->name), record->default_value); 433 break; 434 default: 435 FATAL("error:" STR(__FILE__) ":" STR(__LINE__)": Invalid type.\n"); 436 } 437 prev_line_num = record->line_num; 438 } 439 } 440 441 static res_T 442 parse_line 443 (const char* file_name, 444 const int idx, 445 struct green* green, 446 struct applied_settings_elt* settings, 447 struct darray_logs* logs) 448 { 449 res_T res = RES_OK; 450 int name_count; 451 struct str name; 452 char* tk_name = NULL; 453 char* tk_value = NULL; 454 char* tok_ctx = NULL; 455 const struct str* keep_line; 456 struct str line; 457 ASSERT(file_name && green && settings && logs && idx >= 0); 458 (void)file_name; 459 460 keep_line = darray_str_cdata_get(&green->settings) + idx; 461 str_init(green->allocator, &line); 462 str_init(green->allocator, &name); 463 464 ERR(str_set(&line, str_cget(keep_line))); 465 466 /* At least one name=value, no name without value */ 467 for(name_count = 0; ; name_count++) { 468 struct variable_data* vd; 469 struct log_record* record = NULL; 470 double prev; 471 472 if(res != RES_OK) goto error; 473 474 tk_name = strtok_r(name_count ? NULL : str_get(&line), "=", &tok_ctx); 475 if(tk_name) { 476 char* c; 477 c = tk_name + strlen(tk_name) - 1; 478 while(c >= tk_name && isspace(*tk_name)) tk_name++; 479 while(c >= tk_name && isspace(*c)) { 480 *c = '\0'; 481 c--; 482 } 483 if(c < tk_name) tk_name = NULL; 484 } 485 tk_value = strtok_r(NULL, " \t", &tok_ctx); 486 if((tk_name == NULL) != (tk_value == NULL)) { 487 /* '<name> = <value>' or '' */ 488 ERR(log_record_create(green->allocator, RECORD_INVALID, &record)); 489 if(tk_name) { ERR(str_set(&record->name, tk_name)); } 490 res = RES_BAD_ARG; 491 goto push_record; 492 } 493 if(!tk_name) break; /* End of data */ 494 /* Check name validity */ 495 ERR(str_set(&name, tk_name)); 496 vd = htable_variable_ptr_find(&green->variable_ptrs, &name); 497 if(!vd) { 498 ERR(log_record_create(green->allocator, RECORD_UNKNOWN, &record)); 499 ERR(str_set(&record->name, tk_name)); 500 res = RES_BAD_ARG; 501 goto push_record; 502 } 503 if(!vd->used) { 504 ERR(log_record_create(green->allocator, RECORD_UNUSED, &record)); 505 ERR(str_set(&record->name, tk_name)); 506 goto push_record; 507 } 508 /* Variable exists and is used in Green function: check if really changed */ 509 prev = ((double*)settings)[vd->idx]; 510 ERR(cstr_to_double(tk_value, ((double*)settings) + vd->idx)); 511 if(prev == ((double*)settings)[vd->idx]) { 512 ERR(log_record_create(green->allocator, RECORD_DEFAULT, &record)); 513 ERR(str_set(&record->name, tk_name)); 514 record->default_value = prev; 515 goto push_record; 516 } 517 continue; 518 push_record: 519 #pragma omp critical 520 { 521 res_T tmp_res; 522 ASSERT(record && record->type != RECORD_TYPES_UNDEF__); 523 record->line_num = idx; 524 record->var_num = name_count; 525 record->keep_line = str_cget(keep_line); 526 tmp_res = darray_logs_push_back(logs, &record); 527 if(res == RES_OK) res = tmp_res; 528 } 529 } 530 531 end: 532 str_release(&line); 533 str_release(&name); 534 return res; 535 error: 536 goto end; 537 } 538 539 res_T 540 green_compute 541 (struct green* green, 542 const char* in_name) 543 { 544 res_T res = RES_OK; 545 size_t sz, n; 546 int i; 547 struct applied_settings_elt** settings = NULL; 548 struct darray_logs logs; 549 550 ASSERT(green && in_name); 551 sz = darray_str_size_get(&green->settings); 552 ASSERT(sz < INT_MAX); 553 554 settings = MEM_CALLOC(green->allocator, green->nthreads, sizeof(*settings)); 555 if(settings == NULL) { 556 res = RES_MEM_ERR; 557 goto error; 558 } 559 darray_logs_init(green->allocator, &logs); 560 561 omp_set_num_threads((int)green->nthreads); 562 #pragma omp parallel for schedule(static) 563 for(i = 0; i < (int)sz; i++) { 564 struct result* result = darray_result_data_get(&green->results) + i; 565 res_T tmp_res = RES_OK; 566 int t = omp_get_thread_num(); 567 568 if(res != RES_OK) continue; 569 570 if(settings[t] == NULL) { 571 /* description_count + 1 for AMBIENT */ 572 settings[t] = MEM_CALLOC(green->allocator, 1 + green->header.description_count, 573 sizeof(**settings)); 574 if(settings[t] == NULL) { 575 res = RES_MEM_ERR; 576 continue; 577 } 578 } 579 580 /* Apply settings */ 581 applied_settings_set_defaults(green, settings[t]); 582 tmp_res = parse_line(in_name, i, green, settings[t], &logs); 583 if(tmp_res != RES_OK) { 584 res = tmp_res; 585 continue; 586 } 587 /* Compute */ 588 green_compute_1(green, settings[t], &result->E, &result->STD); 589 } /* Implicit barrier */ 590 if(res != RES_OK) goto error; 591 592 exit: 593 do_log(in_name, green, &logs); 594 if(settings) { 595 FOR_EACH(n, 0, green->nthreads) MEM_RM(green->allocator, settings[n]); 596 MEM_RM(green->allocator, settings); 597 } 598 FOR_EACH(n, 0, darray_logs_size_get(&logs)) { 599 struct log_record* record = darray_logs_data_get(&logs)[n]; 600 log_record_release(record); 601 } 602 darray_logs_release(&logs); 603 return res; 604 error: 605 goto exit; 606 }