solstice_solve.c (22244B)
1 /* Copyright (C) 2018-2026 |Meso|Star> (contact@meso-star.com) 2 * Copyright (C) 2016-2018 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 #include "solstice.h" 18 #include "solstice_c.h" 19 #include "parser/solparser_sun.h" 20 #include <solstice/ssol.h> 21 #include <star/ssp.h> 22 23 /* How many percent of random walk realisations may fail before solve() stops 24 * in a standard solve() invocation. 25 * It is not used when solve() is invoked in order to dump paths. */ 26 #define MAX_PERCENT_FAILURES 0.01 27 28 /******************************************************************************* 29 * Helper function 30 ******************************************************************************/ 31 static void 32 write_mc_global(struct solstice* solstice, struct ssol_estimator* estimator) 33 { 34 struct ssol_mc_global mc_global; 35 struct htable_primary_iterator p_it, p_end; 36 const struct solparser_sun* solparser_sun = NULL; 37 size_t nexperiments, nfailed, nprimary; 38 size_t ircv; 39 double area, potential, irradiance_factor; 40 ASSERT(solstice && estimator); 41 42 /* get global information */ 43 SSOL(estimator_get_mc_global(estimator, &mc_global)); 44 SSOL(estimator_get_realisation_count(estimator, &nexperiments)); 45 SSOL(estimator_get_sampled_count(estimator, &nprimary)); 46 SSOL(estimator_get_failed_count(estimator, &nfailed)); 47 SSOL(estimator_get_sampled_area(estimator, &area)); 48 solparser_sun = solparser_get_sun(solstice->parser); 49 potential = solparser_sun->dni * area; 50 irradiance_factor = 1 / potential; 51 52 /* Counts */ 53 fprintf(solstice->output, "%d %lu %lu %lu %lu\n", 54 7, /* #global results count */ 55 (unsigned long)darray_receiver_size_get(&solstice->rcvs_list), 56 (unsigned long)nprimary, 57 (unsigned long)nexperiments, 58 (unsigned long)nfailed); 59 60 /* Global data */ 61 #define PRINT_MC_GLOBAL(Name) \ 62 fprintf(solstice->output, "%g %g\n", mc_global.Name.E, mc_global.Name.SE) 63 fprintf(solstice->output, "%g %g\n", potential, 0.); 64 PRINT_MC_GLOBAL(absorbed_by_receivers); 65 PRINT_MC_GLOBAL(cos_factor); 66 PRINT_MC_GLOBAL(shadowed); 67 PRINT_MC_GLOBAL(missing); 68 PRINT_MC_GLOBAL(other_absorbed); 69 PRINT_MC_GLOBAL(extinguished_by_atmosphere); 70 #undef PRINT_MC_GLOBAL 71 72 /* Receivers' data */ 73 FOR_EACH(ircv, 0, darray_receiver_size_get(&solstice->rcvs_list)) { 74 struct solstice_receiver* rcv = darray_receiver_data_get 75 (&solstice->rcvs_list) + ircv; 76 struct ssol_instance* inst = rcv->node->instance; 77 struct ssol_mc_receiver front = MC_RCV_NONE__; 78 struct ssol_mc_receiver back = MC_RCV_NONE__; 79 double f_eff_E = -1, f_eff_SE = -1; /* Front efficiency */ 80 double b_eff_E = -1, b_eff_SE = -1; /* Back efficiency */ 81 uint32_t id; 82 83 switch(rcv->side) { 84 case SRCVL_FRONT: 85 SSOL(estimator_get_mc_receiver(estimator, inst, SSOL_FRONT, &front)); 86 f_eff_E = front.absorbed_flux.E * irradiance_factor; 87 f_eff_SE = front.absorbed_flux.SE * irradiance_factor; 88 break; 89 case SRCVL_BACK: 90 SSOL(estimator_get_mc_receiver(estimator, inst, SSOL_BACK, &back)); 91 b_eff_E = back.absorbed_flux.E * irradiance_factor; 92 b_eff_SE = back.absorbed_flux.SE * irradiance_factor; 93 break; 94 case SRCVL_FRONT_AND_BACK: 95 SSOL(estimator_get_mc_receiver(estimator, inst, SSOL_FRONT, &front)); 96 SSOL(estimator_get_mc_receiver(estimator, inst, SSOL_BACK, &back)); 97 f_eff_E = front.absorbed_flux.E * irradiance_factor; 98 f_eff_SE = front.absorbed_flux.SE * irradiance_factor; 99 b_eff_E = back.absorbed_flux.E * irradiance_factor; 100 b_eff_SE = back.absorbed_flux.SE * irradiance_factor; 101 break; 102 default: FATAL("Unreachable code.\n"); break; 103 } 104 SSOL(instance_get_id(inst, &id)); 105 SSOL(instance_get_area(inst, &area)); 106 fprintf(solstice->output, 107 "%s %u %g " 108 "%g %g %g %g %g %g %g %g %g %g %g %g " 109 "%g %g %g %g %g %g %g %g %g %g " 110 "%g %g %g %g %g %g %g %g %g %g %g %g " 111 "%g %g %g %g %g %g %g %g %g %g\n", 112 str_cget(&rcv->name), (unsigned)id, area, 113 front.incoming_flux.E, front.incoming_flux.SE, 114 front.incoming_if_no_field_loss.E, front.incoming_if_no_field_loss.SE, 115 front.incoming_if_no_atm_loss.E, front.incoming_if_no_atm_loss.SE, 116 front.incoming_lost_in_field.E, front.incoming_lost_in_field.SE, 117 front.incoming_lost_in_atmosphere.E, front.incoming_lost_in_atmosphere.SE, 118 front.absorbed_flux.E, front.absorbed_flux.SE, 119 front.absorbed_if_no_field_loss.E, front.absorbed_if_no_field_loss.SE, 120 front.absorbed_if_no_atm_loss.E, front.absorbed_if_no_atm_loss.SE, 121 front.absorbed_lost_in_field.E, front.absorbed_lost_in_field.SE, 122 front.absorbed_lost_in_atmosphere.E, front.absorbed_lost_in_atmosphere.SE, 123 f_eff_E, f_eff_SE, 124 back.incoming_flux.E, back.incoming_flux.SE, 125 back.incoming_if_no_field_loss.E, back.incoming_if_no_field_loss.SE, 126 back.incoming_if_no_atm_loss.E, back.incoming_if_no_atm_loss.SE, 127 back.incoming_lost_in_field.E, back.incoming_lost_in_field.SE, 128 back.incoming_lost_in_atmosphere.E, back.incoming_lost_in_atmosphere.SE, 129 back.absorbed_flux.E, back.absorbed_flux.SE, 130 back.absorbed_if_no_field_loss.E, back.absorbed_if_no_field_loss.SE, 131 back.absorbed_if_no_atm_loss.E, back.absorbed_if_no_atm_loss.SE, 132 back.absorbed_lost_in_field.E, back.absorbed_lost_in_field.SE, 133 back.absorbed_lost_in_atmosphere.E, back.absorbed_lost_in_atmosphere.SE, 134 b_eff_E, b_eff_SE); 135 } 136 137 /* Primary-instances' data */ 138 htable_primary_begin(&solstice->primaries, &p_it); 139 htable_primary_end(&solstice->primaries, &p_end); 140 while(!htable_primary_iterator_eq(&p_it, &p_end)) { 141 const struct str* name = htable_primary_iterator_key_get(&p_it); 142 struct solstice_primary* prim = htable_primary_iterator_data_get(&p_it); 143 struct ssol_mc_sampled sampled; 144 uint32_t id; 145 146 htable_primary_iterator_next(&p_it); 147 SSOL(estimator_get_mc_sampled(estimator, prim->node->instance, &sampled)); 148 SSOL(instance_get_id(prim->node->instance, &id)); 149 SSOL(instance_get_area(prim->node->instance, &area)); 150 fprintf(solstice->output, 151 "%s %u %g %lu " 152 "%g %g %g %g\n", 153 str_cget(name), (unsigned) id, area, (unsigned long)sampled.nb_samples, 154 sampled.cos_factor.E, sampled.cos_factor.SE, 155 sampled.shadowed.E, sampled.shadowed.SE 156 ); 157 } 158 159 /* ReceiverXprimarys' data */ 160 FOR_EACH(ircv, 0, darray_receiver_size_get(&solstice->rcvs_list)) { 161 struct solstice_receiver* rcv = darray_receiver_data_get 162 (&solstice->rcvs_list) + ircv; 163 struct ssol_instance* rcv_inst = rcv->node->instance; 164 uint32_t rcv_id, prim_id; 165 166 SSOL(instance_get_id(rcv_inst, &rcv_id)); 167 htable_primary_begin(&solstice->primaries, &p_it); 168 htable_primary_end(&solstice->primaries, &p_end); 169 while(!htable_primary_iterator_eq(&p_it, &p_end)) { 170 struct solstice_primary* prim = htable_primary_iterator_data_get(&p_it); 171 struct ssol_instance* prim_inst = prim->node->instance; 172 struct ssol_mc_receiver front = MC_RCV_NONE__; 173 struct ssol_mc_receiver back = MC_RCV_NONE__; 174 175 SSOL(instance_get_id(prim_inst, &prim_id)); 176 switch (rcv->side) { 177 case SRCVL_FRONT: 178 SSOL(estimator_get_mc_sampled_x_receiver 179 (estimator, prim_inst, rcv_inst, SSOL_FRONT, &front)); 180 break; 181 case SRCVL_BACK: 182 SSOL(estimator_get_mc_sampled_x_receiver 183 (estimator, prim_inst, rcv_inst, SSOL_BACK, &back)); 184 break; 185 case SRCVL_FRONT_AND_BACK: 186 SSOL(estimator_get_mc_sampled_x_receiver 187 (estimator, prim_inst, rcv_inst, SSOL_FRONT, &front)); 188 SSOL(estimator_get_mc_sampled_x_receiver 189 (estimator, prim_inst, rcv_inst, SSOL_BACK, &back)); 190 break; 191 default: FATAL("Unreachable code.\n"); break; 192 } 193 fprintf(solstice->output, 194 "%u %u " 195 "%g %g %g %g %g %g %g %g %g %g %g %g " 196 "%g %g %g %g %g %g %g %g " 197 "%g %g %g %g %g %g %g %g %g %g %g %g " 198 "%g %g %g %g %g %g %g %g\n", 199 (unsigned)rcv_id, (unsigned) prim_id, 200 front.incoming_flux.E, front.incoming_flux.SE, 201 front.incoming_if_no_field_loss.E, front.incoming_if_no_field_loss.SE, 202 front.incoming_if_no_atm_loss.E, front.incoming_if_no_atm_loss.SE, 203 front.incoming_lost_in_field.E, front.incoming_lost_in_field.SE, 204 front.incoming_lost_in_atmosphere.E, front.incoming_lost_in_atmosphere.SE, 205 front.absorbed_flux.E, front.absorbed_flux.SE, 206 front.absorbed_if_no_field_loss.E, front.absorbed_if_no_field_loss.SE, 207 front.absorbed_if_no_atm_loss.E, front.absorbed_if_no_atm_loss.SE, 208 front.absorbed_lost_in_field.E, front.absorbed_lost_in_field.SE, 209 front.absorbed_lost_in_atmosphere.E, front.absorbed_lost_in_atmosphere.SE, 210 back.incoming_flux.E, back.incoming_flux.SE, 211 back.incoming_if_no_field_loss.E, back.incoming_if_no_field_loss.SE, 212 back.incoming_if_no_atm_loss.E, back.incoming_if_no_atm_loss.SE, 213 back.incoming_lost_in_field.E, back.incoming_lost_in_field.SE, 214 back.incoming_lost_in_atmosphere.E, back.incoming_lost_in_atmosphere.SE, 215 back.absorbed_flux.E, back.absorbed_flux.SE, 216 back.absorbed_if_no_field_loss.E, back.absorbed_if_no_field_loss.SE, 217 back.absorbed_if_no_atm_loss.E, back.absorbed_if_no_atm_loss.SE, 218 back.absorbed_lost_in_field.E, back.absorbed_lost_in_field.SE, 219 back.absorbed_lost_in_atmosphere.E, back.absorbed_lost_in_atmosphere.SE); 220 htable_primary_iterator_next(&p_it); 221 } 222 } 223 } 224 225 static void 226 dump_instantiated_shaded_shape_vertices 227 (struct solstice* solstice, 228 const struct ssol_instantiated_shaded_shape* inst_sshape) 229 { 230 unsigned ivert, nverts; 231 ASSERT(solstice && inst_sshape); 232 233 SSOL(shape_get_vertices_count(inst_sshape->shape, &nverts)); 234 FOR_EACH(ivert, 0, nverts) { 235 double pos[3]; 236 SSOL(instantiated_shaded_shape_get_vertex_attrib 237 (inst_sshape, ivert, SSOL_POSITION, pos)); 238 fprintf(solstice->output, "%f %f %f\n", SPLIT3(pos)); 239 } 240 } 241 242 static void 243 dump_shape_triangle_indices 244 (struct solstice* solstice, 245 const struct ssol_shape* shape, 246 const size_t offset) 247 { 248 unsigned itri, ntris; 249 ASSERT(solstice && shape); 250 251 SSOL(shape_get_triangles_count(shape, &ntris)); 252 FOR_EACH(itri, 0, ntris) { 253 unsigned ids[3]; 254 SSOL(shape_get_triangle_indices(shape, itri, ids)); 255 fprintf(solstice->output, "3 %lu %lu %lu\n", 256 (unsigned long)(ids[0] + offset), 257 (unsigned long)(ids[1] + offset), 258 (unsigned long)(ids[2] + offset)); 259 } 260 } 261 262 static void 263 dump_mc_shape_in 264 (struct solstice* solstice, 265 struct ssol_shape* shape, 266 struct ssol_mc_shape* mc_shape) 267 { 268 unsigned itri, ntris; 269 ASSERT(solstice && shape && mc_shape); 270 271 SSOL(shape_get_triangles_count(shape, &ntris)); 272 FOR_EACH(itri, 0, ntris) { 273 struct ssol_mc_primitive mc_prim; 274 SSOL(mc_shape_get_mc_primitive(mc_shape, itri, &mc_prim)); 275 fprintf(solstice->output, "%g %g\n", 276 mc_prim.incoming_flux.E, 277 mc_prim.incoming_flux.SE); 278 } 279 } 280 281 static void 282 dump_mc_shape_abs 283 (struct solstice* solstice, 284 struct ssol_shape* shape, 285 struct ssol_mc_shape* mc_shape) 286 { 287 unsigned itri, ntris; 288 ASSERT(solstice && shape && mc_shape); 289 290 SSOL(shape_get_triangles_count(shape, &ntris)); 291 FOR_EACH(itri, 0, ntris) { 292 struct ssol_mc_primitive mc_prim; 293 SSOL(mc_shape_get_mc_primitive(mc_shape, itri, &mc_prim)); 294 fprintf(solstice->output, "%g %g\n", 295 mc_prim.absorbed_flux.E, 296 mc_prim.absorbed_flux.SE); 297 } 298 } 299 300 static void 301 dump_per_primitive_mc_estimations 302 (struct solstice* solstice, 303 struct ssol_estimator* estimator, 304 struct ssol_instance* inst, 305 const enum ssol_side_flag side, 306 const enum srcvl_pp_output output) 307 { 308 size_t ishape, nshapes; 309 struct ssol_mc_receiver mc_rcv; 310 const char* name; 311 const char* flux; 312 ASSERT(solstice && estimator && inst); 313 314 SSOL(estimator_get_mc_receiver(estimator, inst, side, &mc_rcv)); 315 316 switch(side) { 317 case SSOL_FRONT: name = "Front_faces"; break; 318 case SSOL_BACK: name = "Back_faces"; break; 319 default: FATAL("Unreachable code.\n"); break; 320 } 321 switch (output) { 322 case SRCVL_PP_INCOMING: flux = "Incoming_flux"; break; 323 case SRCVL_PP_ABSORBED: flux = "Absorbed_flux"; break; 324 default: FATAL("Unreachable code.\n"); break; 325 } 326 327 fprintf(solstice->output, "SCALARS %s_%s double 2\n", name, flux); 328 fprintf(solstice->output, "LOOKUP_TABLE default\n"); 329 330 SSOL(instance_get_shaded_shapes_count(inst, &nshapes)); 331 FOR_EACH(ishape, 0, nshapes) { 332 struct ssol_instantiated_shaded_shape inst_sshape; 333 struct ssol_mc_shape mc_shape; 334 SSOL(instance_get_shaded_shape(inst, ishape, &inst_sshape)); 335 SSOL(mc_receiver_get_mc_shape(&mc_rcv, inst_sshape.shape, &mc_shape)); 336 switch (output) { 337 case SRCVL_PP_INCOMING: 338 dump_mc_shape_in(solstice, inst_sshape.shape, &mc_shape); 339 break; 340 case SRCVL_PP_ABSORBED: 341 dump_mc_shape_abs(solstice, inst_sshape.shape, &mc_shape); 342 break; 343 default: FATAL("Unreachable code.\n"); break; 344 } 345 } 346 } 347 348 static void 349 write_per_receiver_mc_primitive 350 (struct solstice* solstice, struct ssol_estimator* estimator) 351 { 352 size_t ircv; 353 ASSERT(solstice && estimator); 354 355 FOR_EACH(ircv, 0, darray_receiver_size_get(&solstice->rcvs_list)) { 356 struct ssol_instantiated_shaded_shape inst_sshape; 357 struct solstice_receiver* rcv = darray_receiver_data_get 358 (&solstice->rcvs_list) + ircv; 359 struct ssol_instance* inst = rcv->node->instance; 360 size_t ishape, nshapes; 361 size_t nverts, ntris; 362 size_t offset; 363 int mask, prim; 364 365 ASSERT(rcv->side); 366 SSOL(instance_is_receiver(inst, &mask, &prim)); 367 CHK(mask != 0); 368 if(!prim) continue; 369 370 SSOL(instance_get_shaded_shapes_count(inst, &nshapes)); 371 372 /* Write the header */ 373 fprintf(solstice->output, "# vtk DataFile Version 2.0\n"); 374 fprintf(solstice->output, "%s\n", str_cget(&rcv->name)); 375 fprintf(solstice->output, "ASCII\n"); 376 fprintf(solstice->output, "DATASET POLYDATA\n"); 377 378 /* Compute the overall number of vertices & triangles of the receiver */ 379 nverts = ntris = 0; 380 FOR_EACH(ishape, 0, nshapes) { 381 unsigned shape_nverts, shape_ntris; 382 SSOL(instance_get_shaded_shape(inst, ishape, &inst_sshape)); 383 SSOL(shape_get_vertices_count(inst_sshape.shape, &shape_nverts)); 384 SSOL(shape_get_triangles_count(inst_sshape.shape, &shape_ntris)); 385 nverts += shape_nverts; 386 ntris += shape_ntris; 387 } 388 389 /* Write the positions of the receiver shaded shapes */ 390 fprintf(solstice->output, "POINTS %lu float\n", (unsigned long)nverts); 391 FOR_EACH(ishape, 0, nshapes) { 392 SSOL(instance_get_shaded_shape(inst, ishape, &inst_sshape)); 393 dump_instantiated_shaded_shape_vertices(solstice, &inst_sshape); 394 } 395 396 /* Write the triangles of the receiver shade shapes */ 397 offset = 0; 398 fprintf(solstice->output, "POLYGONS %lu %lu\n", 399 (unsigned long)ntris, (unsigned long)ntris*4); 400 FOR_EACH(ishape, 0, nshapes) { 401 unsigned shape_nverts; 402 SSOL(instance_get_shaded_shape(inst, ishape, &inst_sshape)); 403 SSOL(shape_get_vertices_count(inst_sshape.shape, &shape_nverts)); 404 dump_shape_triangle_indices(solstice, inst_sshape.shape, offset); 405 offset += shape_nverts; 406 } 407 408 /* Write front faces MC estimations */ 409 fprintf(solstice->output, "CELL_DATA %lu\n", (unsigned long)ntris); 410 if(rcv->side & SRCVL_FRONT) { 411 if(rcv->per_primitive_output & SRCVL_PP_INCOMING) { 412 dump_per_primitive_mc_estimations(solstice, estimator, inst, 413 SSOL_FRONT, SRCVL_PP_INCOMING); 414 } 415 if(rcv->per_primitive_output & SRCVL_PP_ABSORBED) { 416 dump_per_primitive_mc_estimations(solstice, estimator, inst, 417 SSOL_FRONT, SRCVL_PP_ABSORBED); 418 } 419 } 420 if(rcv->side & SRCVL_BACK) { 421 if(rcv->per_primitive_output & SRCVL_PP_INCOMING) { 422 dump_per_primitive_mc_estimations(solstice, estimator, inst, 423 SSOL_BACK, SRCVL_PP_INCOMING); 424 } 425 if(rcv->per_primitive_output & SRCVL_PP_ABSORBED) { 426 dump_per_primitive_mc_estimations(solstice, estimator, inst, 427 SSOL_BACK, SRCVL_PP_ABSORBED); 428 } 429 } 430 } 431 } 432 433 static void 434 dump_path_positions(struct solstice* solstice, const struct ssol_path* path) 435 { 436 size_t ivert, nverts; 437 ASSERT(solstice && path); 438 439 SSOL(path_get_vertices_count(path, &nverts)); 440 FOR_EACH(ivert, 0, nverts) { 441 struct ssol_path_vertex vertex; 442 SSOL(path_get_vertex(path, ivert, &vertex)); 443 fprintf(solstice->output, "%f %f %f\n", SPLIT3(vertex.pos)); 444 } 445 } 446 447 static void 448 dump_path_segments 449 (struct solstice* solstice, 450 const struct ssol_path* path, 451 const size_t offset) 452 { 453 size_t i, nverts; 454 ASSERT(solstice && path); 455 456 SSOL(path_get_vertices_count(path, &nverts)); 457 ASSERT(nverts); 458 459 fprintf(solstice->output, "%lu", (unsigned long)nverts); 460 FOR_EACH(i, 0, nverts) { 461 fprintf(solstice->output, " %lu", (unsigned long)(i + offset)); 462 } 463 fprintf(solstice->output, "\n"); 464 } 465 466 static void 467 write_paths(struct solstice* solstice, struct ssol_estimator* estimator) 468 { 469 struct ssol_path path; 470 size_t ipath, npaths; 471 size_t nverts; 472 size_t offset; 473 ASSERT(solstice && estimator); 474 475 /* Write the header */ 476 fprintf(solstice->output, "# vtk DataFile Version 2.0\n"); 477 fprintf(solstice->output, "Radiative paths\n"); 478 fprintf(solstice->output, "ASCII\n"); 479 fprintf(solstice->output, "DATASET POLYDATA\n"); 480 481 /* Compute the overall number of path vertices & segments */ 482 nverts = 0; 483 SSOL(estimator_get_tracked_paths_count(estimator, &npaths)); 484 FOR_EACH(ipath, 0, npaths) { 485 size_t n; 486 SSOL(estimator_get_tracked_path(estimator, ipath, &path)); 487 SSOL(path_get_vertices_count(&path, &n)); 488 nverts += n; 489 } 490 491 /* Write the positions of the tracked paths */ 492 fprintf(solstice->output, "POINTS %lu float\n", (unsigned long)nverts); 493 FOR_EACH(ipath, 0, npaths) { 494 SSOL(estimator_get_tracked_path(estimator, ipath, &path)); 495 dump_path_positions(solstice, &path); 496 } 497 498 /* Write the segment of the tracked paths */ 499 offset = 0; 500 fprintf(solstice->output, "LINES %lu %lu\n", 501 (unsigned long)npaths, (unsigned long)(nverts + npaths)); 502 FOR_EACH(ipath, 0, npaths) { 503 size_t n; 504 SSOL(estimator_get_tracked_path(estimator, ipath, &path)); 505 SSOL(path_get_vertices_count(&path, &n)); 506 dump_path_segments(solstice, &path, offset); 507 offset += n; 508 } 509 510 /* Write path type */ 511 fprintf(solstice->output, "CELL_DATA %lu\n", (unsigned long)npaths); 512 fprintf(solstice->output, "SCALARS Radiative_path_type float 1\n"); 513 fprintf(solstice->output, "LOOKUP_TABLE path_type\n"); 514 FOR_EACH(ipath, 0, npaths) { 515 enum ssol_path_type type; 516 SSOL(estimator_get_tracked_path(estimator, ipath, &path)); 517 SSOL(path_get_type(&path, &type)); 518 switch(type) { 519 case SSOL_PATH_ERROR: fprintf(solstice->output, "0.0\n"); break; 520 case SSOL_PATH_SUCCESS: fprintf(solstice->output, "0.5\n"); break; 521 case SSOL_PATH_MISSING: fprintf(solstice->output, "0.75\n"); break; 522 case SSOL_PATH_SHADOW: fprintf(solstice->output, "1.0\n"); break; 523 default: FATAL("Unreachable code.\n"); break; 524 } 525 } 526 fprintf(solstice->output, "LOOKUP_TABLE path_type 5\n"); 527 fprintf(solstice->output, "1.0 0.0 0.0 1.0\n"); /* 0.0 = Red: for error paths */ 528 fprintf(solstice->output, "0.0 1.0 0.0 1.0\n"); /* 0.25 = Green: unused */ 529 fprintf(solstice->output, "0.0 0.0 1.0 1.0\n"); /* 0.5 = Blue: for success paths */ 530 fprintf(solstice->output, "0.0 1.0 1.0 1.0\n"); /* 0.75 = Turquoise: for missing paths */ 531 fprintf(solstice->output, "1.0 1.0 0.0 1.0\n"); /* 1.0 = Yellow: for occluded paths */ 532 } 533 534 /******************************************************************************* 535 * Local functions 536 ******************************************************************************/ 537 res_T 538 solstice_solve(struct solstice* solstice) 539 { 540 struct ssol_estimator* estimator = NULL; 541 struct ssp_rng* rng = NULL; 542 size_t max_failure; 543 res_T res = RES_OK; 544 ASSERT(solstice); 545 546 res = ssp_rng_create(solstice->allocator, SSP_RNG_THREEFRY, &rng); 547 if(res != RES_OK) { 548 fprintf(stderr, "Could not create the Random Number Generator .\n"); 549 goto error; 550 } 551 552 if(solstice->rng_state_input) { 553 rewind(solstice->rng_state_input); 554 res = ssp_rng_read(rng, solstice->rng_state_input); 555 if(res != RES_OK) { 556 fprintf(stderr, "Could not read the input RNG state.\n"); 557 goto error; 558 } 559 } 560 561 max_failure = solstice->dump_paths ? 562 solstice->nexperiments 563 : (size_t)((double)solstice->nexperiments * MAX_PERCENT_FAILURES); 564 565 res = ssol_solve(solstice->scene, rng, solstice->nexperiments, max_failure, 566 solstice->dump_paths ? &solstice->path_tracker : NULL, &estimator); 567 if(res != RES_OK) { 568 fprintf(stderr, "Error in integrating the solar flux.\n"); 569 if(!estimator) goto error; 570 } 571 572 if(solstice->dump_paths) { 573 write_paths(solstice, estimator); 574 } else { 575 write_mc_global(solstice, estimator); 576 write_per_receiver_mc_primitive(solstice, estimator); 577 } 578 579 if(solstice->rng_state_output) { 580 const struct ssp_rng* rng_state = NULL; 581 SSOL(estimator_get_rng_state(estimator, &rng_state)); 582 res = ssp_rng_write(rng_state, solstice->rng_state_output); 583 if(res != RES_OK) { 584 fprintf(stderr, "Could not write the RNG state.\n"); 585 goto error; 586 } 587 } 588 589 exit: 590 if(estimator) SSOL(estimator_ref_put(estimator)); 591 if(rng) SSP(rng_ref_put(rng)); 592 return res; 593 error: 594 goto exit; 595 } 596