solstice-solver

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

test_ssol_solver7.c (9183B)


      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 #include <rsys/mem_allocator.h>
     21 #include <rsys/image.h>
     22 
     23 #define REFLECTIVITY 0.1
     24 #include "test_ssol_materials.h"
     25 
     26 #define TARGET_SZ 2
     27 #define PLANE_NAME TARGET
     28 #define HALF_X (TARGET_SZ/2)
     29 #define HALF_Y (TARGET_SZ/2)
     30 STATIC_ASSERT((HALF_X * 2 == TARGET_SZ), ONLY_ENVEN_VALUES_FOR_SZ);
     31 #include "test_ssol_rect_geometry.h"
     32 
     33 #define HELIOSTAT_SZ 4
     34 #define POLYGON_NAME HELIOSTAT
     35 #define HALF_X (HELIOSTAT_SZ/2)
     36 #define HALF_Y (HELIOSTAT_SZ/2)
     37 STATIC_ASSERT(HALF_X * 2 == HELIOSTAT_SZ, ONLY_ENVEN_VALUES_FOR_SZ);
     38 #include "test_ssol_rect2D_geometry.h"
     39 
     40 #define HYPERBOL_SZ 24
     41 #define POLYGON_NAME HYPERBOL
     42 #define HALF_X (HYPERBOL_SZ/2)
     43 #define HALF_Y (HYPERBOL_SZ/2)
     44 STATIC_ASSERT(HALF_X * 2 == HYPERBOL_SZ, ONLY_ENVEN_VALUES_FOR_SZ);
     45 #include "test_ssol_rect2D_geometry.h"
     46 
     47 #include <rsys/double33.h>
     48 
     49 #include <star/s3d.h>
     50 #include <star/ssp.h>
     51 
     52 static void
     53 get_wlen(const size_t i, double* wlen, double* data, void* ctx)
     54 {
     55   double wavelengths[3] = { 1, 2, 3 };
     56   double intensities[3] = { 1, 0.8, 1 };
     57   CHK(i < 3);
     58   (void)ctx;
     59   *wlen = wavelengths[i];
     60   *data = intensities[i];
     61 }
     62 
     63 int
     64 main(int argc, char** argv)
     65 {
     66   struct mem_allocator allocator;
     67   struct ssol_device* dev;
     68   struct ssp_rng* rng;
     69   struct ssol_scene* scene;
     70   struct ssol_shape* square;
     71   struct ssol_vertex_data attribs[1] = { SSOL_VERTEX_DATA_NULL__ };
     72   struct ssol_carving carving1 = SSOL_CARVING_NULL;
     73   struct ssol_carving carving2 = SSOL_CARVING_NULL;
     74   struct ssol_material *m_mtl, *bck_mtl, *v_mtl;
     75   struct ssol_mirror_shader m_shader = SSOL_MIRROR_SHADER_NULL;
     76   struct ssol_matte_shader bck_shader = SSOL_MATTE_SHADER_NULL;
     77   struct ssol_object* t_object;
     78   struct ssol_instance* target;
     79   struct ssol_sun* sun;
     80   struct ssol_spectrum* spectrum;
     81   struct ssol_estimator* estimator;
     82   struct ssol_mc_global mc_global;
     83   struct ssol_mc_receiver mc_rcv;
     84   double dir[3], area;
     85   double transform[12]; /* 3x4 column major matrix */
     86   /* primary is a parabol */
     87   struct ssol_quadric quadric1 = SSOL_QUADRIC_DEFAULT;
     88   struct ssol_punched_surface punched1 = SSOL_PUNCHED_SURFACE_NULL;
     89   struct ssol_object* m_object1;
     90   struct ssol_shape* quad_square1;
     91   struct ssol_instance* heliostat;
     92   /* secondary is an hyperbol */
     93   struct ssol_quadric quadric2 = SSOL_QUADRIC_DEFAULT;
     94   struct ssol_punched_surface punched2 = SSOL_PUNCHED_SURFACE_NULL;
     95   struct ssol_object* m_object2;
     96   struct ssol_shape* quad_square2;
     97   struct ssol_instance* secondary;
     98   (void) argc, (void) argv;
     99 
    100 #define H 10
    101 
    102   mem_init_proxy_allocator(&allocator, &mem_default_allocator);
    103 
    104   CHK(ssol_device_create(NULL, &allocator, 1, 0, &dev) == RES_OK);
    105 
    106   CHK(ssp_rng_create(&allocator, SSP_RNG_THREEFRY, &rng) == RES_OK);
    107   CHK(ssol_spectrum_create(dev, &spectrum) == RES_OK);
    108   CHK(ssol_spectrum_setup(spectrum, get_wlen, 3, NULL) == RES_OK);
    109   CHK(ssol_sun_create_directional(dev, &sun) == RES_OK);
    110   CHK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)) == RES_OK);
    111   CHK(ssol_sun_set_spectrum(sun, spectrum) == RES_OK);
    112   CHK(ssol_sun_set_dni(sun, 1000) == RES_OK);
    113   CHK(ssol_scene_create(dev, &scene) == RES_OK);
    114   CHK(ssol_scene_attach_sun(scene, sun) == RES_OK);
    115 
    116   /* create scene content */
    117 
    118   CHK(ssol_shape_create_mesh(dev, &square) == RES_OK);
    119   attribs[0].usage = SSOL_POSITION;
    120   attribs[0].get = get_position;
    121   CHK(ssol_mesh_setup(square, TARGET_NTRIS__, get_ids,
    122     TARGET_NVERTS__, attribs, 1, (void*) &TARGET_DESC__) == RES_OK);
    123 
    124   CHK(ssol_material_create_mirror(dev, &m_mtl) == RES_OK);
    125   m_shader.normal = get_shader_normal;
    126   m_shader.reflectivity = get_shader_reflectivity;
    127   m_shader.roughness = get_shader_roughness;
    128   CHK(ssol_mirror_setup(m_mtl, &m_shader, SSOL_MICROFACET_BECKMANN) == RES_OK);
    129   CHK(ssol_material_create_matte(dev, &bck_mtl) == RES_OK);
    130   bck_shader.normal = get_shader_normal;
    131   bck_shader.reflectivity = get_shader_reflectivity_2;
    132   CHK(ssol_matte_setup(bck_mtl, &bck_shader) == RES_OK);
    133   CHK(ssol_material_create_virtual(dev, &v_mtl) == RES_OK);
    134 
    135   carving1.get = get_polygon_vertices;
    136   carving1.operation = SSOL_AND;
    137   carving1.nb_vertices = HELIOSTAT_NVERTS__;
    138   carving1.context = &HELIOSTAT_EDGES__;
    139 
    140   CHK(ssol_shape_create_punched_surface(dev, &quad_square1) == RES_OK);
    141   quadric1.type = SSOL_QUADRIC_PARABOL;
    142   quadric1.data.parabol.focal = 10 * sqrt(2) * H;
    143   punched1.nb_carvings = 1;
    144   punched1.quadric = &quadric1;
    145   punched1.carvings = &carving1;
    146   CHK(ssol_punched_surface_setup(quad_square1, &punched1) == RES_OK);
    147   CHK(ssol_object_create(dev, &m_object1) == RES_OK);
    148   CHK(ssol_object_add_shaded_shape(m_object1, quad_square1, m_mtl, m_mtl) == RES_OK);
    149   CHK(ssol_object_instantiate(m_object1, &heliostat) == RES_OK);
    150   CHK(ssol_scene_attach_instance(scene, heliostat) == RES_OK);
    151   d33_rotation_yaw(transform, -0.25 * PI);
    152   d3_splat(transform + 9, 0);
    153   transform[9] = 10 * H; /* target the img focal point of the hyperbol */
    154   CHK(ssol_instance_set_transform(heliostat, transform) == RES_OK);
    155 
    156   carving2.get = get_polygon_vertices;
    157   carving2.operation = SSOL_AND;
    158   carving2.nb_vertices = HYPERBOL_NVERTS__;
    159   carving2.context = &HYPERBOL_EDGES__;
    160 
    161   CHK(ssol_shape_create_punched_surface(dev, &quad_square2) == RES_OK);
    162   quadric2.type = SSOL_QUADRIC_HYPERBOL;
    163   quadric2.data.hyperbol.real_focal = 9 * H;
    164   quadric2.data.hyperbol.img_focal = 1 * H;
    165   punched2.nb_carvings = 1;
    166   punched2.quadric = &quadric2;
    167   punched2.carvings = &carving2;
    168   CHK(ssol_punched_surface_setup(quad_square2, &punched2) == RES_OK);
    169   CHK(ssol_object_create(dev, &m_object2) == RES_OK);
    170   CHK(ssol_object_add_shaded_shape(m_object2, quad_square2, m_mtl, v_mtl) == RES_OK);
    171   CHK(ssol_object_instantiate(m_object2, &secondary) == RES_OK);
    172   CHK(ssol_scene_attach_instance(scene, secondary) == RES_OK);
    173   d33_set_identity(transform);
    174   d3_splat(transform + 9, 0);
    175   transform[11] = 9 * H;
    176   CHK(ssol_instance_set_transform(secondary, transform) == RES_OK);
    177   CHK(ssol_instance_sample(secondary, 0) == RES_OK);
    178 
    179   CHK(ssol_object_create(dev, &t_object) == RES_OK);
    180   CHK(ssol_object_add_shaded_shape(t_object, square, bck_mtl, v_mtl) == RES_OK);
    181   CHK(ssol_object_instantiate(t_object, &target) == RES_OK);
    182   CHK(ssol_instance_set_receiver(target, SSOL_FRONT, 0) == RES_OK);
    183   CHK(ssol_instance_sample(target, 0) == RES_OK);
    184   CHK(ssol_scene_attach_instance(scene, target) == RES_OK);
    185   d33_set_identity(transform);
    186   d3_splat(transform + 9, 0);
    187   CHK(ssol_instance_set_transform(target, transform) == RES_OK);
    188 
    189 #define N__ 10000
    190 #define DNI_cos (1000 * cos(0))
    191 #define TOTAL (HELIOSTAT_SZ * HELIOSTAT_SZ * DNI_cos)
    192 #define GET_MC_RCV ssol_estimator_get_mc_receiver
    193 
    194   CHK(ssol_solve(scene, rng, N__, 0, NULL, &estimator) == RES_OK);
    195 
    196   CHK(ssol_estimator_get_mc_global(estimator, &mc_global) == RES_OK);
    197   CHK(ssol_estimator_get_sampled_area(estimator, &area) == RES_OK);
    198   printf("Total = %g\n", area * DNI_cos);
    199   CHK(eq_eps(area * DNI_cos, TOTAL, TOTAL * 1e-4) == 1);
    200   print_global(&mc_global);
    201 
    202   CHK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv) == RES_OK);
    203   printf("Abs(target1) = %g +/- %g\n",
    204     mc_rcv.absorbed_flux.E,
    205     mc_rcv.absorbed_flux.SE);
    206   printf("Ir(target1) = %g +/- %g\n",
    207     mc_rcv.incoming_flux.E, mc_rcv.incoming_flux.SE);
    208   CHK(eq_eps(mc_rcv.incoming_flux.E, TOTAL, TOTAL * 1e-4));
    209   CHK(eq_eps(mc_rcv.absorbed_flux.E,
    210     (1 - REFLECTIVITY) * TOTAL, (1 - REFLECTIVITY) *TOTAL * 1e-4));
    211   CHK(eq_eps(mc_rcv.incoming_flux.SE, 0, 1e-2) == 1);
    212 
    213   /* Free data */
    214   CHK(ssol_instance_ref_put(heliostat) == RES_OK);
    215   CHK(ssol_instance_ref_put(secondary) == RES_OK);
    216   CHK(ssol_instance_ref_put(target) == RES_OK);
    217   CHK(ssol_object_ref_put(m_object1) == RES_OK);
    218   CHK(ssol_object_ref_put(m_object2) == RES_OK);
    219   CHK(ssol_object_ref_put(t_object) == RES_OK);
    220   CHK(ssol_shape_ref_put(square) == RES_OK);
    221   CHK(ssol_shape_ref_put(quad_square1) == RES_OK);
    222   CHK(ssol_shape_ref_put(quad_square2) == RES_OK);
    223   CHK(ssol_material_ref_put(m_mtl) == RES_OK);
    224   CHK(ssol_material_ref_put(bck_mtl) == RES_OK);
    225   CHK(ssol_material_ref_put(v_mtl) == RES_OK);
    226   CHK(ssol_estimator_ref_put(estimator) == RES_OK);
    227   CHK(ssol_device_ref_put(dev) == RES_OK);
    228   CHK(ssol_scene_ref_put(scene) == RES_OK);
    229   CHK(ssp_rng_ref_put(rng) == RES_OK);
    230   CHK(ssol_spectrum_ref_put(spectrum) == RES_OK);
    231   CHK(ssol_sun_ref_put(sun) == RES_OK);
    232 
    233   check_memory_allocator(&allocator);
    234   mem_shutdown_proxy_allocator(&allocator);
    235   CHK(mem_allocated_size() == 0);
    236 
    237   return 0;
    238 }
    239