solstice-solver

Solver library of the solstice app
git clone git://git.meso-star.com/solstice-solver.git
Log | Files | Refs | README | LICENSE

commit 391b88739003e4be5c069784988e94fbafa5edd1
parent 39b4e392b048f2717e3414bbe72e2570ca1dc113
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Thu,  7 Sep 2017 14:20:20 +0200

Fix MC weights and add new ones for receivers.

Lot of renames, though.

Diffstat:
Msrc/ssol.h | 52+++++++++++++++++++++++++++++++++++++++++-----------
Msrc/ssol_estimator.c | 20+++++++++++++-------
Msrc/ssol_estimator_c.h | 57+++++++++++++++++++++++++++++++++++++--------------------
Msrc/ssol_material.c | 19+++++++++----------
Msrc/ssol_material_c.h | 3+++
Msrc/ssol_mc_receiver.c | 28++++++++++++++++------------
Msrc/ssol_scene.c | 4++++
Msrc/ssol_scene_c.h | 1+
Msrc/ssol_solver.c | 353+++++++++++++++++++++++++++++++++++++++++++++++++++----------------------------
Msrc/test_ssol_by_receiver_integration.c | 16++++++++--------
Msrc/test_ssol_solver1.c | 96+++++++++++++++++++++++++++----------------------------------------------------
Msrc/test_ssol_solver2.c | 8++++----
Msrc/test_ssol_solver2b.c | 6+++---
Msrc/test_ssol_solver3.c | 6+++---
Msrc/test_ssol_solver4.c | 12++++++------
Msrc/test_ssol_solver5.c | 6+++---
Msrc/test_ssol_solver6.c | 17+++++++----------
Msrc/test_ssol_solver7.c | 18+++++++-----------
Msrc/test_ssol_solver8.c | 6+++---
Msrc/test_ssol_solver9.c | 6+++---
Msrc/test_ssol_utils.h | 43+++++++++++++++++++++++++++++++++++++++++++
21 files changed, 474 insertions(+), 303 deletions(-)

diff --git a/src/ssol.h b/src/ssol.h @@ -388,11 +388,11 @@ static const struct ssol_mc_result SSOL_MC_RESULT_NULL = SSOL_MC_RESULT_NULL__; struct ssol_mc_global { struct ssol_mc_result cos_factor; /* [0 1] */ - struct ssol_mc_result absorbed; /* In W */ + struct ssol_mc_result absorbed_by_receivers; /* In W */ struct ssol_mc_result shadowed; /* In W */ struct ssol_mc_result missing; /* In W */ - struct ssol_mc_result atmosphere; /* In W */ - struct ssol_mc_result reflectivity; /* In W */ + struct ssol_mc_result absorbed_by_atmosphere; /* In W */ + struct ssol_mc_result other_absorbed; /* In W */ }; #define SSOL_MC_GLOBAL_NULL__ { \ SSOL_MC_RESULT_NULL__, \ @@ -405,10 +405,16 @@ struct ssol_mc_global { static const struct ssol_mc_global SSOL_MC_GLOBAL_NULL = SSOL_MC_GLOBAL_NULL__; struct ssol_mc_receiver { - struct ssol_mc_result integrated_irradiance; /* In W */ - struct ssol_mc_result integrated_absorbed_irradiance; /* In W */ - struct ssol_mc_result absorptivity_loss; /* In W */ - struct ssol_mc_result reflectivity_loss; /* In W */ + struct ssol_mc_result incoming_flux; /* In W */ + struct ssol_mc_result incoming_if_no_atm_loss; /* In W */ + struct ssol_mc_result incoming_if_no_field_loss; /* In W */ + struct ssol_mc_result incoming_lost_in_field; /* In W */ + struct ssol_mc_result incoming_lost_in_atmosphere; /* In W */ + struct ssol_mc_result absorbed_flux; /* In W */ + struct ssol_mc_result absorbed_if_no_atm_loss; /* In W */ + struct ssol_mc_result absorbed_if_no_field_loss; /* In W */ + struct ssol_mc_result absorbed_lost_in_field; /* In W */ + struct ssol_mc_result absorbed_lost_in_atmosphere; /* In W */ /* Internal data */ size_t N__; @@ -420,6 +426,12 @@ struct ssol_mc_receiver { SSOL_MC_RESULT_NULL__, \ SSOL_MC_RESULT_NULL__, \ SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ 0, NULL, NULL \ } static const struct ssol_mc_receiver SSOL_MC_RECEIVER_NULL = @@ -430,6 +442,12 @@ static const struct ssol_mc_receiver SSOL_MC_RECEIVER_NULL = { -1, -1, -1 }, \ { -1, -1, -1 }, \ { -1, -1, -1 }, \ + { -1, -1, -1 }, \ + { -1, -1, -1 }, \ + { -1, -1, -1 }, \ + { -1, -1, -1 }, \ + { -1, -1, -1 }, \ + { -1, -1, -1 }, \ 0, NULL, NULL \ } @@ -449,15 +467,27 @@ struct ssol_mc_sampled { }; struct ssol_mc_primitive { - struct ssol_mc_result integrated_irradiance; /* In W */ - struct ssol_mc_result integrated_absorbed_irradiance; /* In W */ - struct ssol_mc_result absorptivity_loss; /* In W */ - struct ssol_mc_result reflectivity_loss; /* In W */ + struct ssol_mc_result incoming_flux; /* In W */ + struct ssol_mc_result incoming_if_no_atm_loss; /* In W */ + struct ssol_mc_result incoming_if_no_field_loss; /* In W */ + struct ssol_mc_result incoming_lost_in_field; /* In W */ + struct ssol_mc_result incoming_lost_in_atmosphere; /* In W */ + struct ssol_mc_result absorbed_flux; /* In W */ + struct ssol_mc_result absorbed_if_no_atm_loss; /* In W */ + struct ssol_mc_result absorbed_if_no_field_loss; /* In W */ + struct ssol_mc_result absorbed_lost_in_field; /* In W */ + struct ssol_mc_result absorbed_lost_in_atmosphere; /* In W */ }; #define SSOL_MC_PRIMITIVE_NULL__ { \ SSOL_MC_RESULT_NULL__, \ SSOL_MC_RESULT_NULL__, \ SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ SSOL_MC_RESULT_NULL__ \ } static const struct ssol_mc_primitive SSOL_MC_PRIMITIVE_NULL = diff --git a/src/ssol_estimator.c b/src/ssol_estimator.c @@ -120,11 +120,11 @@ ssol_estimator_get_mc_global global->Name.SE = sqrt(global->Name.V / N); \ } (void)0 SETUP_MC_RESULT(cos_factor); - SETUP_MC_RESULT(absorbed); + SETUP_MC_RESULT(absorbed_by_receivers); SETUP_MC_RESULT(shadowed); SETUP_MC_RESULT(missing); - SETUP_MC_RESULT(atmosphere); - SETUP_MC_RESULT(reflectivity); + SETUP_MC_RESULT(absorbed_by_atmosphere); + SETUP_MC_RESULT(other_absorbed); #undef SETUP_MC_RESULT return RES_OK; } @@ -171,10 +171,16 @@ ssol_estimator_get_mc_sampled_x_receiver rcv->Name.V = rcv->Name.V > 0 ? rcv->Name.V : 0; \ rcv->Name.SE = sqrt(rcv->Name.V / N); \ } (void)0 - SETUP_MC_RESULT(integrated_irradiance); - SETUP_MC_RESULT(integrated_absorbed_irradiance); - SETUP_MC_RESULT(absorptivity_loss); - SETUP_MC_RESULT(reflectivity_loss); + SETUP_MC_RESULT(incoming_flux); + SETUP_MC_RESULT(incoming_if_no_atm_loss); + SETUP_MC_RESULT(incoming_if_no_field_loss); + SETUP_MC_RESULT(incoming_lost_in_field); + SETUP_MC_RESULT(incoming_lost_in_atmosphere); + SETUP_MC_RESULT(absorbed_flux); + SETUP_MC_RESULT(absorbed_if_no_atm_loss); + SETUP_MC_RESULT(absorbed_if_no_field_loss); + SETUP_MC_RESULT(absorbed_lost_in_field); + SETUP_MC_RESULT(absorbed_lost_in_atmosphere); #undef SETUP_MC_RESULT rcv->mc__ = mc_rcv1; rcv->N__ = mc_samp->nb_samples; diff --git a/src/ssol_estimator_c.h b/src/ssol_estimator_c.h @@ -35,16 +35,28 @@ struct mc_data { #define MC_DATA_NULL__ { 0, 0 } static const struct mc_data MC_DATA_NULL = MC_DATA_NULL__; -#define MC_RECEIVER_DATA \ - struct mc_data integrated_irradiance; /* In W */ \ - struct mc_data integrated_absorbed_irradiance; /* In W */ \ - struct mc_data absorptivity_loss; /* In W */ \ - struct mc_data reflectivity_loss; /* In W */ +#define MC_RECEIVER_DATA \ + struct mc_data incoming_flux; /* In W */ \ + struct mc_data incoming_if_no_atm_loss; /* In W */ \ + struct mc_data incoming_if_no_field_loss; /* In W */ \ + struct mc_data incoming_lost_in_field; /* In W */ \ + struct mc_data incoming_lost_in_atmosphere; /* In W */ \ + struct mc_data absorbed_flux; /* In W */ \ + struct mc_data absorbed_if_no_atm_loss; /* In W */ \ + struct mc_data absorbed_if_no_field_loss; /* In W */ \ + struct mc_data absorbed_lost_in_field; /* In W */ \ + struct mc_data absorbed_lost_in_atmosphere; /* In W */ #define MC_RECEIVER_DATA_NULL__ \ MC_DATA_NULL__, \ MC_DATA_NULL__, \ MC_DATA_NULL__, \ + MC_DATA_NULL__, \ + MC_DATA_NULL__, \ + MC_DATA_NULL__, \ + MC_DATA_NULL__, \ + MC_DATA_NULL__, \ + MC_DATA_NULL__, \ MC_DATA_NULL__ /******************************************************************************* @@ -142,15 +154,26 @@ struct mc_receiver_1side { struct htable_shape2mc shape2mc; }; +#define COPY_LOCAL_MC_WEIGHTS(Dst, Src) { \ + (Dst)->incoming_flux = (Src)->incoming_flux; \ + (Dst)->incoming_if_no_atm_loss = (Src)->incoming_if_no_atm_loss; \ + (Dst)->incoming_if_no_field_loss = (Src)->incoming_if_no_field_loss; \ + (Dst)->incoming_lost_in_atmosphere = (Src)->incoming_lost_in_atmosphere; \ + (Dst)->incoming_lost_in_field = (Src)->incoming_lost_in_field; \ + (Dst)->absorbed_flux = (Src)->absorbed_flux; \ + (Dst)->absorbed_if_no_atm_loss = (Src)->absorbed_if_no_atm_loss; \ + (Dst)->absorbed_if_no_field_loss = (Src)->absorbed_if_no_field_loss; \ + (Dst)->absorbed_lost_in_atmosphere = (Src)->absorbed_lost_in_atmosphere; \ + (Dst)->absorbed_lost_in_field = (Src)->absorbed_lost_in_field; \ +} (void)0 + static INLINE void mc_receiver_1side_init (struct mem_allocator* allocator, struct mc_receiver_1side* mc) { + static const struct mc_receiver_1side r_null = { MC_RECEIVER_DATA_NULL__ }; ASSERT(mc); - mc->integrated_irradiance = MC_DATA_NULL; - mc->integrated_absorbed_irradiance = MC_DATA_NULL; - mc->absorptivity_loss = MC_DATA_NULL; - mc->reflectivity_loss = MC_DATA_NULL; + COPY_LOCAL_MC_WEIGHTS(mc, &r_null); htable_shape2mc_init(allocator, &mc->shape2mc); } @@ -166,10 +189,7 @@ mc_receiver_1side_copy (struct mc_receiver_1side* dst, const struct mc_receiver_1side* src) { ASSERT(dst && src); - dst->integrated_irradiance = src->integrated_irradiance; - dst->integrated_absorbed_irradiance = src->integrated_absorbed_irradiance; - dst->absorptivity_loss = src->absorptivity_loss; - dst->reflectivity_loss = src->reflectivity_loss; + COPY_LOCAL_MC_WEIGHTS(dst, src); return htable_shape2mc_copy(&dst->shape2mc, &src->shape2mc); } @@ -178,10 +198,7 @@ mc_receiver_1side_copy_and_release (struct mc_receiver_1side* dst, struct mc_receiver_1side* src) { ASSERT(dst && src); - dst->integrated_irradiance = src->integrated_irradiance; - dst->integrated_absorbed_irradiance = src->integrated_absorbed_irradiance; - dst->absorptivity_loss = src->absorptivity_loss; - dst->reflectivity_loss = src->reflectivity_loss; + COPY_LOCAL_MC_WEIGHTS(dst, src); return htable_shape2mc_copy_and_release(&dst->shape2mc, &src->shape2mc); } @@ -448,11 +465,11 @@ struct ssol_estimator { /* Implicit MC computations */ struct mc_data cos_factor; - struct mc_data absorbed; + struct mc_data absorbed_by_receivers; struct mc_data shadowed; struct mc_data missing; - struct mc_data atmosphere; /* TODO rename in absorptivity_loss */ - struct mc_data reflectivity; /* TODO rename in reflectivity_loss */ + struct mc_data absorbed_by_atmosphere; + struct mc_data other_absorbed; struct htable_receiver mc_receivers; /* Per receiver MC */ struct htable_sampled mc_sampled; /* Per sampled instance MC */ diff --git a/src/ssol_material.c b/src/ssol_material.c @@ -64,16 +64,6 @@ ssol_data_ceq(const struct ssol_data* a, const struct ssol_data* b) return i; } -/* Define if the submitted media are *certainly* equals. Refer to the - * check_ssol_data for more details. */ -static INLINE int -media_ceq(const struct ssol_medium* a, const struct ssol_medium* b) -{ - ASSERT(a && b); - return ssol_data_ceq(&a->refractive_index, &b->refractive_index) - && ssol_data_ceq(&a->absorption, &b->absorption); -} - static void shade_normal_default (struct ssol_device* dev, @@ -738,3 +728,12 @@ material_get_next_medium return RES_OK; } +/* Define if the submitted media are *certainly* equals. Refer to the +* check_ssol_data for more details. */ +int +media_ceq(const struct ssol_medium* a, const struct ssol_medium* b) +{ + ASSERT(a && b); + return ssol_data_ceq(&a->refractive_index, &b->refractive_index) + && ssol_data_ceq(&a->absorption, &b->absorption); +} diff --git a/src/ssol_material_c.h b/src/ssol_material_c.h @@ -92,5 +92,8 @@ material_get_next_medium const struct ssol_medium* medium, /* Current mediu */ struct ssol_medium* next_medium); +extern LOCAL_SYM int +media_ceq(const struct ssol_medium* a, const struct ssol_medium* b); + #endif /* SSOL_MATERIAL_C_H */ diff --git a/src/ssol_mc_receiver.c b/src/ssol_mc_receiver.c @@ -59,10 +59,19 @@ ssol_estimator_get_mc_receiver rcv->Name.V = rcv->Name.V > 0 ? rcv->Name.V : 0; \ rcv->Name.SE = sqrt(rcv->Name.V / N); \ } (void)0 - SETUP_MC_RESULT(integrated_irradiance); - SETUP_MC_RESULT(integrated_absorbed_irradiance); - SETUP_MC_RESULT(absorptivity_loss); - SETUP_MC_RESULT(reflectivity_loss); + #define MC_SETUP_ALL { \ + SETUP_MC_RESULT(incoming_flux); \ + SETUP_MC_RESULT(incoming_if_no_atm_loss); \ + SETUP_MC_RESULT(incoming_if_no_field_loss); \ + SETUP_MC_RESULT(incoming_lost_in_atmosphere); \ + SETUP_MC_RESULT(incoming_lost_in_field); \ + SETUP_MC_RESULT(absorbed_flux); \ + SETUP_MC_RESULT(absorbed_if_no_atm_loss); \ + SETUP_MC_RESULT(absorbed_if_no_field_loss); \ + SETUP_MC_RESULT(absorbed_lost_in_atmosphere); \ + SETUP_MC_RESULT(absorbed_lost_in_field); \ + } (void)0 + MC_SETUP_ALL; #undef SETUP_MC_RESULT rcv->mc__ = mc_rcv1; rcv->N__ = estimator->realisation_count; @@ -109,10 +118,7 @@ ssol_mc_shape_get_mc_primitive prim->Name.V = 0; \ prim->Name.SE = 0; \ } (void)0 - SETUP_MC_RESULT(integrated_irradiance); - SETUP_MC_RESULT(integrated_absorbed_irradiance); - SETUP_MC_RESULT(absorptivity_loss); - SETUP_MC_RESULT(reflectivity_loss); + MC_SETUP_ALL; #undef SETUP_MC_RESULT } else { struct s3d_attrib attr; @@ -153,11 +159,9 @@ ssol_mc_shape_get_mc_primitive prim->Name.V /= area*area; \ prim->Name.SE /= area; \ } (void)0 - SETUP_MC_RESULT(integrated_irradiance); - SETUP_MC_RESULT(integrated_absorbed_irradiance); - SETUP_MC_RESULT(absorptivity_loss); - SETUP_MC_RESULT(reflectivity_loss); + MC_SETUP_ALL; #undef SETUP_MC_RESULT + #undef MC_SETUP_ALL } return RES_OK; diff --git a/src/ssol_scene.c b/src/ssol_scene.c @@ -85,6 +85,9 @@ ssol_scene_create res = s3d_scene_create(dev->s3d, &scene->scn_samp); if(res != RES_OK) goto error; + /* default air medium */ + ssol_medium_copy(&scene->air, &SSOL_MEDIUM_VACUUM); + exit: if(out_scene) *out_scene = scene; return res; @@ -231,6 +234,7 @@ ssol_scene_clear(struct ssol_scene* scene) htable_instance_clear(&scene->instances_samp); S3D(scene_clear(scene->scn_rt)); S3D(scene_clear(scene->scn_samp)); + ssol_medium_clear(&scene->air); if(scene->sun) SSOL(scene_detach_sun(scene, scene->sun)); if(scene->atmosphere) SSOL(scene_detach_atmosphere(scene, scene->atmosphere)); return RES_OK; diff --git a/src/ssol_scene_c.h b/src/ssol_scene_c.h @@ -45,6 +45,7 @@ struct ssol_scene { struct ssol_sun* sun; /* Sun of the scene */ struct ssol_atmosphere* atmosphere; /* Atmosphere of the scene */ + struct ssol_medium air; /* Defined according to atmosphere's properties */ struct ssol_device* dev; ref_T ref; diff --git a/src/ssol_solver.c b/src/ssol_solver.c @@ -53,11 +53,11 @@ struct thread_context { struct ssp_rng* rng; struct ssf_bsdf* bsdf; struct mc_data cos_factor; - struct mc_data absorbed; + struct mc_data absorbed_by_receivers; struct mc_data shadowed; struct mc_data missing; - struct mc_data atmosphere; /* TODO rename in absorptivity_loss */ - struct mc_data reflectivity; /* TODO rename in reflectivity_loss */ + struct mc_data absorbed_by_atmosphere; + struct mc_data other_absorbed; struct htable_receiver mc_rcvs; struct htable_sampled mc_samps; struct darray_path paths; /* paths */ @@ -107,11 +107,11 @@ thread_context_copy dst->rng = src->rng; dst->bsdf = src->bsdf; dst->cos_factor = src->cos_factor; - dst->absorbed = src->absorbed; + dst->absorbed_by_receivers = src->absorbed_by_receivers; dst->shadowed = src->shadowed; dst->missing = src->missing; - dst->atmosphere = src->atmosphere; - dst->reflectivity = src->reflectivity; + dst->absorbed_by_atmosphere = src->absorbed_by_atmosphere; + dst->other_absorbed = src->other_absorbed; res = htable_receiver_copy(&dst->mc_rcvs, &src->mc_rcvs); if(res != RES_OK) return res; res = htable_sampled_copy(&dst->mc_samps, &src->mc_samps); @@ -170,12 +170,26 @@ struct point { double dir[3]; float uv[2]; double wl; /* Sampled wavelength */ - /* MC weights, before and after hit */ - double incoming_weight, weight; + const struct ssol_material* material; + /* tmp quantities to compute weights */ + double kabs_at_pt; + double partial_atm, partial_recv, partial_other; + /* MC weights */ + /* Set once */ + double initial_flux; /* the initial flux*/ double cos_factor; /* local cos at the starting point */ - double absorbed_irradiance; /* current hit only */ - double absorptivity_loss; - double reflectivity_loss_before, reflectivity_loss; + /* outgoing weights at previous hit */ + double prev_outgoing_flux; + double prev_outgoing_if_no_atm_loss; + double prev_outgoing_if_no_field_loss; + /* incoming weights at current hit */ + double incoming_flux; + double incoming_if_no_atm_loss; + double incoming_if_no_field_loss; + /* outgoing weights at current hit */ + double outgoing_flux; + double outgoing_if_no_atm_loss; + double outgoing_if_no_field_loss; enum ssol_side_flag side; }; @@ -189,7 +203,9 @@ struct point { {0, 0, 0}, /* Direction */ \ {0, 0}, /* UV */ \ 0, /* Wavelength */ \ - 0, 0, 0, 0, 0, 0, 0,/* MC weights */ \ + NULL, /* Material */ \ + 0, 0, 0, 0, /* tmp values */ \ + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, /* MC weights */ \ SSOL_FRONT /* Side */ \ } static const struct point POINT_NULL = POINT_NULL__; @@ -211,6 +227,7 @@ point_init struct ranst_sun_dir* ran_sun_dir, struct ranst_sun_wl* ran_sun_wl, struct ssp_rng* rng, + struct ssol_medium* current_medium, int* is_lit) { struct ssol_surface_fragment frag; @@ -219,6 +236,7 @@ point_init struct ray_data ray_data = RAY_DATA_NULL; struct ssol_material* mtl; double N[3], Np[3]; + double w0; float dir[3], pos[3], range[2] = { 0, FLT_MAX }; size_t id; res_T res = RES_OK; @@ -277,18 +295,23 @@ point_init /* Initialise the Monte Carlo weight */ if(pt->sshape->shape->type != SHAPE_PUNCHED) { double surface_sun_cos = fabs(d3_dot(Np, pt->dir)); - pt->weight = scn->sun->dni * sampled_area_proxy * surface_sun_cos; + w0 = scn->sun->dni * sampled_area_proxy * surface_sun_cos; pt->cos_factor = surface_sun_cos; } else { double cos_ratio, surface_proxy_cos, surface_sun_cos; surface_proxy_cos = fabs(d3_dot(pt->N, Np)); surface_sun_cos = fabs(d3_dot(Np, pt->dir)); cos_ratio = surface_sun_cos / surface_proxy_cos; - pt->weight = scn->sun->dni * sampled_area_proxy * cos_ratio; + w0 = scn->sun->dni * sampled_area_proxy * cos_ratio; pt->cos_factor = surface_sun_cos; } + + pt->partial_atm = pt->partial_recv = pt->partial_other = 0; + pt->initial_flux = w0; + pt->prev_outgoing_flux = w0; + pt->prev_outgoing_if_no_atm_loss = w0; + pt->prev_outgoing_if_no_field_loss = w0; d3_set(pt->N, N); - pt->absorptivity_loss = pt->reflectivity_loss = 0; ASSERT(d3_dot(pt->N, pt->dir) <= 0); /* Store sampled entity related weights */ @@ -296,6 +319,24 @@ point_init if(res != RES_OK) goto error; pt->mc_samp->nb_samples++; + /* Define the medium in which the sampled point lies */ + pt->material = point_get_material(pt); + switch (pt->material->type) { + case SSOL_MATERIAL_DIELECTRIC: + case SSOL_MATERIAL_THIN_DIELECTRIC: + /* TODO: check sampled face role!!! */ + ssol_medium_copy(current_medium, + (pt->side == SSOL_FRONT) ? + &pt->material->in_medium : &pt->material->out_medium); + break; + case SSOL_MATERIAL_MATTE: + case SSOL_MATERIAL_MIRROR: + case SSOL_MATERIAL_VIRTUAL: + ssol_medium_copy(current_medium, &scn->air); + break; + default: FATAL("Unreachable code\n"); break; + } + /* Initialise the ray data to avoid self intersection */ ray_data.scn = scn; ray_data.prim_from = pt->prim; @@ -369,23 +410,33 @@ point_update_from_hit pt->side = SSOL_BACK; d3_minus(pt->N, pt->N); /* Force the normal to look forward dir */ } + + /* Update material */ + pt->material = point_get_material(pt); +} + +static FINLINE int +point_is_receiver(const struct point* pt) +{ + return (pt->inst->receiver_mask & (int)pt->side) != 0; } static FINLINE res_T point_shade (struct point* pt, struct ssf_bsdf* bsdf, - struct ssol_medium* medium, + const struct ssol_medium* in_medium, + struct ssol_medium* out_medium, struct ssp_rng* rng, double dir[3]) { struct ssol_material* mtl; struct ssol_surface_fragment frag; - double r = 1; + double propagated = 0; double wi[3], N[3], pdf; int type = 0; res_T res; - ASSERT(pt && bsdf && medium && rng && dir); + ASSERT(pt && bsdf && in_medium && out_medium && rng && dir); /* TODO ensure that if `prim' was sampled, then the surface fragment setup * remains valid in *all* situations, i.e. even though the point primitive @@ -403,7 +454,7 @@ point_shade /* Shade the surface fragment */ mtl = point_get_material(pt); SSF(bsdf_clear(bsdf)); - res = material_setup_bsdf(mtl, &frag, pt->wl, medium, 0, bsdf); + res = material_setup_bsdf(mtl, &frag, pt->wl, in_medium, 0, bsdf); if(res != RES_OK) return res; /* Perturbe the normal */ @@ -414,42 +465,45 @@ point_shade d3_minus(wi, pt->dir); if(d3_dot(wi, N) <= 0) { - r = 0; + propagated = 0; } else { double cos_dir_Ng; - r = ssf_bsdf_sample(bsdf, rng, wi, N, dir, &type, &pdf); - ASSERT(0 <= r && r <= 1); + propagated = ssf_bsdf_sample(bsdf, rng, wi, N, dir, &type, &pdf); + ASSERT(0 <= propagated && propagated <= 1); /* Due to the perturbed normal, the sampled direction may point in the * wrong direction wrt the sampled BSDF component. */ cos_dir_Ng = d3_dot(frag.Ng, dir); if((cos_dir_Ng > 0 && (type & SSF_TRANSMISSION)) || (cos_dir_Ng < 0 && (type & SSF_REFLECTION))) { - r = 0; + propagated = 0; } } - pt->incoming_weight = pt->weight; - pt->reflectivity_loss_before = pt->reflectivity_loss; - pt->absorbed_irradiance = (1 - r) * pt->weight; - pt->reflectivity_loss += pt->absorbed_irradiance; - pt->weight = pt->incoming_weight - pt->absorbed_irradiance; - - if(type & SSF_TRANSMISSION) material_get_next_medium(mtl, medium, medium); + pt->kabs_at_pt = (1 - propagated); + pt->outgoing_flux = pt->incoming_flux * propagated; + pt->outgoing_if_no_atm_loss = pt->incoming_if_no_atm_loss * propagated; + pt->outgoing_if_no_field_loss = point_is_receiver(pt) ? + pt->incoming_if_no_field_loss * propagated : pt->incoming_if_no_field_loss; + + if(type & SSF_TRANSMISSION) { + material_get_next_medium(mtl, in_medium, out_medium); + } else { + ssol_medium_copy(out_medium, in_medium); + } return RES_OK; } static FINLINE void -point_hit_virtual(struct point* pt) -{ - pt->absorbed_irradiance = 0; - pt->incoming_weight = pt->weight; - pt->reflectivity_loss_before = pt->reflectivity_loss; -} - -static FINLINE int -point_is_receiver(const struct point* pt) +point_hit_virtual + (struct point* pt, + const struct ssol_medium* in_medium, + struct ssol_medium* out_medium) { - return (pt->inst->receiver_mask & (int)pt->side) != 0; + pt->kabs_at_pt = 0; + pt->outgoing_flux = pt->incoming_flux; + pt->outgoing_if_no_atm_loss = pt->incoming_if_no_atm_loss; + pt->outgoing_if_no_field_loss = pt->incoming_if_no_field_loss; + ssol_medium_copy(out_medium, in_medium); } static FINLINE int32_t @@ -480,7 +534,7 @@ point_dump f3_set_d3(out.in_dir, pt->dir); f3_set_d3(out.normal, pt->N); f2_set(out.uv, pt->uv); - out.weight = pt->weight; + out.weight = pt->outgoing_flux; n = fwrite(&out, sizeof(out), 1, stream); return n != 1 ? RES_IO_ERR : RES_OK; } @@ -531,10 +585,19 @@ accum_mc_receivers_1side dst->Name.weight += src->Name.weight; \ dst->Name.sqr_weight += src->Name.sqr_weight; \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance); - ACCUM_WEIGHT(integrated_absorbed_irradiance); - ACCUM_WEIGHT(absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss); + #define ACCUM_ALL { \ + ACCUM_WEIGHT(incoming_flux); \ + ACCUM_WEIGHT(incoming_if_no_atm_loss); \ + ACCUM_WEIGHT(incoming_lost_in_field); \ + ACCUM_WEIGHT(incoming_lost_in_atmosphere); \ + ACCUM_WEIGHT(incoming_if_no_field_loss); \ + ACCUM_WEIGHT(absorbed_flux); \ + ACCUM_WEIGHT(absorbed_if_no_atm_loss); \ + ACCUM_WEIGHT(absorbed_if_no_field_loss); \ + ACCUM_WEIGHT(absorbed_lost_in_field); \ + ACCUM_WEIGHT(absorbed_lost_in_atmosphere); \ + } (void)0 + ACCUM_ALL; #undef ACCUM_WEIGHT /* Merge the per shape MC */ @@ -568,11 +631,9 @@ accum_mc_receivers_1side mc_prim1_dst->Name.weight += mc_prim1_src->Name.weight; \ mc_prim1_dst->Name.sqr_weight += mc_prim1_src->Name.sqr_weight; \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance); - ACCUM_WEIGHT(integrated_absorbed_irradiance); - ACCUM_WEIGHT(absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss); + ACCUM_ALL; #undef ACCUM_WEIGHT + #undef ACCUM_ALL htable_prim2mc_iterator_next(&it_prim); } @@ -644,7 +705,7 @@ error: */ static res_T update_mc - (const struct point* pt, + (struct point* pt, const size_t irealisation, const size_t ibounce, struct thread_context* thread_ctx, @@ -658,13 +719,7 @@ update_mc res = point_dump(pt, irealisation, ibounce, output); if(res != RES_OK) goto error; - /* Global MC accumulation */ - #define ACCUM_WEIGHT(Res, W) { \ - Res.weight += (W); \ - Res.sqr_weight += (W)*(W); \ - } (void)0 - ACCUM_WEIGHT(thread_ctx->absorbed, pt->absorbed_irradiance); - #undef ACCUM_WEIGHT + pt->partial_recv += pt->incoming_flux - pt->outgoing_flux; /* Per receiver MC accumulation */ res = get_mc_receiver_1side(&thread_ctx->mc_rcvs, pt->inst, pt->side, &mc_rcv1); @@ -674,10 +729,25 @@ update_mc mc_rcv1->Name.weight += (W); \ mc_rcv1->Name.sqr_weight += (W)*(W); \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance, pt->incoming_weight); - ACCUM_WEIGHT(integrated_absorbed_irradiance, pt->absorbed_irradiance); - ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss_before); + #define ACCUM_ALL { \ + ACCUM_WEIGHT(incoming_flux, pt->incoming_flux); \ + ACCUM_WEIGHT(incoming_if_no_atm_loss, pt->incoming_if_no_atm_loss); \ + ACCUM_WEIGHT(incoming_if_no_field_loss, pt->incoming_if_no_field_loss); \ + ACCUM_WEIGHT(incoming_lost_in_field, \ + pt->incoming_if_no_field_loss - pt->incoming_flux); \ + ACCUM_WEIGHT(incoming_lost_in_atmosphere, \ + pt->incoming_if_no_atm_loss - pt->incoming_flux); \ + ACCUM_WEIGHT(absorbed_flux, pt->incoming_flux * pt->kabs_at_pt); \ + ACCUM_WEIGHT(absorbed_if_no_atm_loss, \ + pt->incoming_if_no_atm_loss * pt->kabs_at_pt); \ + ACCUM_WEIGHT(absorbed_if_no_field_loss, \ + pt->incoming_if_no_field_loss * pt->kabs_at_pt); \ + ACCUM_WEIGHT(absorbed_lost_in_field, \ + (pt->incoming_if_no_field_loss - pt->incoming_flux) * pt->kabs_at_pt); \ + ACCUM_WEIGHT(absorbed_lost_in_atmosphere, \ + (pt->incoming_if_no_atm_loss - pt->incoming_flux) * pt->kabs_at_pt); \ + } (void)0 + ACCUM_ALL; #undef ACCUM_WEIGHT /* Per-sampled/receiver MC accumulation */ @@ -689,10 +759,7 @@ update_mc mc_samp_x_rcv1->Name.weight += (W); \ mc_samp_x_rcv1->Name.sqr_weight += (W)*(W); \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance, pt->incoming_weight); - ACCUM_WEIGHT(integrated_absorbed_irradiance, pt->absorbed_irradiance); - ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss_before); + ACCUM_ALL; #undef ACCUM_WEIGHT /* Per primitive receiver MC accumulation */ @@ -710,11 +777,9 @@ update_mc mc_prim1->Name.weight += (W); \ mc_prim1->Name.sqr_weight += (W)*(W); \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance, pt->incoming_weight); - ACCUM_WEIGHT(integrated_absorbed_irradiance, pt->absorbed_irradiance); - ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss_before); + ACCUM_ALL; #undef ACCUM_WEIGHT + #undef ACCUM_ALL } exit: @@ -737,7 +802,8 @@ trace_radiative_path FILE* output) /* May be NULL */ { struct path path; - struct ssol_medium medium = SSOL_MEDIUM_VACUUM; + struct ssol_medium in_medium = SSOL_MEDIUM_VACUUM; + struct ssol_medium out_medium = SSOL_MEDIUM_VACUUM; struct s3d_hit hit = S3D_HIT_NULL; struct point pt; float org[3], dir[3], range[2] = { 0, FLT_MAX }; @@ -751,7 +817,8 @@ trace_radiative_path /* Find a new starting point of the radiative random walk */ res = point_init(&pt, sampled_area_proxy, scn, &thread_ctx->mc_samps, - view_samp, view_rt, ran_sun_dir, ran_sun_wl, thread_ctx->rng, &is_lit); + view_samp, view_rt, ran_sun_dir, ran_sun_wl, thread_ctx->rng, + &in_medium, &is_lit); if(res != RES_OK) goto error; if(tracker) { @@ -765,9 +832,8 @@ trace_radiative_path if(res != RES_OK) goto error; } - /* Register the init position onto the sampled geometry */ - res = path_add_vertex(&path, pt.pos, pt.weight); + res = path_add_vertex(&path, pt.pos, pt.initial_flux); if(res != RES_OK) goto error; } @@ -775,60 +841,67 @@ trace_radiative_path Res.weight += (W); \ Res.sqr_weight += (W)*(W); \ } (void)0 - ACCUM_WEIGHT(thread_ctx->cos_factor, pt.cos_factor); - ACCUM_WEIGHT(pt.mc_samp->cos_factor, pt.cos_factor); if(!is_lit) { /* The starting point is not lit */ - ACCUM_WEIGHT(pt.mc_samp->shadowed, pt.weight); - ACCUM_WEIGHT(thread_ctx->shadowed, pt.weight); + ACCUM_WEIGHT(pt.mc_samp->shadowed, pt.initial_flux); + ACCUM_WEIGHT(thread_ctx->shadowed, pt.initial_flux); if(tracker) path.type = SSOL_PATH_SHADOW; } else { - struct ssol_material* mtl = point_get_material(&pt); - - /* Define the medium in which the sampled point lies */ - switch(mtl->type) { - case SSOL_MATERIAL_DIELECTRIC: - case SSOL_MATERIAL_THIN_DIELECTRIC: - ssol_medium_copy(&medium, &mtl->out_medium); - break; - default: - if(scn->atmosphere) - ssol_data_copy(&medium.absorption, &scn->atmosphere->absorption); - break; - } - /* Setup the ray as if it starts from the current point position in order * to handle the points that start from a virtual material */ f3_set_d3(org, pt.pos); f3_set_d3(dir, pt.dir); - hit.distance = 0; + hit.distance = 0; /* first loop has no atmospheric absorption */ for(;;) { /* Here we go for the radiative random walk */ + const int in_atm = media_ceq(&in_medium, &scn->air); + const int hit_receiver = point_is_receiver(&pt); + const int hit_virtual = pt.material->type == SSOL_MATERIAL_VIRTUAL; + int last_segment = 0; struct ray_data ray_data = RAY_DATA_NULL; - double absorption; + double trans = 1; + + /* Compute medium absorption along the incoming segment. */ + if(hit.distance > 0) { + const double kabs = ssol_data_get_value(&in_medium.absorption, pt.wl); + ASSERT(0 <= kabs && kabs <= 1); + if(kabs > 0) { + trans = exp(-kabs * hit.distance); + } + } + pt.incoming_flux = pt.prev_outgoing_flux * trans; + pt.incoming_if_no_atm_loss = in_atm ? + pt.prev_outgoing_if_no_atm_loss : pt.prev_outgoing_if_no_atm_loss * trans; + pt.incoming_if_no_field_loss = (!in_atm) ? + pt.prev_outgoing_if_no_field_loss : pt.prev_outgoing_if_no_field_loss * trans; /* Compute interaction with material */ - mtl = point_get_material(&pt); - if(mtl->type == SSOL_MATERIAL_VIRTUAL) { - point_hit_virtual(&pt); + if(hit_virtual) { + point_hit_virtual(&pt, &in_medium, &out_medium); } else { - /* Modulate the point weight wrt its scattering functions and generate - * an outgoing direction */ - res = point_shade(&pt, thread_ctx->bsdf, &medium, thread_ctx->rng, pt.dir); + /* Modulate the point weights wrt its scattering functions and generate + * an outgoing direction and set out_medium accordingly */ + res = point_shade(&pt, thread_ctx->bsdf, &in_medium, &out_medium, + thread_ctx->rng, pt.dir); if(res != RES_OK) goto error; } - if(point_is_receiver(&pt)) { + /* If receiver update MC results */ + if(hit_receiver) { hit_a_receiver = 1; res = update_mc(&pt, path_id, depth, thread_ctx, output); if(res != RES_OK) goto error; + } else { + pt.partial_other += pt.incoming_flux * pt.kabs_at_pt; } - /* Stop the radiative random walk */ - if(pt.weight == 0) break; + /* Stop the radiative random walk if no more flux */ + if(!pt.outgoing_flux) { + break; + } /* Setup new ray parameters */ - if(mtl->type == SSOL_MATERIAL_VIRTUAL) { + if(hit_virtual) { /* Note that for Virtual materials, the ray parameters 'org' & 'dir' * are not updated to ensure that it pursues its traversal without any * accuracy issue */ @@ -851,42 +924,63 @@ trace_radiative_path ray_data.range_min = range[0]; ray_data.dst = FLT_MAX; S3D(scene_view_trace_ray(view_rt, org, dir, range, &ray_data, &hit)); - if(S3D_HIT_NONE(&hit)) { + if(S3D_HIT_NONE(&hit)) { /* The ray is lost! */ /* Add the point of the last path segment going to the infinite */ - if (tracker && tracker->infinite_ray_length > 0) { + if(tracker && tracker->infinite_ray_length > 0) { double pos[3], wi[3]; d3_set_f3(wi, dir); d3_add(pos, pt.pos, d3_muld(wi, wi, tracker->infinite_ray_length)); - res = path_add_vertex(&path, pos, pt.weight); + res = path_add_vertex(&path, pos, pt.outgoing_flux); if (res != RES_OK) goto error; } - break; + last_segment = 1; /* Path reached its last segment */ + ASSERT(in_atm); } - depth += mtl->type != SSOL_MATERIAL_VIRTUAL; - - /* Take into account the medium attenuation */ - absorption = ssol_data_get_value(&medium.absorption, pt.wl); - if(absorption > 0 && hit.distance > 0) { - const double transmissivity = exp(-absorption * hit.distance); - ASSERT(0 < transmissivity && transmissivity <= 1); - pt.absorptivity_loss += (1 - transmissivity) * pt.weight; - pt.weight *= transmissivity; + /* Don't change prev_outgoing weigths nor record segment absorption until + * a non-virtual material is hit or this segment is the last one. + * This is because propagation is restarted from the same origin until + * a non-virtual material is hit or no further hit can be found. */ + if(last_segment || !hit_virtual) { + if(in_atm) { + pt.partial_atm += pt.prev_outgoing_flux - pt.incoming_flux; + } else { + pt.partial_other += pt.prev_outgoing_flux - pt.incoming_flux; + } + if(last_segment) { + break; + } + pt.prev_outgoing_flux = pt.outgoing_flux; + pt.prev_outgoing_if_no_atm_loss = pt.outgoing_if_no_atm_loss; + pt.prev_outgoing_if_no_field_loss = pt.outgoing_if_no_field_loss; } + depth += !hit_virtual; + ASSERT(depth < 100); /* FIXME: create a true cancel path for this MC sample */ + /* Update the point */ point_update_from_hit(&pt, scn, org, dir, &hit, &ray_data); if(tracker) { - res = path_add_vertex(&path, pt.pos, pt.weight); + res = path_add_vertex(&path, pt.pos, pt.outgoing_flux); if (res != RES_OK) goto error; } + + ssol_medium_copy(&in_medium, &out_medium); } - /* Handle the overall absorptivity losses along the radiative path, the - * total of reflectivity losses and register the remaining power as lost */ - ACCUM_WEIGHT(thread_ctx->atmosphere, pt.absorptivity_loss); - ACCUM_WEIGHT(thread_ctx->reflectivity, pt.reflectivity_loss); - ACCUM_WEIGHT(thread_ctx->missing, pt.weight); + /* Check conservation of energy */ + ASSERT(depth * pt.initial_flux * DBL_EPSILON >= + fabs(pt.initial_flux - + (pt.outgoing_flux + pt.partial_recv + pt.partial_atm + pt.partial_other))); + + /* Now that the sample ends successfully, + * record MC weights and register the remaining power as missing */ + ACCUM_WEIGHT(pt.mc_samp->cos_factor, pt.cos_factor); + ACCUM_WEIGHT(thread_ctx->cos_factor, pt.cos_factor); + ACCUM_WEIGHT(thread_ctx->missing, pt.outgoing_flux); + ACCUM_WEIGHT(thread_ctx->absorbed_by_receivers, pt.partial_recv); + ACCUM_WEIGHT(thread_ctx->absorbed_by_atmosphere, pt.partial_atm); + ACCUM_WEIGHT(thread_ctx->other_absorbed, pt.partial_other); #undef ACCUM_WEIGHT if(tracker) { @@ -899,7 +993,8 @@ trace_radiative_path if(res != RES_OK) goto error; } exit: - ssol_medium_clear(&medium); + ssol_medium_clear(&in_medium); + ssol_medium_clear(&out_medium); if(tracker) path_release(&path); return res; error: @@ -957,6 +1052,12 @@ ssol_solve res = scene_check(scn, FUNC_NAME); if(res != RES_OK) goto error; + + /* init air properties */ + if(scn->atmosphere) + ssol_data_copy(&scn->air.absorption, &scn->atmosphere->absorption); + else + ssol_data_copy(&scn->air.absorption, &SSOL_MEDIUM_VACUUM.absorption); /* Create data structures shared by all threads */ res = scene_create_s3d_views(scn, &view_rt, &view_samp, &sampled_area, @@ -1034,11 +1135,11 @@ ssol_solve estimator->Name.sqr_weight += thread_ctx->Name.sqr_weight; \ } (void)0 ACCUM_WEIGHT(cos_factor); - ACCUM_WEIGHT(absorbed); + ACCUM_WEIGHT(absorbed_by_receivers); ACCUM_WEIGHT(shadowed); ACCUM_WEIGHT(missing); - ACCUM_WEIGHT(atmosphere); - ACCUM_WEIGHT(reflectivity); + ACCUM_WEIGHT(absorbed_by_atmosphere); + ACCUM_WEIGHT(other_absorbed); estimator->realisation_count += thread_ctx->realisation_count; #undef ACCUM_WEIGHT } diff --git a/src/test_ssol_by_receiver_integration.c b/src/test_ssol_by_receiver_integration.c @@ -137,24 +137,24 @@ main(int argc, char** argv) CHECK(ssol_solve(scene, rng, N__, 0, NULL, &estimator1), RES_OK); CHECK(GET_MC_RCV(estimator1, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); CHECK(ssol_instance_set_receiver(heliostat, SSOL_FRONT, 0), RES_OK); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S_DNI_cos, S_DNI_cos * 2e-1), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.E, S_DNI_cos, S_DNI_cos * 2e-1), 1); CHECK(ssol_solve(scene, rng, 8 * N__, 0, NULL, &estimator2), RES_OK); CHECK(GET_MC_RCV(estimator2, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S_DNI_cos, S_DNI_cos * 5e-2), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, S_DNI_cos, S_DNI_cos * 5e-2), 1); CHECK(ssol_estimator_ref_put(estimator1), RES_OK); CHECK(ssol_solve(scene, rng, 3 * N__, 0, NULL, &estimator1), RES_OK); CHECK(GET_MC_RCV(estimator1, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S_DNI_cos, S_DNI_cos * 1e-1), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, S_DNI_cos, S_DNI_cos * 1e-1), 1); CHECK(GET_MC_SAMP_X_RCV(estimator1, heliostat, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(heliostat=>target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S_DNI_cos, S_DNI_cos * 1e-1), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, S_DNI_cos, S_DNI_cos * 1e-1), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat), RES_OK); diff --git a/src/test_ssol_solver1.c b/src/test_ssol_solver1.c @@ -118,7 +118,8 @@ main(int argc, char** argv) CHECK(ssol_sun_create_directional(dev, &sun), RES_OK); CHECK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)), RES_OK); CHECK(ssol_sun_set_spectrum(sun, spectrum), RES_OK); - CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); +#define DNI 1000 + CHECK(ssol_sun_set_dni(sun, DNI), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_solve(NULL, rng, 10, 0, NULL, &estimator), RES_BAD_ARG); @@ -263,7 +264,7 @@ main(int argc, char** argv) /* Sun with no spectrum */ CHECK(ssol_sun_create_directional(dev, &sun), RES_OK); CHECK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)), RES_OK); - CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); + CHECK(ssol_sun_set_dni(sun, DNI), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG); CHECK(ssol_scene_detach_sun(scene, sun), RES_OK); @@ -275,7 +276,7 @@ main(int argc, char** argv) CHECK(ssol_sun_set_spectrum(sun, spectrum), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG); - CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); + CHECK(ssol_sun_set_dni(sun, DNI), RES_OK); /* No receiver in scene */ CHECK(ssol_instance_set_receiver(heliostat, 0, 0), RES_OK); @@ -302,7 +303,6 @@ main(int argc, char** argv) CHECK(fcount, 0); printf("Ir = %g +/- %g; ", m, std); #define COS cos(PI / 4) -#define DNI 1000 #define DNI_cos (DNI * COS) CHECK(eq_eps(m, 4 * DNI_cos, MMAX(4 * DNI_cos * 1e-2, 2*std)), 1); #define SQR(x) ((x)*(x)) @@ -310,9 +310,7 @@ main(int argc, char** argv) CHECK(eq_eps(std, dbl, dbl*1e-2), 1); /* Target was sampled but shadowed by secondary */ CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); - printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); - printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); + PRINT_GLOBAL(mc_global); CHECK(eq_eps(mc_global.shadowed.E, m, 2 * dbl), 1); CHECK(eq_eps(mc_global.missing.E, 2*m, 2*mc_global.missing.SE), 1); CHECK(GET_MC_RCV(NULL, NULL, SSOL_BACK, NULL), RES_BAD_ARG); @@ -332,9 +330,9 @@ main(int argc, char** argv) CHECK(GET_MC_RCV(NULL, target, SSOL_FRONT, &mc_rcv), RES_BAD_ARG); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); CHECK(ssol_estimator_ref_put(estimator), RES_OK); /* Sample primary mirror only; variance is low */ @@ -351,17 +349,15 @@ main(int argc, char** argv) CHECK(eq_eps(m, 4 * DNI_cos, MMAX(4 * DNI_cos * 1e-2, std)), 1); CHECK(eq_eps(std, 0, 1e-4), 1); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); - printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); - printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); + PRINT_GLOBAL(mc_global); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); CHECK(eq_eps(mc_global.missing.E, m, 1e-4), 1); CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); CHECK(ssol_estimator_ref_put(estimator), RES_OK); /* Check atmosphere model; with no absorption result is unchanged */ @@ -383,17 +379,15 @@ main(int argc, char** argv) CHECK(ssol_scene_detach_atmosphere(scene, atm), RES_OK); CHECK(ssol_atmosphere_ref_put(atm), RES_OK); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); - printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); - printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); + PRINT_GLOBAL(mc_global); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); CHECK(eq_eps(mc_global.missing.E, m, 1e-4), 1); CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(ssol_estimator_ref_put(estimator), RES_OK); @@ -433,46 +427,23 @@ main(int argc, char** argv) MMAX(4 * K * DNI_cos * 1e-1, a_std)), 1); CHECK(eq_eps(a_std, 0, 1e-4), 1); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); - printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); - printf("Atmosphere = %g +/- %g; ", mc_global.atmosphere.E, mc_global.atmosphere.SE); - printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); + PRINT_GLOBAL(mc_global); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E + mc_global.atmosphere.E + mc_global.absorbed.E, - m, 1e-4), 1); + CHECK(eq_eps( + mc_global.missing.E + mc_global.shadowed.E + mc_global.absorbed_by_receivers.E + + mc_global.absorbed_by_atmosphere.E + mc_global.other_absorbed.E, + m, + 1e-4), 1); CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); - printf - ("\tIr(target) = %g +/- %g (%.2g %%)\n", - mc_rcv.integrated_irradiance.E, - mc_rcv.integrated_irradiance.SE, - 100 * mc_rcv.integrated_irradiance.E / m); - printf - ("\tAtmospheric Loss(target) = %g +/- %g (%.2g %%)\n", - mc_rcv.absorptivity_loss.E, - mc_rcv.absorptivity_loss.SE, - 100 * mc_rcv.absorptivity_loss.E / m); - printf - ("\tReflectivity Loss(target) = %g +/- %g (%.2g %%)\n", - mc_rcv.reflectivity_loss.E, - mc_rcv.reflectivity_loss.SE, - 100 * mc_rcv.reflectivity_loss.E / m); - + PRINT_RCV(mc_rcv); CHECK(ssol_estimator_get_sampled_count(estimator, &scount), RES_OK); CHECK(ssol_estimator_get_mc_sampled(estimator, heliostat, &sampled), RES_BAD_ARG); CHECK(ssol_estimator_get_mc_sampled(estimator, heliostat2, &sampled), RES_OK); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, a_m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, a_std, 1e-4), 1); - CHECK(eq_eps - ( mc_rcv.integrated_irradiance.E - + mc_rcv.absorptivity_loss.E - + mc_rcv.reflectivity_loss.E, m, 1e-8), 1); - CHECK(eq_eps - ( mc_rcv.integrated_irradiance.E - + mc_rcv.absorptivity_loss.E - + mc_rcv.reflectivity_loss.E - + (1 - mc_global.cos_factor.E) * 4 * DNI, 4 * DNI, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.E, a_m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, a_std, 1e-4), 1); + CHECK(ssol_mc_receiver_get_mc_shape(NULL, NULL, NULL), RES_BAD_ARG); CHECK(ssol_mc_receiver_get_mc_shape(&mc_rcv, NULL, NULL), RES_BAD_ARG); CHECK(ssol_mc_receiver_get_mc_shape(NULL, square, NULL), RES_BAD_ARG); @@ -505,7 +476,7 @@ main(int argc, char** argv) area = d3_len(d3_cross(N, d3_sub(E1, v1, v0), d3_sub(E2, v2, v0))) * 0.5; CHECK(ssol_mc_shape_get_mc_primitive(&mc_shape, (unsigned)i, &mc_prim), RES_OK); - dbl += mc_prim.integrated_irradiance.E * area; + dbl += mc_prim.incoming_flux.E * area; } CHECK(eq_eps(dbl, a_m, 1.e-6), 1); @@ -523,7 +494,7 @@ main(int argc, char** argv) CHECK(ssol_sun_create_directional(dev, &sun_mono), RES_OK); CHECK(ssol_sun_set_direction(sun_mono, d3(dir, 1, 0, -1)), RES_OK); CHECK(ssol_sun_set_spectrum(sun_mono, spectrum), RES_OK); - CHECK(ssol_sun_set_dni(sun_mono, 1000), RES_OK); + CHECK(ssol_sun_set_dni(sun_mono, DNI), RES_OK); CHECK(ssol_scene_detach_sun(scene, sun), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun_mono), RES_OK); ka[1] = 0.2; ka[0] = ka[2] = 0.1; @@ -547,17 +518,16 @@ main(int argc, char** argv) CHECK(eq_eps(m, 4 * K2 * DNI_cos, MMAX(4 * K2 * DNI_cos * 1e-4, std)), 1); CHECK(eq_eps(std, 0, 1e-4), 1); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); - printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); - printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); + PRINT_GLOBAL(mc_global); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); CHECK(eq_eps(mc_global.missing.E, m, 1e-4), 1); CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + PRINT_RCV(mc_rcv); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat2), RES_OK); diff --git a/src/test_ssol_solver2.c b/src/test_ssol_solver2.c @@ -202,12 +202,12 @@ main(int argc, char** argv) CHECK(GET_MC_RCV(estimator, heliostat1, SSOL_BACK, &mc_rcv), RES_BAD_ARG); CHECK(GET_MC_RCV(estimator, secondary, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(secondary) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat1), RES_OK); diff --git a/src/test_ssol_solver2b.c b/src/test_ssol_solver2b.c @@ -206,9 +206,9 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_receiver (estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK); CHECK(count, 0); diff --git a/src/test_ssol_solver3.c b/src/test_ssol_solver3.c @@ -161,9 +161,9 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_receiver (estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK); CHECK(count, 0); CHECK(ssol_instance_get_area(heliostat, &area), RES_OK); diff --git a/src/test_ssol_solver4.c b/src/test_ssol_solver4.c @@ -170,14 +170,14 @@ main(int argc, char** argv) CHECK(eq_eps(mc_global.missing.E, 400 * DNI_cos, 1e-2), 1); /* nothing absorbed */ CHECK(GET_MC_RCV(estimator, target1, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target1) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m1, 1e-2), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std1, 1e-2), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m1, 1e-2), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std1, 1e-2), 1); CHECK(GET_MC_RCV(estimator, target2, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target2) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m2, 1e-2), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std2, 1e-2), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m2, 1e-2), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std2, 1e-2), 1); CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK); CHECK(count, 0); diff --git a/src/test_ssol_solver5.c b/src/test_ssol_solver5.c @@ -158,9 +158,9 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_receiver (estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, m, 1e-8), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, std, 1e-4), 1); CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK); CHECK(count, 0); diff --git a/src/test_ssol_solver6.c b/src/test_ssol_solver6.c @@ -181,24 +181,21 @@ main(int argc, char** argv) CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK); CHECK(fclose(tmp), 0); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); - printf("Absorbed = %g +/- %g\n", mc_global.absorbed.E, mc_global.absorbed.SE); - printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE); - printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); + PRINT_GLOBAL(mc_global); CHECK(eq_eps(mc_global.shadowed.E, 100000, 2 * 100000/sqrt(N__)), 1); CHECK(eq_eps(mc_global.missing.E, 0, 0), 1); CHECK(GET_MC_RCV(estimator, target1, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target1) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, 100000, 2*100000/sqrt(N__)), 1); - CHECK(mc_rcv.integrated_irradiance.E, mc_rcv.integrated_absorbed_irradiance.E); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, 100000, 2*100000/sqrt(N__)), 1); + CHECK(mc_rcv.incoming_flux.E, mc_rcv.absorbed_flux.E); CHECK(GET_MC_RCV(estimator, target2, SSOL_BACK, &mc_rcv), RES_OK); printf("Ir(target2) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, 0, 1), 1); - CHECK(mc_rcv.integrated_irradiance.E, mc_rcv.integrated_absorbed_irradiance.E); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, 0, 1), 1); + CHECK(mc_rcv.incoming_flux.E, mc_rcv.absorbed_flux.E); /* Free data */ CHECK(ssol_instance_ref_put(heliostat1), RES_OK); diff --git a/src/test_ssol_solver7.c b/src/test_ssol_solver7.c @@ -200,22 +200,18 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_sampled_area(estimator, &area), RES_OK); printf("Total = %g\n", area * DNI_cos); CHECK(eq_eps(area * DNI_cos, TOTAL, TOTAL * 1e-4), 1); - printf("Absorbed = %g +/- %g\n", mc_global.absorbed.E, mc_global.absorbed.SE); - printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); - printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE); - CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE); + PRINT_GLOBAL(mc_global); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Abs(target1) = %g +/- %g\n", - mc_rcv.integrated_absorbed_irradiance.E, - mc_rcv.integrated_absorbed_irradiance.SE); + mc_rcv.absorbed_flux.E, + mc_rcv.absorbed_flux.SE); printf("Ir(target1) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, TOTAL, TOTAL * 1e-4), 1); - CHECK(eq_eps(mc_rcv.integrated_absorbed_irradiance.E, + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, TOTAL, TOTAL * 1e-4), 1); + CHECK(eq_eps(mc_rcv.absorbed_flux.E, (1 - REFLECTIVITY) * TOTAL, (1 - REFLECTIVITY) *TOTAL * 1e-4), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, 0, 1e-2), 1); + CHECK(eq_eps(mc_rcv.incoming_flux.SE, 0, 1e-2), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat), RES_OK); diff --git a/src/test_ssol_solver8.c b/src/test_ssol_solver8.c @@ -159,9 +159,9 @@ main(int argc, char** argv) 3 * mc_global.missing.SE), 1); /* nothing absorbed */ CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target1) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S * DNI, - 2 * mc_rcv.integrated_irradiance.SE), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, S * DNI, + 2 * mc_rcv.incoming_flux.SE), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat), RES_OK); diff --git a/src/test_ssol_solver9.c b/src/test_ssol_solver9.c @@ -154,9 +154,9 @@ main(int argc, char** argv) 3 * mc_global.missing.SE), 1); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target1) = %g +/- %g\n", - mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, MMIN(DNI_S, DNI_TGT_S), - 3 * mc_rcv.integrated_irradiance.SE), 1); + mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE); + CHECK(eq_eps(mc_rcv.incoming_flux.E, MMIN(DNI_S, DNI_TGT_S), + 3 * mc_rcv.incoming_flux.SE), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat), RES_OK); diff --git a/src/test_ssol_utils.h b/src/test_ssol_utils.h @@ -38,4 +38,47 @@ pp_sum double* mean, double* std); +#define PRINT_GLOBAL(Mc) { \ + printf("Shadows = %g +/- %g; ", (Mc).shadowed.E, (Mc).shadowed.SE); \ + printf("Missing = %g +/- %g; ", (Mc).missing.E, (Mc).missing.SE); \ + printf("Receivers = %g +/- %g; ", \ + (Mc).absorbed_by_receivers.E, (Mc).absorbed_by_receivers.SE); \ + printf("Atmosphere = %g +/- %g; ", \ + (Mc).absorbed_by_atmosphere.E, (Mc).absorbed_by_atmosphere.SE); \ + printf("Other absorbed = %g +/- %g; ", \ + (Mc).other_absorbed.E, (Mc).other_absorbed.SE); \ + printf("Cos = %g +/- %g\n", (Mc).cos_factor.E, (Mc).cos_factor.SE); \ +} (void)0 + +#define PRINT_RCV(Rcv) { \ + printf("\tIncoming flux(target) = %g +/- %g \n", \ + (Rcv).incoming_flux.E, (Rcv).incoming_flux.SE); \ + printf("\tIncoming flux wo Atmosphere(target) = %g +/- %g (%.2g %%)\n", \ + (Rcv).incoming_if_no_atm_loss.E, (Rcv).incoming_if_no_atm_loss.SE, \ + 100 * (Rcv).incoming_if_no_atm_loss.E / (Rcv).incoming_flux.E); \ + printf("\tIncoming flux wo Field Loss(target) = %g +/- %g (%.2g %%)\n", \ + (Rcv).incoming_if_no_field_loss.E, (Rcv).incoming_if_no_field_loss.SE, \ + 100 * (Rcv).incoming_if_no_field_loss.E / (Rcv).incoming_flux.E); \ + printf("\tAtmospheric Loss on Incoming(target) = %g +/- %g (%.2g %%)\n", \ + (Rcv).incoming_lost_in_atmosphere.E, (Rcv).incoming_lost_in_atmosphere.SE,\ + 100 * (Rcv).incoming_lost_in_atmosphere.E / (Rcv).incoming_flux.E); \ + printf("\tOptical Field Loss(target) on Incoming = %g +/- %g (%.2g %%)\n", \ + (Rcv).incoming_lost_in_field.E, (Rcv).incoming_lost_in_field.SE, \ + 100 * (Rcv).incoming_lost_in_field.E / (Rcv).incoming_flux.E); \ + printf("\tAbsorbed flux(target) = %g +/- %g \n", \ + (Rcv).absorbed_flux.E, (Rcv).absorbed_flux.SE); \ + printf("\tAbsorbed flux wo Atmosphere(target) = %g +/- %g (%.2g %%)\n", \ + (Rcv).absorbed_if_no_atm_loss.E, (Rcv).absorbed_if_no_atm_loss.SE, \ + 100 * (Rcv).absorbed_if_no_atm_loss.E / (Rcv).absorbed_flux.E); \ + printf("\tAbsorbed flux wo Field Loss(target) = %g +/- %g (%.2g %%)\n", \ + (Rcv).absorbed_if_no_field_loss.E, (Rcv).absorbed_if_no_field_loss.SE, \ + 100 * (Rcv).absorbed_if_no_field_loss.E / (Rcv).absorbed_flux.E); \ + printf("\tAtmospheric Loss(target) on Absorbed = %g +/- %g (%.2g %%)\n", \ + (Rcv).absorbed_lost_in_atmosphere.E, (Rcv).absorbed_lost_in_atmosphere.SE,\ + 100 * (Rcv).absorbed_lost_in_atmosphere.E / (Rcv).absorbed_flux.E); \ + printf("\tOptical Field Loss(target) on Absorbed = %g +/- %g (%.2g %%)\n", \ + (Rcv).absorbed_lost_in_field.E, (Rcv).absorbed_lost_in_field.SE, \ + 100 * (Rcv).absorbed_lost_in_field.E / (Rcv).absorbed_flux.E); \ +} (void)0 + #endif /* TEST_SSOL_UTILS_H */