test_ssol_solver1.c (24697B)
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 "ssol.h" 18 #include "test_ssol_utils.h" 19 20 #define REFLECTIVITY 0.87 21 #include "test_ssol_materials.h" 22 23 #define HALF_X 1 24 #define HALF_Y 1 25 #define PLANE_NAME SQUARE 26 #include "test_ssol_rect_geometry.h" 27 28 #include <rsys/double33.h> 29 30 #include <star/s3d.h> 31 #include <star/ssp.h> 32 33 struct spectrum_desc { 34 const double* wavelengths; 35 const double* intensities; 36 size_t count; 37 }; 38 39 static void 40 get_wlen(const size_t i, double* wlen, double* data, void* ctx) 41 { 42 const struct spectrum_desc* desc = ctx; 43 CHK(i < desc->count); 44 *wlen = desc->wavelengths[i]; 45 *data = desc->intensities[i]; 46 } 47 48 int 49 main(int argc, char** argv) 50 { 51 struct spectrum_desc desc = {0}; 52 struct mem_allocator allocator; 53 FILE* stream; 54 struct ssol_device* dev; 55 struct ssp_rng* rng; 56 struct ssp_rng* rng2; 57 const struct ssp_rng* rng_state; 58 enum ssp_rng_type rng_type0; 59 enum ssp_rng_type rng_type1; 60 struct ssol_scene* scene; 61 struct ssol_shape* dummy; 62 struct ssol_shape* square; 63 struct ssol_vertex_data attribs[1] = { SSOL_VERTEX_DATA_NULL__ }; 64 struct ssol_material *m_mtl, *m_mtl2; 65 struct ssol_material* v_mtl; 66 struct ssol_mirror_shader shader = SSOL_MIRROR_SHADER_NULL; 67 struct ssol_object *m_object, *m_object2; 68 struct ssol_object* t_object; 69 struct ssol_instance *heliostat, *heliostat2; 70 struct ssol_instance* secondary; 71 struct ssol_instance* target; 72 struct ssol_sun* sun; 73 struct ssol_sun* sun_mono; 74 struct ssol_spectrum* spectrum; 75 struct ssol_spectrum* abs_spectrum; 76 struct ssol_data extinction; 77 struct ssol_atmosphere* atm; 78 struct ssol_estimator* estimator; 79 struct ssol_mc_sampled sampled; 80 struct ssol_mc_global mc_global; 81 struct ssol_mc_receiver mc_rcv; 82 struct ssol_mc_shape mc_shape; 83 struct ssol_mc_primitive mc_prim; 84 struct ssol_path path; 85 struct ssol_path_vertex vertex; 86 double dir[3]; 87 double wavelengths[3] = { 1, 2, 3 }; 88 double intensities[3] = { 1, 0.8, 1 }; 89 double ka[3] = { 0, 0, 0 }; 90 double mono = 1.21; 91 double transform1[12]; /* 3x4 column major matrix */ 92 double transform2[12]; /* 3x4 column major matrix */ 93 double dbl; 94 size_t i, count, fcount, scount; 95 double m, std; 96 double a_m, a_std; 97 unsigned ntris; 98 (void) argc, (void) argv; 99 100 d33_splat(transform1, 0); 101 d3_splat(transform1 + 9, 0); 102 d33_rotation_pitch(transform1, PI); /* flip faces: invert normal */ 103 transform1[9] = 2; /* +2 offset along X axis */ 104 transform1[11] = 2; /* +2 offset along Z axis */ 105 106 d33_set_identity(transform2); 107 d3_splat(transform2 + 9, 0); 108 transform2[9] = 4; /* +4 offset along X axis */ 109 110 mem_init_proxy_allocator(&allocator, &mem_default_allocator); 111 112 CHK(ssol_device_create 113 (NULL, &allocator, SSOL_NTHREADS_DEFAULT, 0, &dev) == RES_OK); 114 115 CHK(ssp_rng_create(&allocator, SSP_RNG_THREEFRY, &rng) == RES_OK); 116 117 desc.wavelengths = wavelengths; 118 desc.intensities = intensities; 119 desc.count = 3; 120 CHK(ssol_spectrum_create(dev, &spectrum) == RES_OK); 121 CHK(ssol_spectrum_setup(spectrum, get_wlen, 3, &desc) == RES_OK); 122 CHK(ssol_sun_create_directional(dev, &sun) == RES_OK); 123 CHK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)) == RES_OK); 124 CHK(ssol_sun_set_spectrum(sun, spectrum) == RES_OK); 125 #define DNI 1000 126 CHK(ssol_sun_set_dni(sun, DNI) == RES_OK); 127 CHK(ssol_scene_create(dev, &scene) == RES_OK); 128 129 CHK(ssol_solve(NULL, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 130 CHK(ssol_solve(scene, NULL, 10, 0, NULL, &estimator) == RES_BAD_ARG); 131 CHK(ssol_solve(scene, rng, 0, 0, NULL, &estimator) == RES_BAD_ARG); 132 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 133 CHK(ssol_solve(scene, rng, 10, 0, NULL, NULL) == RES_BAD_ARG); 134 135 /* No geometry */ 136 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 137 138 /* Create scene content */ 139 CHK(ssol_shape_create_mesh(dev, &dummy) == RES_OK); 140 CHK(ssol_shape_create_mesh(dev, &square) == RES_OK); 141 attribs[0].usage = SSOL_POSITION; 142 attribs[0].get = get_position; 143 CHK(ssol_mesh_setup(square, SQUARE_NTRIS__, get_ids, 144 SQUARE_NVERTS__, attribs, 1, (void*)&SQUARE_DESC__) == RES_OK); 145 146 CHK(ssol_material_create_mirror(dev, &m_mtl) == RES_OK); 147 shader.normal = get_shader_normal; 148 shader.reflectivity = get_shader_reflectivity; 149 shader.roughness = get_shader_roughness; 150 CHK(ssol_mirror_setup(m_mtl, &shader, SSOL_MICROFACET_BECKMANN) == RES_OK); 151 CHK(ssol_material_create_virtual(dev, &v_mtl) == RES_OK); 152 153 CHK(ssol_object_create(dev, &m_object) == RES_OK); 154 CHK(ssol_object_add_shaded_shape(m_object, square, m_mtl, m_mtl) == RES_OK); 155 CHK(ssol_object_instantiate(m_object, &heliostat) == RES_OK); 156 CHK(ssol_instance_set_receiver(heliostat, SSOL_FRONT, 0) == RES_OK); 157 CHK(ssol_scene_attach_instance(scene, heliostat) == RES_OK); 158 159 CHK(ssol_object_instantiate(m_object, &secondary) == RES_OK); 160 CHK(ssol_instance_set_receiver(secondary, SSOL_FRONT, 0) == RES_OK); 161 CHK(ssol_instance_set_transform(secondary, transform1) == RES_OK); 162 CHK(ssol_scene_attach_instance(scene, secondary) == RES_OK); 163 164 CHK(ssol_object_create(dev, &t_object) == RES_OK); 165 CHK(ssol_object_add_shaded_shape(t_object, square, v_mtl, v_mtl) == RES_OK); 166 CHK(ssol_object_instantiate(t_object, &target) == RES_OK); 167 CHK(ssol_instance_set_transform(target, transform2) == RES_OK); 168 CHK(ssol_instance_set_receiver(target, SSOL_FRONT, 0) == RES_OK); 169 CHK(ssol_scene_attach_instance(scene, target) == RES_OK); 170 171 /* No sun */ 172 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 173 174 CHK(ssol_scene_attach_sun(scene, sun) == RES_OK); 175 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_OK); 176 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 177 178 CHK(ssol_solve 179 (scene, rng, 1, 0, &SSOL_PATH_TRACKER_DEFAULT, &estimator) == RES_OK); 180 181 CHK(ssol_estimator_get_tracked_paths_count(NULL, NULL) == RES_BAD_ARG); 182 CHK(ssol_estimator_get_tracked_paths_count(estimator, NULL) == RES_BAD_ARG); 183 CHK(ssol_estimator_get_tracked_paths_count(NULL, &count) == RES_BAD_ARG); 184 CHK(ssol_estimator_get_tracked_paths_count(estimator, &count) == RES_OK); 185 CHK(count == 1); 186 187 CHK(ssol_estimator_get_tracked_path(NULL, count, NULL) == RES_BAD_ARG); 188 CHK(ssol_estimator_get_tracked_path(estimator, count, NULL) == RES_BAD_ARG); 189 CHK(ssol_estimator_get_tracked_path(NULL, 0, NULL) == RES_BAD_ARG); 190 CHK(ssol_estimator_get_tracked_path(estimator, 0, NULL) == RES_BAD_ARG); 191 CHK(ssol_estimator_get_tracked_path(NULL, count, &path) == RES_BAD_ARG); 192 CHK(ssol_estimator_get_tracked_path(estimator, count, &path) == RES_BAD_ARG); 193 CHK(ssol_estimator_get_tracked_path(NULL, 0, &path) == RES_BAD_ARG); 194 CHK(ssol_estimator_get_tracked_path(estimator, 0, &path) == RES_OK); 195 196 CHK(ssol_path_get_vertices_count(NULL, NULL) == RES_BAD_ARG); 197 CHK(ssol_path_get_vertices_count(&path, NULL) == RES_BAD_ARG); 198 CHK(ssol_path_get_vertices_count(NULL, &count) == RES_BAD_ARG); 199 CHK(ssol_path_get_vertices_count(&path, &count) == RES_OK); 200 CHK(count != 0); 201 202 CHK(ssol_path_get_vertex(NULL, count, NULL) == RES_BAD_ARG); 203 CHK(ssol_path_get_vertex(&path, count, NULL) == RES_BAD_ARG); 204 CHK(ssol_path_get_vertex(NULL, 0, NULL) == RES_BAD_ARG); 205 CHK(ssol_path_get_vertex(&path, 0, NULL) == RES_BAD_ARG); 206 CHK(ssol_path_get_vertex(NULL, count, &vertex) == RES_BAD_ARG); 207 CHK(ssol_path_get_vertex(&path, count, &vertex) == RES_BAD_ARG); 208 CHK(ssol_path_get_vertex(NULL, 0, &vertex) == RES_BAD_ARG); 209 FOR_EACH(i, 0, count) { 210 CHK(ssol_path_get_vertex(&path, i, &vertex) == RES_OK); 211 } 212 213 CHK(ssol_estimator_get_sampled_area(NULL, NULL) == RES_BAD_ARG); 214 CHK(ssol_estimator_get_sampled_area(estimator, NULL) == RES_BAD_ARG); 215 CHK(ssol_estimator_get_sampled_area(NULL, &dbl) == RES_BAD_ARG); 216 CHK(ssol_estimator_get_sampled_area(estimator, &dbl) == RES_OK); 217 CHK(eq_eps(dbl, 12, DBL_EPSILON) == 1); 218 219 CHK(ssol_estimator_get_realisation_count(NULL, NULL) == RES_BAD_ARG); 220 CHK(ssol_estimator_get_realisation_count(estimator, NULL) == RES_BAD_ARG); 221 CHK(ssol_estimator_get_realisation_count(NULL, &count) == RES_BAD_ARG); 222 CHK(ssol_estimator_get_realisation_count(estimator, &count) == RES_OK); 223 CHK(count == 1); 224 225 CHK(ssol_estimator_get_failed_count(NULL, NULL) == RES_BAD_ARG); 226 CHK(ssol_estimator_get_failed_count(estimator, NULL) == RES_BAD_ARG); 227 CHK(ssol_estimator_get_failed_count(NULL, &count) == RES_BAD_ARG); 228 CHK(ssol_estimator_get_failed_count(estimator, &fcount) == RES_OK); 229 CHK(fcount == 0); 230 231 CHK(ssol_estimator_get_sampled_count(NULL, NULL) == RES_BAD_ARG); 232 CHK(ssol_estimator_get_sampled_count(estimator, NULL) == RES_BAD_ARG); 233 CHK(ssol_estimator_get_sampled_count(NULL, &scount) == RES_BAD_ARG); 234 CHK(ssol_estimator_get_sampled_count(estimator, &scount) == RES_OK); 235 CHK(scount == 3); 236 237 CHK(ssol_estimator_get_mc_sampled(NULL, heliostat, &sampled) == RES_BAD_ARG); 238 CHK(ssol_estimator_get_mc_sampled(estimator, NULL, &sampled) == RES_BAD_ARG); 239 CHK(ssol_estimator_get_mc_sampled(estimator, heliostat, NULL) == RES_BAD_ARG); 240 CHK(ssol_estimator_get_mc_sampled(estimator, heliostat, &sampled) == RES_OK); 241 242 CHK(ssol_estimator_get_mc_global(NULL, &mc_global) == RES_BAD_ARG); 243 CHK(ssol_estimator_get_mc_global(estimator, NULL) == RES_BAD_ARG); 244 CHK(ssol_estimator_get_mc_global(estimator, &mc_global) == RES_OK); 245 246 CHK(ssol_estimator_get_rng_state(NULL, &rng_state) == RES_BAD_ARG); 247 CHK(ssol_estimator_get_rng_state(estimator, NULL) == RES_BAD_ARG); 248 CHK(ssol_estimator_get_rng_state(estimator, &rng_state) == RES_OK); 249 CHK(ssp_rng_get_type(rng_state, &rng_type0) == RES_OK); 250 CHK(ssp_rng_get_type(rng, &rng_type1) == RES_OK); 251 CHK(rng_type0 == rng_type1); 252 253 /* Clone the rng_state */ 254 CHK(stream = tmpfile()); 255 CHK(ssp_rng_create(&allocator, SSP_RNG_THREEFRY, &rng2) == RES_OK); 256 CHK(ssp_rng_write(rng_state, stream) == RES_OK); 257 rewind(stream); 258 CHK(ssp_rng_read(rng2, stream) == RES_OK); 259 CHK(fclose(stream) == 0); 260 CHK(ssp_rng_get(rng2) != ssp_rng_get(rng)); 261 CHK(ssp_rng_ref_put(rng2) == RES_OK); 262 263 CHK(ssol_estimator_ref_get(NULL) == RES_BAD_ARG); 264 CHK(ssol_estimator_ref_get(estimator) == RES_OK); 265 CHK(ssol_estimator_ref_put(NULL) == RES_BAD_ARG); 266 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 267 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 268 269 /* No geometry to sample */ 270 CHK(ssol_instance_sample(target, 0) == RES_OK); 271 CHK(ssol_instance_sample(secondary, 0) == RES_OK); 272 CHK(ssol_instance_sample(heliostat, 0) == RES_OK); 273 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 274 CHK(ssol_estimator_get_mc_sampled(estimator, heliostat, &sampled) == RES_BAD_ARG); 275 276 CHK(ssol_instance_sample(target, 1) == RES_OK); 277 CHK(ssol_instance_sample(secondary, 1) == RES_OK); 278 CHK(ssol_instance_sample(heliostat, 1) == RES_OK); 279 280 /* No attached sun */ 281 CHK(ssol_scene_detach_sun(scene, sun) == RES_OK); 282 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 283 CHK(ssol_sun_ref_put(sun) == RES_OK); 284 285 /* Sun with no spectrum */ 286 CHK(ssol_sun_create_directional(dev, &sun) == RES_OK); 287 CHK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)) == RES_OK); 288 CHK(ssol_sun_set_dni(sun, DNI) == RES_OK); 289 CHK(ssol_scene_attach_sun(scene, sun) == RES_OK); 290 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 291 CHK(ssol_scene_detach_sun(scene, sun) == RES_OK); 292 CHK(ssol_sun_ref_put(sun) == RES_OK); 293 294 /* Sun with undefined DNI */ 295 CHK(ssol_sun_create_directional(dev, &sun) == RES_OK); 296 CHK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)) == RES_OK); 297 CHK(ssol_sun_set_spectrum(sun, spectrum) == RES_OK); 298 CHK(ssol_scene_attach_sun(scene, sun) == RES_OK); 299 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_BAD_ARG); 300 CHK(ssol_sun_set_dni(sun, DNI) == RES_OK); 301 302 /* No receiver in scene */ 303 CHK(ssol_instance_set_receiver(heliostat, 0, 0) == RES_OK); 304 CHK(ssol_instance_set_receiver(secondary, 0, 0) == RES_OK); 305 CHK(ssol_instance_set_receiver(target, 0, 0) == RES_OK); 306 CHK(ssol_solve(scene, rng, 10, 0, NULL, &estimator) == RES_OK); 307 CHK(ssol_instance_set_receiver(heliostat, SSOL_FRONT, 0) == RES_OK); 308 CHK(ssol_instance_set_receiver(secondary, SSOL_FRONT, 0) == RES_OK); 309 CHK(ssol_instance_set_receiver(target, SSOL_FRONT, 0) == RES_OK); 310 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 311 312 /* Can sample any geometry; variance is high */ 313 #define N__ 10000 314 #define GET_MC_RCV ssol_estimator_get_mc_receiver 315 #define GET_MC_GLOBAL ssol_estimator_get_mc_global 316 CHK(ssol_solve(scene, rng, N__, 0, NULL, &estimator) == RES_OK); 317 CHK(ssol_estimator_get_realisation_count(estimator, &count) == RES_OK); 318 CHK(count == N__); 319 CHK(ssol_estimator_get_failed_count(estimator, &fcount) == RES_OK); 320 CHK(fcount == 0); 321 #define COS cos(PI / 4) 322 #define DNI_cos (DNI * COS) 323 m = 4 * DNI_cos; 324 #define SQR(x) ((x)*(x)) 325 dbl = sqrt((SQR(12 * DNI_cos) / 3 - SQR(4 * DNI_cos)) / (double)count); 326 std = dbl; 327 /* Target was sampled but shadowed by secondary */ 328 CHK(ssol_estimator_get_mc_global(estimator, &mc_global) == RES_OK); 329 print_global(&mc_global); 330 CHK(eq_eps(mc_global.shadowed.E, m, 2 * dbl) == 1); 331 CHK(eq_eps(mc_global.missing.E, 2*m, 2*mc_global.missing.SE) == 1); 332 CHK(GET_MC_RCV(NULL, NULL, SSOL_BACK, NULL) == RES_BAD_ARG); 333 CHK(GET_MC_RCV(estimator, NULL, SSOL_BACK, NULL) == RES_BAD_ARG); 334 CHK(GET_MC_RCV(NULL, target, SSOL_BACK, NULL) == RES_BAD_ARG); 335 CHK(GET_MC_RCV(estimator, target, SSOL_BACK, NULL) == RES_BAD_ARG); 336 CHK(GET_MC_RCV(NULL, NULL, SSOL_BACK, &mc_rcv) == RES_BAD_ARG); 337 CHK(GET_MC_RCV(estimator, NULL, SSOL_BACK, &mc_rcv) == RES_BAD_ARG); 338 CHK(GET_MC_RCV(NULL, target, SSOL_BACK, &mc_rcv) == RES_BAD_ARG); 339 CHK(GET_MC_RCV(estimator, target, SSOL_BACK, &mc_rcv) == RES_BAD_ARG); 340 CHK(GET_MC_RCV(NULL, NULL, SSOL_FRONT, NULL) == RES_BAD_ARG); 341 CHK(GET_MC_RCV(estimator, NULL, SSOL_FRONT, NULL) == RES_BAD_ARG); 342 CHK(GET_MC_RCV(NULL, target, SSOL_FRONT, NULL) == RES_BAD_ARG); 343 CHK(GET_MC_RCV(estimator, target, SSOL_FRONT, NULL) == RES_BAD_ARG); 344 CHK(GET_MC_RCV(NULL, NULL, SSOL_FRONT, &mc_rcv) == RES_BAD_ARG); 345 CHK(GET_MC_RCV(estimator, NULL, SSOL_FRONT, &mc_rcv) == RES_BAD_ARG); 346 CHK(GET_MC_RCV(NULL, target, SSOL_FRONT, &mc_rcv) == RES_BAD_ARG); 347 CHK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv) == RES_OK); 348 printf("Ir(target) = %g +/- %g\n", 349 mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); 350 CHK(eq_eps(mc_rcv.incoming_flux.E, m, 3 * std) == 1); 351 CHK(eq_eps(mc_rcv.incoming_flux.SE, std, std*1e-2) == 1); 352 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 353 354 /* Sample primary mirror only; variance is low */ 355 CHK(ssol_instance_sample(target, 0) == RES_OK); 356 CHK(ssol_instance_sample(secondary, 0) == RES_OK); 357 358 CHK(ssol_solve(scene, rng, N__, 0, NULL, &estimator) == RES_OK); 359 CHK(ssol_estimator_get_realisation_count(estimator, &count) == RES_OK); 360 CHK(count == N__); 361 m = 4 * DNI_cos; 362 std = 0; 363 CHK(ssol_estimator_get_mc_global(estimator, &mc_global) == RES_OK); 364 print_global(&mc_global); 365 CHK(eq_eps(mc_global.shadowed.E, 0, 1e-4) == 1); 366 CHK(eq_eps(mc_global.missing.E, m, 1e-4) == 1); 367 CHK(eq_eps(mc_global.cos_factor.E, COS, 1e-4) == 1); 368 CHK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv) == RES_OK); 369 printf("Ir(target) = %g +/- %g\n", 370 mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); 371 CHK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8) == 1); 372 CHK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4) == 1); 373 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 374 375 /* Check atmosphere model; with no extinction result is unchanged */ 376 CHK(ssol_atmosphere_create(dev, &atm) == RES_OK); 377 extinction.type = SSOL_DATA_REAL; 378 extinction.value.real = 0; 379 CHK(ssol_atmosphere_set_extinction(atm, &extinction) == RES_OK); 380 CHK(ssol_scene_attach_atmosphere(scene, atm) == RES_OK); 381 382 CHK(ssol_solve(scene, rng, N__, 0, NULL, &estimator) == RES_OK); 383 CHK(ssol_estimator_get_realisation_count(estimator, &count) == RES_OK); 384 CHK(count == N__); 385 m = 4 * DNI_cos; 386 std = 0; 387 CHK(ssol_scene_detach_atmosphere(scene, atm) == RES_OK); 388 CHK(ssol_atmosphere_ref_put(atm) == RES_OK); 389 CHK(ssol_estimator_get_mc_global(estimator, &mc_global) == RES_OK); 390 print_global(&mc_global); 391 CHK(eq_eps(mc_global.shadowed.E, 0, 1e-4) == 1); 392 CHK(eq_eps(mc_global.missing.E, m, 1e-4) == 1); 393 CHK(eq_eps(mc_global.cos_factor.E, COS, 1e-4) == 1); 394 CHK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv) == RES_OK); 395 printf("Ir(target) = %g +/- %g\n", 396 mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); 397 CHK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8) == 1); 398 CHK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4) == 1); 399 CHK(eq_eps(mc_global.cos_factor.E, COS, 1e-4) == 1); 400 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 401 402 /* Check atmosphere model and imperfect mirror: there are losses */ 403 CHK(ssol_scene_detach_instance(scene, heliostat) == RES_OK); 404 405 CHK(ssol_material_create_mirror(dev, &m_mtl2) == RES_OK); 406 shader.normal = get_shader_normal; 407 shader.reflectivity = get_shader_reflectivity_2; 408 shader.roughness = get_shader_roughness; 409 CHK(ssol_mirror_setup(m_mtl2, &shader, SSOL_MICROFACET_BECKMANN) == RES_OK); 410 411 CHK(ssol_object_create(dev, &m_object2) == RES_OK); 412 CHK(ssol_object_add_shaded_shape(m_object2, square, m_mtl2, m_mtl2) == RES_OK); 413 CHK(ssol_object_instantiate(m_object2, &heliostat2) == RES_OK); 414 CHK(ssol_instance_set_receiver(heliostat2, SSOL_FRONT, 0) == RES_OK); 415 CHK(ssol_scene_attach_instance(scene, heliostat2) == RES_OK); 416 417 #define KA 0.03 418 extinction.value.real = KA; 419 CHK(ssol_spectrum_create(dev, &abs_spectrum) == RES_OK); 420 CHK(ssol_spectrum_setup(abs_spectrum, get_wlen, 3, &desc) == RES_OK); 421 CHK(ssol_atmosphere_create(dev, &atm) == RES_OK); 422 CHK(ssol_atmosphere_set_extinction(atm, &extinction) == RES_OK); 423 CHK(ssol_scene_attach_atmosphere(scene, atm) == RES_OK); 424 CHK(ssol_instance_set_receiver(target, SSOL_FRONT, 1) == RES_OK); 425 426 CHK(ssol_solve(scene, rng, N__, 0, NULL, &estimator) == RES_OK); 427 CHK(ssol_estimator_get_realisation_count(estimator, &count) == RES_OK); 428 CHK(count == N__); 429 #define K (exp(-KA * 4 * sqrt(2))) 430 a_m = REFLECTIVITY * 4 * K * DNI_cos; 431 a_std = 0; 432 CHK(ssol_estimator_get_mc_global(estimator, &mc_global) == RES_OK); 433 print_global(&mc_global); 434 CHK(eq_eps(mc_global.shadowed.E, 0, 1e-4) == 1); 435 CHK(eq_eps( 436 mc_global.missing.E + mc_global.shadowed.E + mc_global.absorbed_by_receivers.E 437 + mc_global.extinguished_by_atmosphere.E + mc_global.other_absorbed.E, 438 m, 1e-4)); 439 CHK(eq_eps(mc_global.cos_factor.E, COS, 1e-4) == 1); 440 CHK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv) == RES_OK); 441 print_rcv(&mc_rcv); 442 CHK(ssol_estimator_get_sampled_count(estimator, &scount) == RES_OK); 443 CHK(ssol_estimator_get_mc_sampled(estimator, heliostat, &sampled) == RES_BAD_ARG); 444 CHK(ssol_estimator_get_mc_sampled(estimator, heliostat2, &sampled) == RES_OK); 445 446 CHK(eq_eps(mc_rcv.incoming_flux.E, a_m, 1e-4) == 1); 447 CHK(eq_eps(mc_rcv.incoming_flux.SE, a_std, 1e-4) == 1); 448 449 CHK(ssol_mc_receiver_get_mc_shape(NULL, NULL, NULL) == RES_BAD_ARG); 450 CHK(ssol_mc_receiver_get_mc_shape(&mc_rcv, NULL, NULL) == RES_BAD_ARG); 451 CHK(ssol_mc_receiver_get_mc_shape(NULL, square, NULL) == RES_BAD_ARG); 452 CHK(ssol_mc_receiver_get_mc_shape(&mc_rcv, square, NULL) == RES_BAD_ARG); 453 CHK(ssol_mc_receiver_get_mc_shape(NULL, NULL, &mc_shape) == RES_BAD_ARG); 454 CHK(ssol_mc_receiver_get_mc_shape(&mc_rcv, NULL, &mc_shape) == RES_BAD_ARG); 455 CHK(ssol_mc_receiver_get_mc_shape(NULL, square, &mc_shape) == RES_BAD_ARG); 456 CHK(ssol_mc_receiver_get_mc_shape(&mc_rcv, dummy, &mc_shape) == RES_BAD_ARG); 457 CHK(ssol_mc_receiver_get_mc_shape(&mc_rcv, square, &mc_shape) == RES_OK); 458 459 CHK(ssol_shape_get_triangles_count(square, &ntris) == RES_OK); 460 CHK(ntris != 0); 461 462 CHK(ssol_mc_shape_get_mc_primitive(NULL, ntris, NULL) == RES_BAD_ARG); 463 CHK(ssol_mc_shape_get_mc_primitive(&mc_shape, ntris, NULL) == RES_BAD_ARG); 464 CHK(ssol_mc_shape_get_mc_primitive(NULL, 0, NULL) == RES_BAD_ARG); 465 CHK(ssol_mc_shape_get_mc_primitive(&mc_shape, 0, NULL) == RES_BAD_ARG); 466 CHK(ssol_mc_shape_get_mc_primitive(NULL, ntris, &mc_prim) == RES_BAD_ARG); 467 CHK(ssol_mc_shape_get_mc_primitive(&mc_shape, ntris, &mc_prim) == RES_BAD_ARG); 468 CHK(ssol_mc_shape_get_mc_primitive(NULL, 0, &mc_prim) == RES_BAD_ARG); 469 470 dbl = 0; 471 FOR_EACH(i, 0, ntris) { 472 double v0[3], v1[3], v2[3], E1[3], E2[3], N[3], area; 473 unsigned ids[3]; 474 CHK(ssol_shape_get_triangle_indices(square, (unsigned)i, ids) == RES_OK); 475 CHK(ssol_shape_get_vertex_attrib(square, ids[0], SSOL_POSITION, v0) == RES_OK); 476 CHK(ssol_shape_get_vertex_attrib(square, ids[1], SSOL_POSITION, v1) == RES_OK); 477 CHK(ssol_shape_get_vertex_attrib(square, ids[2], SSOL_POSITION, v2) == RES_OK); 478 area = d3_len(d3_cross(N, d3_sub(E1, v1, v0), d3_sub(E2, v2, v0))) * 0.5; 479 480 CHK(ssol_mc_shape_get_mc_primitive(&mc_shape, (unsigned)i, &mc_prim) == RES_OK); 481 dbl += mc_prim.incoming_flux.E * area; 482 } 483 484 CHK(eq_eps(dbl, a_m, 1e-4) == 1); 485 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 486 487 CHK(ssol_scene_detach_instance(scene, heliostat2) == RES_OK); 488 CHK(ssol_scene_attach_instance(scene, heliostat) == RES_OK); 489 CHK(ssol_instance_set_receiver(target, SSOL_FRONT, 0) == RES_OK); 490 491 /* Check a monochromatic sun */ 492 desc.wavelengths = &mono; 493 desc.intensities = intensities; 494 desc.count = 1; 495 CHK(ssol_spectrum_setup(spectrum, get_wlen, 1, &desc) == RES_OK); 496 CHK(ssol_sun_create_directional(dev, &sun_mono) == RES_OK); 497 CHK(ssol_sun_set_direction(sun_mono, d3(dir, 1, 0, -1)) == RES_OK); 498 CHK(ssol_sun_set_spectrum(sun_mono, spectrum) == RES_OK); 499 CHK(ssol_sun_set_dni(sun_mono, DNI) == RES_OK); 500 CHK(ssol_scene_detach_sun(scene, sun) == RES_OK); 501 CHK(ssol_scene_attach_sun(scene, sun_mono) == RES_OK); 502 ka[1] = 0.2; ka[0] = ka[2] = 0.1; 503 desc.wavelengths = wavelengths; 504 desc.intensities = ka; 505 desc.count = 3; 506 CHK(ssol_spectrum_setup(abs_spectrum, get_wlen, 3, &desc) == RES_OK); 507 CHK(ssol_spectrum_setup(abs_spectrum, get_wlen, 2, &desc) == RES_OK); 508 extinction.type = SSOL_DATA_SPECTRUM; 509 extinction.value.spectrum = abs_spectrum; 510 CHK(ssol_atmosphere_set_extinction(atm, &extinction) == RES_OK); 511 512 CHK(ssol_solve(scene, rng, N__, 0, NULL, &estimator) == RES_OK); 513 CHK(ssol_estimator_get_realisation_count(estimator, &count) == RES_OK); 514 CHK(count == N__); 515 #define K2 (exp(-0.121 * 4 * sqrt(2))) 516 m = 4 * K2 * DNI_cos; 517 std = 0; 518 CHK(ssol_estimator_get_mc_global(estimator, &mc_global) == RES_OK); 519 print_global(&mc_global); 520 CHK(eq_eps(mc_global.shadowed.E, 0, 1e-4) == 1); 521 CHK(eq_eps(mc_global.missing.E, m, 1e-4) == 1); 522 CHK(eq_eps(mc_global.cos_factor.E, COS, 1e-4) == 1); 523 CHK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv) == RES_OK); 524 printf("Ir(target) = %g +/- %g\n", 525 mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); 526 print_rcv(&mc_rcv); 527 CHK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-4) == 1); 528 CHK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4) == 1); 529 530 /* Free data */ 531 CHK(ssol_instance_ref_put(heliostat2) == RES_OK); 532 CHK(ssol_object_ref_put(m_object2) == RES_OK); 533 CHK(ssol_material_ref_put(m_mtl2) == RES_OK); 534 CHK(ssol_instance_ref_put(heliostat) == RES_OK); 535 CHK(ssol_instance_ref_put(secondary) == RES_OK); 536 CHK(ssol_instance_ref_put(target) == RES_OK); 537 CHK(ssol_object_ref_put(m_object) == RES_OK); 538 CHK(ssol_object_ref_put(t_object) == RES_OK); 539 CHK(ssol_shape_ref_put(dummy) == RES_OK); 540 CHK(ssol_shape_ref_put(square) == RES_OK); 541 CHK(ssol_material_ref_put(m_mtl) == RES_OK); 542 CHK(ssol_material_ref_put(v_mtl) == RES_OK); 543 CHK(ssol_device_ref_put(dev) == RES_OK); 544 CHK(ssol_scene_ref_put(scene) == RES_OK); 545 CHK(ssol_atmosphere_ref_put(atm) == RES_OK); 546 CHK(ssol_estimator_ref_put(estimator) == RES_OK); 547 CHK(ssol_spectrum_ref_put(spectrum) == RES_OK); 548 CHK(ssol_spectrum_ref_put(abs_spectrum) == RES_OK); 549 CHK(ssol_sun_ref_put(sun) == RES_OK); 550 CHK(ssol_sun_ref_put(sun_mono) == RES_OK); 551 CHK(ssp_rng_ref_put(rng) == RES_OK); 552 553 check_memory_allocator(&allocator); 554 mem_shutdown_proxy_allocator(&allocator); 555 CHK(mem_allocated_size() == 0); 556 557 return 0; 558 }