sdis_solve_probe_Xd.h (23478B)
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_c.h" 17 #include "sdis_device_c.h" 18 #include "sdis_estimator_c.h" 19 #include "sdis_estimator_buffer_c.h" 20 #include "sdis_log.h" 21 #include "sdis_green.h" 22 #include "sdis_misc.h" 23 #include "sdis_realisation.h" 24 #include "sdis_scene_c.h" 25 26 #include <rsys/clock_time.h> 27 #include <star/ssp.h> 28 #include <omp.h> 29 30 #include "sdis_Xd_begin.h" 31 32 /******************************************************************************* 33 * Helper function 34 ******************************************************************************/ 35 #ifndef SDIS_SOLVE_PROBE_XD_H 36 #define SDIS_SOLVE_PROBE_XD_H 37 38 static INLINE res_T 39 check_solve_probe_args(const struct sdis_solve_probe_args* args) 40 { 41 if(!args) return RES_BAD_ARG; 42 43 /* Check #realisations */ 44 if(!args->nrealisations || args->nrealisations > INT64_MAX) { 45 return RES_BAD_ARG; 46 } 47 48 /* Check time range */ 49 if(args->time_range[0] < 0 || args->time_range[1] < args->time_range[0]) { 50 return RES_BAD_ARG; 51 } 52 if(args->time_range[1] > DBL_MAX 53 && args->time_range[0] != args->time_range[1]) { 54 return RES_BAD_ARG; 55 } 56 57 /* Check picard order */ 58 if(args->picard_order < 1) { 59 return RES_BAD_ARG; 60 } 61 62 /* Check the RNG type */ 63 if(!args->rng_state && args->rng_type >= SSP_RNG_TYPES_COUNT__) { 64 return RES_BAD_ARG; 65 } 66 67 /* Check the diffusion algorithm */ 68 if((unsigned)args->diff_algo >= SDIS_DIFFUSION_ALGORITHMS_COUNT__) { 69 return RES_BAD_ARG; 70 } 71 72 return RES_OK; 73 } 74 75 static INLINE res_T 76 check_solve_probe_list_args 77 (struct sdis_device* dev, 78 const struct sdis_solve_probe_list_args* args) 79 { 80 size_t iprobe = 0; 81 82 if(!args) return RES_BAD_ARG; 83 84 /* Check the list of probes */ 85 if(!args->probes || !args->nprobes) { 86 return RES_BAD_ARG; 87 } 88 89 /* Check the RNG type */ 90 if(!args->rng_state && args->rng_type >= SSP_RNG_TYPES_COUNT__) { 91 return RES_BAD_ARG; 92 } 93 94 FOR_EACH(iprobe, 0, args->nprobes) { 95 const res_T res = check_solve_probe_args(args->probes+iprobe); 96 if(res != RES_OK) return res; 97 98 if(args->probes[iprobe].register_paths != SDIS_HEAT_PATH_NONE) { 99 log_warn(dev, 100 "Unable to save paths for probe %lu. " 101 "Saving path is not supported when solving multiple probes\n", 102 (unsigned long)iprobe); 103 } 104 } 105 106 return RES_OK; 107 } 108 109 #endif /* SDIS_SOLVE_PROBE_XD_H */ 110 111 static res_T 112 XD(solve_one_probe) 113 (struct sdis_scene* scn, 114 struct ssp_rng* rng, 115 const struct sdis_solve_probe_args* args, 116 struct accum* acc_temp, 117 struct accum* acc_time) 118 { 119 /* Enclosure in which the probe lies */ 120 unsigned enc_id = ENCLOSURE_ID_NULL; 121 122 size_t irealisation = 0; 123 res_T res = RES_OK; 124 ASSERT(scn && rng && check_solve_probe_args(args) == RES_OK); 125 ASSERT(acc_temp && acc_time); 126 127 *acc_temp = ACCUM_NULL; 128 *acc_time = ACCUM_NULL; 129 130 /* Retrieve the medium in which the submitted position lies */ 131 res = scene_get_enclosure_id(scn, args->position, &enc_id); 132 if(res != RES_OK) goto error; 133 134 FOR_EACH(irealisation, 0, args->nrealisations) { 135 struct probe_realisation_args realis_args = PROBE_REALISATION_ARGS_NULL; 136 double w = NaN; /* MC weight */ 137 double usec = 0; /* Time of a realisation */ 138 double time = 0; /* Sampled observation time */ 139 struct time t0, t1; /* Register the time spent solving a realisation */ 140 141 /* Begin time registration of the realisation */ 142 time_current(&t0); 143 144 /* Sample observation time */ 145 time = sample_time(rng, args->time_range); 146 147 /* Run a realisation */ 148 realis_args.rng = rng; 149 realis_args.enc_id = enc_id; 150 realis_args.time = time; 151 realis_args.picard_order = args->picard_order; 152 realis_args.irealisation = irealisation; 153 realis_args.diff_algo = args->diff_algo; 154 dX(set)(realis_args.position, args->position); 155 res = XD(probe_realisation)(scn, &realis_args, &w); 156 if(res != RES_OK && res != RES_BAD_OP) goto error; 157 158 switch(res) { 159 /* Reject the realisation */ 160 case RES_BAD_OP: 161 res = RES_OK; 162 break; 163 164 /* Update the accumulators */ 165 case RES_OK: 166 /* Stop time registration */ 167 time_sub(&t0, time_current(&t1), &t0); 168 usec = (double)time_val(&t0, TIME_NSEC) * 0.001; 169 170 /* Update MC weights */ 171 acc_temp->sum += w; 172 acc_temp->sum2 += w*w; 173 acc_temp->count += 1; 174 acc_time->sum += usec; 175 acc_time->sum2 += usec*usec; 176 acc_time->count += 1; 177 break; 178 179 default: FATAL("Unreachable code\n"); break; 180 } 181 } 182 183 exit: 184 return res; 185 error: 186 goto exit; 187 } 188 189 /******************************************************************************* 190 * Generic solve function 191 ******************************************************************************/ 192 static res_T 193 XD(solve_probe) 194 (struct sdis_scene* scn, 195 const struct sdis_solve_probe_args* args, 196 struct sdis_green_function** out_green, /* May be NULL <=> No green func */ 197 struct sdis_estimator** out_estimator) /* May be NULL <=> No estimator */ 198 { 199 /* Time registration */ 200 struct time time0, time1; 201 char buf[128]; /* Temporary buffer used to store formated time */ 202 203 /* Device variables */ 204 struct mem_allocator* allocator = NULL; 205 size_t nthreads = 0; 206 207 /* Stardis variables */ 208 struct sdis_estimator* estimator = NULL; 209 struct sdis_green_function* green = NULL; 210 struct sdis_green_function** per_thread_green = NULL; 211 212 /* Random Number generator */ 213 struct ssp_rng_proxy* rng_proxy = NULL; 214 struct ssp_rng** per_thread_rng = NULL; 215 216 /* Enclosure in which the probe lies */ 217 const struct enclosure* enc = NULL; 218 unsigned enc_id = ENCLOSURE_ID_NULL; 219 220 /* Miscellaneous */ 221 struct accum* per_thread_acc_temp = NULL; 222 struct accum* per_thread_acc_time = NULL; 223 size_t nrealisations = 0; 224 int64_t irealisation = 0; 225 int32_t* progress = NULL; /* Per process progress bar */ 226 int pcent_progress = 1; /* Percentage requiring progress update */ 227 int register_paths = SDIS_HEAT_PATH_NONE; 228 int is_master_process = 1; 229 ATOMIC nsolved_realisations = 0; 230 ATOMIC res = RES_OK; 231 232 if(!scn) { res = RES_BAD_ARG; goto error; } 233 if(!out_estimator && !out_green) { res = RES_BAD_ARG; goto error; } 234 res = check_solve_probe_args(args); 235 if(res != RES_OK) goto error; 236 res = XD(scene_check_dimensionality)(scn); 237 if(res != RES_OK) goto error; 238 239 if(out_green && args->picard_order != 1) { 240 log_err(scn->dev, "%s: the evaluation of the green function does not make " 241 "sense when dealing with the non-linearities of the system; i.e. picard " 242 "order must be set to 1 while it is currently set to %lu.\n", 243 FUNC_NAME, (unsigned long)args->picard_order); 244 res = RES_BAD_ARG; 245 goto error; 246 } 247 248 #ifdef SDIS_ENABLE_MPI 249 is_master_process = !scn->dev->use_mpi || scn->dev->mpi_rank == 0; 250 #endif 251 252 nthreads = scn->dev->nthreads; 253 allocator = scn->dev->allocator; 254 255 /* Update the progress bar every percent if escape sequences are allowed in 256 * log messages or only every 10 percent when only plain text is allowed. 257 * This reduces the number of lines of plain text printed */ 258 pcent_progress = scn->dev->no_escape_sequence ? 10 : 1; 259 260 /* Create the per thread RNGs */ 261 res = create_per_thread_rng 262 (scn->dev, args->rng_state, args->rng_type, &rng_proxy, &per_thread_rng); 263 if(res != RES_OK) goto error; 264 265 /* Allocate the per process progress status */ 266 res = alloc_process_progress(scn->dev, &progress); 267 if(res != RES_OK) goto error; 268 269 /* Create the per thread accumulators */ 270 per_thread_acc_temp = MEM_CALLOC(allocator, nthreads, sizeof(struct accum)); 271 per_thread_acc_time = MEM_CALLOC(allocator, nthreads, sizeof(struct accum)); 272 if(!per_thread_acc_temp) { res = RES_MEM_ERR; goto error; } 273 if(!per_thread_acc_time) { res = RES_MEM_ERR; goto error; } 274 275 /* Retrieve the enclosure in which the submitted position lies */ 276 res = scene_get_enclosure_id(scn, args->position, &enc_id); 277 if(res != RES_OK) goto error; 278 279 /* Check that the enclosure does not contain multiple materials */ 280 enc = scene_get_enclosure(scn, enc_id); 281 if(enc->medium_id == MEDIUM_ID_MULTI) { 282 log_err(scn->dev, "%s: probe is in an enclosure with several media " 283 "-- pos=("FORMAT_VECX")\n", FUNC_NAME, SPLITX(args->position)); 284 res = RES_BAD_OP; 285 goto error; 286 } 287 288 /* Create the per thread green function */ 289 if(out_green) { 290 res = create_per_thread_green_function 291 (scn, args->signature, &per_thread_green); 292 if(res != RES_OK) goto error; 293 } 294 295 /* Create the estimator on the master process only. No estimator is needed 296 * for non master process */ 297 if(out_estimator && is_master_process) { 298 res = estimator_create(scn->dev, SDIS_ESTIMATOR_TEMPERATURE, &estimator); 299 if(res != RES_OK) goto error; 300 } 301 302 /* Synchronise the processes */ 303 process_barrier(scn->dev); 304 305 #define PROGRESS_MSG "Solving probe temperature: " 306 print_progress(scn->dev, progress, PROGRESS_MSG); 307 308 /* Begin time registration of the computation */ 309 time_current(&time0); 310 311 /* Here we go! Launch the Monte Carlo estimation */ 312 nrealisations = compute_process_realisations_count(scn->dev, args->nrealisations); 313 register_paths = out_estimator && is_master_process 314 ? args->register_paths : SDIS_HEAT_PATH_NONE; 315 omp_set_num_threads((int)scn->dev->nthreads); 316 #pragma omp parallel for schedule(static) 317 for(irealisation = 0; irealisation < (int64_t)nrealisations; ++irealisation) { 318 struct probe_realisation_args realis_args = PROBE_REALISATION_ARGS_NULL; 319 struct time t0, t1; 320 const int ithread = omp_get_thread_num(); 321 struct ssp_rng* rng = per_thread_rng[ithread]; 322 struct accum* acc_temp = &per_thread_acc_temp[ithread]; 323 struct accum* acc_time = &per_thread_acc_time[ithread]; 324 struct green_path_handle* pgreen_path = NULL; 325 struct green_path_handle green_path = GREEN_PATH_HANDLE_NULL; 326 struct sdis_heat_path* pheat_path = NULL; 327 struct sdis_heat_path heat_path; 328 double w = NaN; 329 double time; 330 size_t n; 331 int pcent; 332 res_T res_local = RES_OK; 333 res_T res_simul = RES_OK; 334 335 if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occurred */ 336 337 /* Begin time registration of the realisation */ 338 time_current(&t0); 339 340 time = sample_time(rng, args->time_range); 341 342 if(out_green) { 343 res_local = green_function_create_path 344 (per_thread_green[ithread], &green_path); 345 if(res_local != RES_OK) { ATOMIC_SET(&res, res_local); goto error_it; } 346 pgreen_path = &green_path; 347 } 348 349 if(register_paths) { 350 heat_path_init(scn->dev->allocator, &heat_path); 351 pheat_path = &heat_path; 352 } 353 354 /* Invoke the probe realisation */ 355 realis_args.rng = rng; 356 realis_args.enc_id = enc_id; 357 realis_args.time = time; 358 realis_args.picard_order = args->picard_order; 359 realis_args.green_path = pgreen_path; 360 realis_args.heat_path = pheat_path; 361 realis_args.irealisation = (size_t)irealisation; 362 realis_args.diff_algo = args->diff_algo; 363 dX(set)(realis_args.position, args->position); 364 res_simul = XD(probe_realisation)(scn, &realis_args, &w); 365 366 /* Handle fatal error */ 367 if(res_simul != RES_OK && res_simul != RES_BAD_OP) { 368 ATOMIC_SET(&res, res_simul); 369 goto error_it; 370 } 371 372 if(pheat_path) { 373 pheat_path->status = res_simul == RES_OK 374 ? SDIS_HEAT_PATH_SUCCESS 375 : SDIS_HEAT_PATH_FAILURE; 376 377 /* Check if the path must be saved regarding the register_paths mask */ 378 if(!(register_paths & (int)pheat_path->status)) { 379 heat_path_release(pheat_path); 380 pheat_path = NULL; 381 } else { /* Register the sampled path */ 382 res_local = estimator_add_and_release_heat_path(estimator, pheat_path); 383 if(res_local != RES_OK) { 384 ATOMIC_SET(&res, res_local); 385 goto error_it; 386 } 387 pheat_path = NULL; 388 } 389 } 390 391 /* Stop time registration */ 392 time_sub(&t0, time_current(&t1), &t0); 393 394 /* Update accumulators */ 395 if(res_simul == RES_OK) { 396 const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; 397 acc_temp->sum += w; acc_temp->sum2 += w*w; ++acc_temp->count; 398 acc_time->sum += usec; acc_time->sum2 += usec*usec; ++acc_time->count; 399 } 400 401 /* Update the progress status */ 402 n = (size_t)ATOMIC_INCR(&nsolved_realisations); 403 pcent = (int)((double)n * 100.0 / (double)nrealisations + 0.5/*round*/); 404 405 #pragma omp critical 406 if(pcent/pcent_progress > progress[0]/pcent_progress) { 407 progress[0] = pcent; 408 print_progress_update(scn->dev, progress, PROGRESS_MSG); 409 } 410 411 exit_it: 412 if(pheat_path) heat_path_release(pheat_path); 413 continue; 414 error_it: 415 goto exit_it; 416 } 417 418 /* Synchronise processes */ 419 process_barrier(scn->dev); 420 421 res = gather_res_T(scn->dev, (res_T)res); 422 if(res != RES_OK) goto error; 423 424 print_progress_completion(scn->dev, progress, PROGRESS_MSG); 425 #undef PROGRESS_MSG 426 427 /* Report computation time */ 428 time_sub(&time0, time_current(&time1), &time0); 429 time_dump(&time0, TIME_ALL, NULL, buf, sizeof(buf)); 430 log_info(scn->dev, "Probe temperature solved in %s.\n", buf); 431 432 /* Gather the RNG proxy sequence IDs and ensure that the RNG proxy state of 433 * the master process is greater than the RNG proxy state of all other 434 * processes */ 435 res = gather_rng_proxy_sequence_id(scn->dev, rng_proxy); 436 if(res != RES_OK) goto error; 437 438 /* Setup the estimated values */ 439 if(out_estimator) { 440 struct accum acc_temp; 441 struct accum acc_time; 442 443 time_current(&time0); 444 445 #define GATHER_ACCUMS(Msg, Acc) { \ 446 res = gather_accumulators(scn->dev, Msg, per_thread_##Acc, &Acc); \ 447 if(res != RES_OK) goto error; \ 448 } (void)0 449 GATHER_ACCUMS(MPI_SDIS_MSG_ACCUM_TEMP, acc_temp); 450 GATHER_ACCUMS(MPI_SDIS_MSG_ACCUM_TIME, acc_time); 451 #undef GATHER_ACCUMS 452 453 time_sub(&time0, time_current(&time1), &time0); 454 time_dump(&time0, TIME_ALL, NULL, buf, sizeof(buf)); 455 log_info(scn->dev, "Accumulators gathered in %s.\n", buf); 456 457 /* Return an estimator only on master process */ 458 if(is_master_process) { 459 ASSERT(acc_temp.count == acc_time.count); 460 estimator_setup_realisations_count(estimator, args->nrealisations, acc_temp.count); 461 estimator_setup_temperature(estimator, acc_temp.sum, acc_temp.sum2); 462 estimator_setup_realisation_time(estimator, acc_time.sum, acc_time.sum2); 463 res = estimator_save_rng_state(estimator, rng_proxy); 464 if(res != RES_OK) goto error; 465 } 466 } 467 468 /* Setup the green function */ 469 if(out_green) { 470 time_current(&time0); 471 472 res = gather_green_functions 473 (scn, rng_proxy, per_thread_green, per_thread_acc_time, &green); 474 if(res != RES_OK) goto error; 475 476 time_sub(&time0, time_current(&time1), &time0); 477 time_dump(&time0, TIME_ALL, NULL, buf, sizeof(buf)); 478 log_info(scn->dev, "Green functions gathered in %s.\n", buf); 479 480 /* Return a green function only on master process */ 481 if(!is_master_process) { 482 SDIS(green_function_ref_put(green)); 483 green = NULL; 484 } 485 } 486 487 exit: 488 if(per_thread_rng) release_per_thread_rng(scn->dev, per_thread_rng); 489 if(per_thread_green) release_per_thread_green_function(scn, per_thread_green); 490 if(progress) free_process_progress(scn->dev, progress); 491 if(per_thread_acc_temp) MEM_RM(scn->dev->allocator, per_thread_acc_temp); 492 if(per_thread_acc_time) MEM_RM(scn->dev->allocator, per_thread_acc_time); 493 if(rng_proxy) SSP(rng_proxy_ref_put(rng_proxy)); 494 if(out_green) *out_green = green; 495 if(out_estimator) *out_estimator = estimator; 496 return (res_T)res; 497 error: 498 if(estimator) { SDIS(estimator_ref_put(estimator)); estimator = NULL; } 499 if(green) { SDIS(green_function_ref_put(green)); green = NULL; } 500 goto exit; 501 } 502 503 static res_T 504 XD(solve_probe_list) 505 (struct sdis_scene* scn, 506 const struct sdis_solve_probe_list_args* args, 507 struct sdis_estimator_buffer** out_estim_buf) 508 { 509 /* Time registration */ 510 struct time time0, time1; 511 char buf[128]; /* Temporary buffer used to store formated time */ 512 513 /* Device variables */ 514 struct mem_allocator* allocator = NULL; 515 516 /* Stardis variables */ 517 struct sdis_estimator_buffer* estim_buf = NULL; 518 519 /* Random Number generator */ 520 struct ssp_rng_proxy* rng_proxy = NULL; 521 struct ssp_rng** per_thread_rng = NULL; 522 523 /* Probe variables */ 524 size_t process_probes[2] = {0, 0}; /* Probes range managed by the process */ 525 size_t process_nprobes = 0; /* Number of probes managed by the process */ 526 size_t nprobes = 0; 527 struct accum* per_probe_acc_temp = NULL; 528 struct accum* per_probe_acc_time = NULL; 529 530 /* Miscellaneous */ 531 int32_t* progress = NULL; /* Per process progress bar */ 532 int pcent_progress = 1; /* Percentage requiring progress update */ 533 int is_master_process = 1; 534 int64_t i = 0; 535 ATOMIC nsolved_probes = 0; 536 ATOMIC res = RES_OK; 537 538 /* Check input arguments */ 539 if(!scn || !out_estim_buf) { res = RES_BAD_ARG; goto error; } 540 res = check_solve_probe_list_args(scn->dev, args); 541 if(res != RES_OK) goto error; 542 res = XD(scene_check_dimensionality)(scn); 543 if(res != RES_OK) goto error; 544 545 #ifdef SDIS_ENABLE_MPI 546 is_master_process = !scn->dev->use_mpi || scn->dev->mpi_rank == 0; 547 #endif 548 549 allocator = scn->dev->allocator; 550 551 /* Update the progress bar every percent if escape sequences are allowed in 552 * log messages or only every 10 percent when only plain text is allowed. 553 * This reduces the number of lines of plain text printed */ 554 pcent_progress = scn->dev->no_escape_sequence ? 10 : 1; 555 556 /* Create the per thread RNGs */ 557 res = create_per_thread_rng 558 (scn->dev, args->rng_state, args->rng_type, &rng_proxy, &per_thread_rng); 559 if(res != RES_OK) goto error; 560 561 /* Allocate the per process progress status */ 562 res = alloc_process_progress(scn->dev, &progress); 563 if(res != RES_OK) goto error; 564 565 /* Synchronise the processes */ 566 process_barrier(scn->dev); 567 568 /* Define the range of probes to manage in this process */ 569 process_nprobes = compute_process_index_range 570 (scn->dev, args->nprobes, process_probes); 571 572 #define PROGRESS_MSG "Solving probes: " 573 print_progress(scn->dev, progress, PROGRESS_MSG); 574 575 /* If there is no work to be done on this process (i.e. no probe to 576 * calculate), simply print its completion and go straight to the 577 * synchronization barrier.*/ 578 if(process_nprobes == 0) { 579 progress[0] = 100; 580 print_progress_update(scn->dev, progress, PROGRESS_MSG); 581 goto post_sync; 582 } 583 584 /* Allocate the list of accumulators per probe. On the master process, 585 * allocate a complete list in which the accumulators of all processes will be 586 * stored. */ 587 nprobes = is_master_process ? args->nprobes : process_nprobes; 588 per_probe_acc_temp = MEM_CALLOC(allocator, nprobes, sizeof(*per_probe_acc_temp)); 589 per_probe_acc_time = MEM_CALLOC(allocator, nprobes, sizeof(*per_probe_acc_time)); 590 if(!per_probe_acc_temp || !per_probe_acc_time) { 591 log_err(scn->dev, "Unable to allocate the list of accumulators per probe.\n"); 592 res = RES_MEM_ERR; 593 goto error; 594 } 595 596 /* Begin time registration of the computation */ 597 time_current(&time0); 598 599 /* Here we go! Calculation of probe list */ 600 omp_set_num_threads((int)scn->dev->nthreads); 601 #pragma omp parallel for schedule(static) 602 for(i = 0; i < (int64_t)process_nprobes; ++i) { 603 /* Thread */ 604 const int ithread = omp_get_thread_num(); /* Thread ID */ 605 struct ssp_rng* rng = per_thread_rng[ithread]; /* RNG of the thread */ 606 607 /* Probe */ 608 struct accum* probe_acc_temp = NULL; 609 struct accum* probe_acc_time = NULL; 610 const struct sdis_solve_probe_args* probe_args = NULL; /* Solve args */ 611 const size_t iprobe = process_probes[0] + (size_t)i; /* Probe ID */ 612 613 /* Misc */ 614 size_t n = 0; /* Number of solved probes */ 615 int pcent = 0; /* Current progress */ 616 res_T res_local = RES_OK; 617 618 if(ATOMIC_GET(&res) != RES_OK) continue; /* An error occured */ 619 620 /* Retrieve the probe arguments */ 621 probe_args = &args->probes[iprobe]; 622 623 /* Retrieve the probe accumulators */ 624 probe_acc_temp = &per_probe_acc_temp[i]; 625 probe_acc_time = &per_probe_acc_time[i]; 626 627 res_local = XD(solve_one_probe) 628 (scn, rng, probe_args, probe_acc_temp, probe_acc_time); 629 if(res_local != RES_OK) { 630 ATOMIC_SET(&res, res_local); 631 continue; 632 } 633 634 /* Update progress */ 635 n = (size_t)ATOMIC_INCR(&nsolved_probes); 636 pcent = (int)((double)n * 100.0 / (double)process_nprobes + 0.5/*round*/); 637 638 #pragma omp critical 639 if(pcent/pcent_progress > progress[0]/pcent_progress) { 640 progress[0] = pcent; 641 print_progress_update(scn->dev, progress, PROGRESS_MSG); 642 } 643 } 644 645 post_sync: 646 /* Synchronise processes */ 647 process_barrier(scn->dev); 648 649 res = gather_res_T(scn->dev, (res_T)res); 650 if(res != RES_OK) goto error; 651 652 print_progress_completion(scn->dev, progress, PROGRESS_MSG); 653 #undef PROGRESS_MSG 654 655 /* Report computation time */ 656 time_sub(&time0, time_current(&time1), &time0); 657 time_dump(&time0, TIME_ALL, NULL, buf, sizeof(buf)); 658 log_info(scn->dev, "%lu probes solved in %s.\n", 659 (unsigned long)args->nprobes, buf); 660 661 /* Gather the RNG proxy sequence IDs and ensure that the RNG proxy state of 662 * the master process is greater than the RNG proxy state of all other 663 * processes */ 664 res = gather_rng_proxy_sequence_id(scn->dev, rng_proxy); 665 if(res != RES_OK) goto error; 666 667 time_current(&time0); 668 669 /* Gather the list of accumulators */ 670 #define GATHER_ACCUMS_LIST(Msg, Acc) { \ 671 res = gather_accumulators_list \ 672 (scn->dev, Msg, args->nprobes, process_probes, per_probe_acc_##Acc); \ 673 if(res != RES_OK) goto error; \ 674 } (void)0 675 GATHER_ACCUMS_LIST(MPI_SDIS_MSG_ACCUM_TEMP, temp); 676 GATHER_ACCUMS_LIST(MPI_SDIS_MSG_ACCUM_TIME, time); 677 #undef GATHER_ACCUMS_LIST 678 679 time_sub(&time0, time_current(&time1), &time0); 680 time_dump(&time0, TIME_ALL, NULL, buf, sizeof(buf)); 681 log_info(scn->dev, "Probes accumulator gathered in %s.\n", buf); 682 683 if(is_master_process) { 684 res = estimator_buffer_create_from_observable_list_probe 685 (scn->dev, rng_proxy, args->probes, per_probe_acc_temp, 686 per_probe_acc_time, args->nprobes, &estim_buf); 687 if(res != RES_OK) goto error; 688 } 689 690 exit: 691 if(per_thread_rng) release_per_thread_rng(scn->dev, per_thread_rng); 692 if(rng_proxy) SSP(rng_proxy_ref_put(rng_proxy)); 693 if(per_probe_acc_temp) MEM_RM(allocator, per_probe_acc_temp); 694 if(per_probe_acc_time) MEM_RM(allocator, per_probe_acc_time); 695 if(progress) free_process_progress(scn->dev, progress); 696 if(out_estim_buf) *out_estim_buf = estim_buf; 697 return (res_T)res; 698 error: 699 if(estim_buf) { 700 SDIS(estimator_buffer_ref_put(estim_buf)); 701 estim_buf = NULL; 702 } 703 goto exit; 704 } 705 706 #include "sdis_Xd_end.h"