test_sgf_square.c (6398B)
1 /* Copyright (C) 2021, 2024 |Meso|Star> (contact@meso-star.com) 2 * Copyright (C) 2015-2018 EDF S.A., France (syrthes-support@edf.fr) 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 "sgf.h" 18 #include "test_sgf_utils.h" 19 20 #include <rsys/logger.h> 21 #include <star/s2d.h> 22 #include <star/ssp.h> 23 24 #include <limits.h> 25 26 #define NSTEPS 100000 27 28 static const float vertices[] = { 29 1.f, 0.f, 30 0.f, 0.f, 31 0.f, 1.f, 32 1.f, 1.f 33 }; 34 const unsigned nverts = sizeof(vertices)/(2*sizeof(float)); 35 36 const unsigned indices[] = { 37 0, 1, /* Bottom */ 38 1, 2, /* Left */ 39 2, 3, /* Top */ 40 3, 0 /* Right */ 41 }; 42 const unsigned nsegs = sizeof(indices)/(2*sizeof(unsigned)); 43 44 static const double emissivity[] = { 45 1.0, /* Bottom */ 46 1.0, /* Left */ 47 1.0, /* Top */ 48 1.0 /* Right */ 49 }; 50 51 /* Emissivity used to simulate 2 infinite planes */ 52 static const double emissivity_inf_bottom_top[] = { 53 1, /* Bottom */ 54 0, /* Left */ 55 1, /* Top */ 56 0 /* Right */ 57 }; 58 59 static const double specularity[] = { 60 0.0, /* Bottom */ 61 0.0, /* Left */ 62 0.0, /* Top */ 63 0.0 /* Right */ 64 }; 65 66 /* Specularity used to simulate 2 infinite segments */ 67 static const double specularity_inf_bottom_top[] = { 68 0.0, /* Bottom */ 69 1.0, /* Left */ 70 0.0, /* Top */ 71 1.0 /* Right */ 72 }; 73 74 /* Check the estimation of the bottom/top & bottom/medium Gebhart factors */ 75 static void 76 check_bottom_top_medium_gf 77 (struct sgf_scene* scn, 78 struct ssp_rng* rng, 79 const double gf_bottom_top, /* Ref of the bottom/top Gebhart factor */ 80 const double gf_bottom_medium) /* Ref of the bottom/medium Gebhart Factor */ 81 { 82 /* Indices of the top/bottom primitives */ 83 const size_t TOP = 2; 84 const size_t BOTTOM = 0; 85 86 struct sgf_estimator* estimator; 87 struct sgf_status status[2]; 88 89 CHK(sgf_scene_begin_integration(scn) == RES_OK); 90 91 /* Estimate the Gebhart factors for the bottom segment */ 92 CHK(sgf_integrate(scn, BOTTOM, rng, NSTEPS, &estimator) == RES_OK); 93 CHK(sgf_estimator_get_status(estimator, TOP, 0, status + 0) == RES_OK); 94 CHK(sgf_estimator_get_status_medium(estimator, 0, status + 1) == RES_OK); 95 CHK(sgf_estimator_ref_put(estimator) == RES_OK); 96 97 /* Check the Gebhart factor between the bottom and the top plane. */ 98 CHK(eq_eps(status[0].E, gf_bottom_top, status[0].SE) == 1); 99 100 /* Check the Gebhart factor between the bottom plane and the medium */ 101 CHK(eq_eps(status[1].E, gf_bottom_medium, status[1].SE) == 1); 102 103 CHK(sgf_scene_end_integration(scn) == RES_OK); 104 } 105 106 int 107 main(int argc, char** argv) 108 { 109 struct mem_allocator allocator; 110 struct shape_context shape; 111 struct sgf_device* sgf; 112 struct sgf_scene* scn; 113 struct sgf_scene_desc desc = SGF_SCENE_DESC_NULL; 114 struct sgf_status status; 115 struct sgf_estimator* estimator; 116 struct ssp_rng* rng; 117 double ka[4]; 118 size_t iprim, i; 119 (void)argc, (void)argv; 120 121 mem_init_proxy_allocator(&allocator, &mem_default_allocator); 122 123 CHK(ssp_rng_create(&allocator, SSP_RNG_THREEFRY, &rng) == RES_OK); 124 CHK(sgf_device_create(NULL, &allocator, 1, &sgf) == RES_OK); 125 CHK(sgf_scene_create(sgf, &scn) == RES_OK); 126 127 shape.vertices = vertices; 128 shape.nvertices = nverts; 129 shape.indices = indices; 130 shape.nprimitives = nsegs; 131 shape.emissivity = emissivity; 132 shape.specularity = specularity; 133 shape.dim = 2; 134 135 desc.get_position = get_position; 136 desc.get_indices = get_indices; 137 desc.get_emissivity = get_emissivity; 138 desc.get_specularity = get_specularity; 139 desc.get_reflectivity = get_reflectivity; 140 desc.nprims = (unsigned)nsegs; 141 desc.nverts = (unsigned)nverts; 142 desc.nbands = 1; 143 desc.context = &shape; 144 145 CHK(sgf_scene_setup_2d(scn, &desc) == RES_OK); 146 CHK(sgf_integrate(scn, 0, rng, NSTEPS, &estimator) == RES_BAD_OP); 147 148 CHK(sgf_scene_begin_integration(scn) == RES_OK); 149 150 FOR_EACH(iprim, 0, 3) { 151 CHK(sgf_integrate(scn, iprim, rng, NSTEPS, &estimator) == RES_OK); 152 FOR_EACH(i, 0, 3) { 153 CHK(sgf_estimator_get_status(estimator, i, 0, &status) == RES_OK); 154 155 if(i == iprim) { 156 CHK(eq_eps(status.E, 0, status.SE) == 1); 157 } else if(i == iprim + 2 || i == iprim - 2) { /* parallel faces */ 158 CHK(eq_eps(status.E, sqrt(2) - 1, 2*status.SE) == 1); 159 } else { /* Orthogonal faces */ 160 CHK(iprim == i+1 || iprim == i-1); 161 CHK(eq_eps(status.E, 1 - sqrt(2)/2, 2*status.SE) == 1); 162 } 163 } 164 CHK(sgf_estimator_get_status_infinity(estimator, 0, &status) == RES_OK); 165 CHK(eq_eps(status.E, 0, status.SE) == 1); 166 167 CHK(sgf_estimator_ref_put(estimator) == RES_OK); 168 } 169 170 CHK(sgf_integrate(scn, 4, rng, NSTEPS, &estimator) == RES_BAD_ARG); 171 172 CHK(sgf_scene_end_integration(scn) == RES_OK); 173 174 /* 175 * Check medium attenuation with 2 parallel infinite segments. To simulate 176 * this configuration, the top and bottom segments of the square are fully 177 * emissive. The other ones are fully specular and have no emissivity 178 */ 179 180 shape.absorption = ka; 181 shape.emissivity = emissivity_inf_bottom_top; 182 shape.specularity = specularity_inf_bottom_top; 183 desc.get_absorption = get_absorption; 184 185 FOR_EACH(i, 0, 4) ka[i] = 0; 186 CHK(sgf_scene_setup_2d(scn, &desc) == RES_OK); 187 check_bottom_top_medium_gf(scn, rng, 1, 0); 188 189 FOR_EACH(i, 0, 4) ka[i] = 0.1; 190 CHK(sgf_scene_setup_2d(scn, &desc) == RES_OK); 191 check_bottom_top_medium_gf(scn, rng, 0.832583, 0.167417); 192 193 FOR_EACH(i, 0, 4) ka[i] = 1; 194 CHK(sgf_scene_setup_2d(scn, &desc) == RES_OK); 195 check_bottom_top_medium_gf(scn, rng, 0.219384, 0.780616); 196 197 FOR_EACH(i, 0, 4) ka[i] = 10; 198 CHK(sgf_scene_setup_2d(scn, &desc) == RES_OK); 199 check_bottom_top_medium_gf(scn, rng, 7.0975e-6, 0.999992902); 200 201 CHK(ssp_rng_ref_put(rng) == RES_OK); 202 CHK(sgf_device_ref_put(sgf) == RES_OK); 203 CHK(sgf_scene_ref_put(scn) == RES_OK); 204 205 check_memory_allocator(&allocator); 206 mem_shutdown_proxy_allocator(&allocator); 207 CHK(mem_allocated_size() == 0); 208 return 0; 209 } 210