star-line

Structure for accelerating line importance sampling
git clone git://git.meso-star.fr/star-line.git
Log | Files | Refs | README | LICENSE

commit 1c0ac0cce49ddb62c08888af958782db8987f5d8
parent ff95468f9c6fa72f95bdbf598df3dd7a73a7a0e1
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Fri, 22 Apr 2022 15:55:41 +0200

Test the sln_tree_create function

Diffstat:
Mcmake/CMakeLists.txt | 1+
Asrc/test_sln_tree.c | 325+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 files changed, 326 insertions(+), 0 deletions(-)

diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -89,6 +89,7 @@ if(NOT NO_TEST) endfunction() new_test(test_sln_device) + new_test(test_sln_tree StarHITRAN) endif() ################################################################################ diff --git a/src/test_sln_tree.c b/src/test_sln_tree.c @@ -0,0 +1,325 @@ +/* Copyright (C) 2022 CNRS - LMD + * Copyright (C) 2022 |Meso|Star> (contact@meso-star.com) + * Copyright (C) 2022 Université Paul Sabatier - IRIT + * Copyright (C) 2022 Université Paul Sabatier - Laplace + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. */ + +#include "sln.h" +#include <star/shtr.h> +#include <rsys/math.h> +#include <rsys/mem_allocator.h> + +/******************************************************************************* + * Input data + ******************************************************************************/ +static const struct shtr_isotope g_H2O_isotopes[] = { + {9.97317E-01, 1.7458E+02, 18.010565, 0, 1, 161}, + {1.99983E-03, 1.7605E+02, 20.014811, 0, 1, 181}, + {3.71884E-04, 1.0521E+03, 19.014780, 0, 6, 171}, + {3.10693E-04, 8.6474E+02, 19.016740, 0, 6, 162}, + {6.23003E-07, 8.7557E+02, 21.020985, 0, 6, 182}, + {1.15853E-07, 5.2268E+03, 20.020956, 0, 36, 172}, + {2.41974E-08, 1.0278E+03, 20.022915, 0, 1, 262} +}; +static const struct shtr_molecule g_H2O = { + "H2O", sizeof(g_H2O_isotopes)/sizeof(struct shtr_isotope), g_H2O_isotopes, 1, +}; + +static const struct shtr_isotope g_CO2_isotopes[] = { + {9.84204E-01, 2.8609E+02, 43.989830, 1, 1, 626}, + {1.10574E-02, 5.7664E+02, 44.993185, 1, 2, 636}, + {3.94707E-03, 6.0781E+02, 45.994076, 1, 1, 628}, + {7.33989E-04, 3.5426E+03, 44.994045, 1, 6, 627}, + {4.43446E-05, 1.2255E+03, 46.997431, 1, 2, 638}, + {8.24623E-06, 7.1413E+03, 45.997400, 1, 12, 637}, + {3.95734E-06, 3.2342E+02, 47.998320, 1, 1, 828}, + {1.47180E-06, 3.7666E+03, 46.998291, 1, 6, 827}, + {1.36847E-07, 1.0972E+04, 45.998262, 1, 1, 727}, + {4.44600E-08, 6.5224E+02, 49.001675, 1, 2, 838}, + {1.65354E-08, 7.5950E+03, 48.001646, 1, 12, 837}, + {1.53745E-09, 2.2120E+04, 47.001618, 1, 2, 737} +}; +static const struct shtr_molecule g_CO2 = { + "CO2", sizeof(g_CO2_isotopes)/sizeof(struct shtr_isotope), g_CO2_isotopes, 2, +}; + +static const struct shtr_isotope g_O3_isotopes[] = { + {9.92901E-01, 3.4750E+03, 47.984745, 2, 1, 666}, + {3.98194E-03, 7.5846E+03, 49.988991, 2, 1, 668}, + {1.99097E-03, 3.7030E+03, 49.988991, 2, 1, 686}, + {7.40475E-04, 4.4044E+04, 48.988960, 2, 6, 667}, + {3.70237E-04, 2.1742E+04, 48.988960, 2, 6, 676} +}; +static const struct shtr_molecule g_O3 = { + "O3", sizeof(g_O3_isotopes)/sizeof(struct shtr_isotope), g_O3_isotopes, 3 +}; + +static const struct shtr_line g_lines[] = { + {0.00156, 2.685e-36, 0.0695, 0.428, 399.7263, 0.79, 0.000240, 1, 5}, + {0.03406, 5.898e-37, 0.0883, 0.410, 1711.6192, 0.79, 0.000000, 1, 5}, + {0.03445, 5.263e-33, 0.0504, 0.329, 1570.0616, 0.79, 0.001940, 1, 3}, + {0.03628, 9.296e-29, 0.0704, 0.095, 522.5576, 0.81, 0.000000, 3, 0}, + {0.03671, 2.555e-37, 0.0389, 0.304, 2337.5190, 0.79, 0.002660, 1, 4}, + {0.08403, 8.190e-33, 0.0753, 0.394, 568.0017, 0.79, 0.002180, 1, 4}, + {0.08653, 5.376e-29, 0.0691, 0.083, 1061.6864, 0.76, 0.000000, 3, 0}, + {0.09642, 6.675e-37, 0.0570, 0.351, 2516.3150, 0.79, 0.000000, 1, 4}, + {0.16772, 2.456e-28, 0.0823, 0.105, 751.2922, 0.77, 0.000000, 3, 0}, + {0.18582, 5.338e-32, 0.0925, 0.428, 1717.3356, 0.79, 0.004100, 1, 3}, + {0.19368, 2.553e-32, 0.0901, 0.428, 293.8010, 0.79, 0.001260, 1, 5}, + {0.19688, 1.447e-31, 0.0901, 0.428, 292.3093, 0.79, 0.001260, 1, 4}, + {0.19757, 6.063e-29, 0.0681, 0.078, 1401.6146, 0.76, 0.000000, 3, 0}, + {0.21281, 8.238e-31, 0.0780, 0.103, 127.4733, 0.78, 0.000000, 3, 4}, + {0.21282, 5.999e-31, 0.0780, 0.103, 127.4733, 0.78, 0.000000, 3, 4}, + {0.21283, 7.737e-31, 0.0780, 0.103, 127.4733, 0.78, 0.000000, 3, 4}, + {0.21283, 6.394e-31, 0.0780, 0.103, 127.4733, 0.78, 0.000000, 3, 4}, + {0.21284, 7.260e-31, 0.0780, 0.103, 127.4732, 0.78, 0.000000, 3, 4}, + {0.21284, 6.813e-31, 0.0780, 0.103, 127.4733, 0.78, 0.000000, 3, 4}, + {0.21728, 5.928e-29, 0.0823, 0.105, 1153.4099, 0.77, 0.000000, 3, 0}, + {0.25818, 3.110e-32, 0.0539, 0.311, 1958.1248, 0.79, 0.008600, 1, 3}, + {0.26618, 1.468e-32, 0.0801, 0.378, 2118.9452, 0.79, 0.003100, 1, 3}, + {0.27091, 1.362e-35, 0.0539, 0.311, 1949.2032, 0.79, 0.008600, 1, 5}, + {0.28910, 2.058e-30, 0.0878, 0.106, 7.8611, 0.76, 0.000000, 3, 3}, + {0.29412, 8.666e-33, 0.0795, 0.378, 679.8760, 0.79, 0.004170, 1, 5}, + {0.29477, 1.457e-30, 0.0660, 0.340, 1238.7943, 0.79, 0.007260, 1, 3}, + {0.29673, 1.609e-30, 0.0773, 0.103, 138.7916, 0.78, 0.000000, 3, 3}, + {0.29673, 1.195e-30, 0.0773, 0.103, 138.7916, 0.78, 0.000000, 3, 3}, + {0.29676, 1.268e-30, 0.0773, 0.103, 138.7916, 0.78, 0.000000, 3, 3}, + {1.46280, 7.205e-26, 0.0696, 0.090, 702.3167, 0.82, 0.000000, 3, 0}, + {1.46899, 7.847e-28, 0.0704, 0.095, 1557.6290, 0.81, 0.000000, 3, 0}, + {1.47237, 1.700e-34, 0.0633, 0.297, 2389.2994, 0.79, 0.013980, 1, 4}, + {1.47273, 1.131e-32, 0.0903, 0.123, 0.7753, 0.69, 0.000814, 2, 2}, + {1.47863, 5.318e-33, 0.0925, 0.428, 2877.6872, 0.79, 0.006200, 1, 3}, + {1.48860, 1.058e-27, 0.0760, 0.102, 163.7760, 0.78, 0.000000, 3, 1}, + {1.49043, 4.787e-31, 0.0773, 0.103, 139.7879, 0.78, 0.000000, 3, 4}, + {1.49044, 3.139e-29, 0.0773, 0.103, 139.7879, 0.78, 0.000000, 3, 4}, + {1.49044, 3.332e-29, 0.0773, 0.103, 139.7880, 0.78, 0.000000, 3, 4}, + {1.49045, 2.785e-29, 0.0773, 0.103, 139.7880, 0.78, 0.000000, 3, 4}, + {1.49045, 2.625e-29, 0.0773, 0.103, 139.7880, 0.78, 0.000000, 3, 4}, + {1.49057, 5.941e-27, 0.0909, 0.106, 703.4398, 0.76, 0.000000, 3, 0}, + {1.49110, 1.537e-29, 0.0714, 0.098, 381.8992, 0.80, 0.000000, 3, 3}, + {1.49111, 1.828e-29, 0.0714, 0.098, 381.8992, 0.80, 0.000000, 3, 3}, + {1.49111, 1.648e-29, 0.0714, 0.098, 381.8992, 0.80, 0.000000, 3, 3}, + {1.49111, 1.766e-29, 0.0714, 0.098, 381.8991, 0.80, 0.000000, 3, 3}, + {1.49111, 1.706e-29, 0.0714, 0.098, 381.8991, 0.80, 0.000000, 3, 3}, + {1.49674, 1.208e-25, 0.0719, 0.102, 190.2125, 0.79, 0.000000, 3, 0}, + {1.50216, 4.971e-28, 0.0696, 0.090, 1734.9796, 0.82, 0.000000, 3, 0}, + {1.51178, 3.047e-27, 0.0766, 0.103, 1195.5580, 0.78, 0.000000, 3, 0}, + {1.51399, 1.986e-27, 0.0765, 0.104, 1129.0675, 0.77, 0.000000, 3, 0}, + {1.51442, 2.601e-34, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51443, 3.981e-34, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51443, 6.637e-35, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51444, 3.902e-34, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51445, 1.394e-33, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51446, 7.166e-34, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51446, 2.787e-34, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51446, 5.096e-34, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51448, 1.672e-34, 0.0903, 0.123, 0.7572, 0.69, 0.000814, 2, 3}, + {1.51597, 4.443e-30, 0.0704, 0.097, 463.6953, 0.81, 0.000000, 3, 3}, + {1.51600, 5.172e-30, 0.0704, 0.097, 463.6952, 0.81, 0.000000, 3, 3}, + {1.51603, 4.580e-30, 0.0704, 0.097, 463.6952, 0.81, 0.000000, 3, 3}, + {1.51605, 5.019e-30, 0.0704, 0.097, 463.6951, 0.81, 0.000000, 3, 3} +}; +const size_t g_nlines = sizeof(g_lines) / sizeof(struct shtr_line); + +static void +write_shtr_isotope(FILE* fp, const struct shtr_isotope* isotope) +{ + CHK(fp && isotope); + fprintf(fp, " %d %.5E %.4E %d %.6f\n", + isotope->id, + isotope->abundance, + isotope->Q296K, + isotope->gj, + isotope->molar_mass); +} + +static void +write_shtr_molecule(FILE* fp, const struct shtr_molecule* molecule) +{ + size_t i; + CHK(fp && molecule); + + fprintf(fp, " %s (%d)\n", molecule->name, molecule->id); + FOR_EACH(i, 0, molecule->nisotopes) { + write_shtr_isotope(fp, molecule->isotopes+i); + } +} + +static void +write_shtr_lines + (FILE* fp, + const struct shtr_line* lines, + const size_t nlines) +{ + size_t i; + + CHK(fp && (!nlines || lines)); + FOR_EACH(i, 0, nlines) { + fprintf(fp, + "%2d%1d%12.6f%10.3e 0.000E-00.%04d%5.3f%10.4f%4.2f%8.6f" + /* Dummy remaining data */ + " 0 0 0" /* Global upper quanta */ + " 0 0 0" /* Global upper quanta */ + " 5 5 0 " /* Local upper quanta */ + " 5 5 1 " /* Local lower quanta */ + "562220" /* Error indices */ + "5041 7833348" /* References */ + " " /* Line mixing flag */ + " 66.0" /* g' */ + " 66.0" /* g'' */ + "\n", + lines[i].molecule_id, + lines[i].isotope_id_local == 9 ? 0 : lines[i].isotope_id_local+1, + lines[i].wavenumber, + lines[i].intensity, + (int)(lines[i].gamma_air*10000+0.5/*round*/), + lines[i].gamma_self, + lines[i].lower_state_energy, + lines[i].n_air, + lines[i].delta_air); + } +} + +/******************************************************************************* + * Helper function + ******************************************************************************/ +static void +test_tree + (struct sln_device* sln, + struct shtr_isotope_metadata* metadata, + struct shtr_lines_list* lines_list) +{ + struct sln_tree_create_args tree_args = SLN_TREE_CREATE_ARGS_DEFAULT; + struct sln_tree* tree = NULL; + + tree_args.metadata = metadata; + tree_args.lines = lines_list; + tree_args.molecules[0].nisotopes = 0; /* Handle all isotopes */ + tree_args.molecules[0].concentration = 1.0/3.0; + tree_args.molecules[0].cutoff = 25; + tree_args.molecules[0].id = 1; /* H2O */ + tree_args.molecules[1].nisotopes = 0; /* Handle all isotopes */ + tree_args.molecules[1].concentration = 1.0/3.0; + tree_args.molecules[1].cutoff = 50; + tree_args.molecules[1].id = 2; /* CO2 */ + tree_args.molecules[2].nisotopes = 0; /* Handle all isotopes */ + tree_args.molecules[2].concentration = 1.0/3.0; + tree_args.molecules[2].cutoff = 25; + tree_args.molecules[2].id = 3; /* CO2 */ + tree_args.nmolecules = 3; + tree_args.wavenumber_range[0] = 0; + tree_args.wavenumber_range[1] = INF; + tree_args.pressure = 1; + tree_args.temperature = 296; + tree_args.max_nlines_per_leaf = 1; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_OK); + + CHK(sln_tree_ref_get(NULL) == RES_BAD_ARG); + CHK(sln_tree_ref_get(tree) == RES_OK); + CHK(sln_tree_ref_put(NULL) == RES_BAD_ARG); + CHK(sln_tree_ref_put(tree) == RES_OK); + CHK(sln_tree_ref_put(tree) == RES_OK); + + tree_args.metadata = NULL; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.metadata = metadata; + tree_args.lines = NULL; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.metadata = metadata; + tree_args.lines = NULL; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.lines = lines_list; + tree_args.molecules[0].concentration = 1; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.molecules[0].concentration = -1; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.molecules[0].concentration = 1.0/3.0; + tree_args.molecules[0].cutoff = 0; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.molecules[0].cutoff = 25; + tree_args.pressure = -1; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.pressure = 1; + tree_args.temperature = -1; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.temperature = 296; + tree_args.max_nlines_per_leaf = 0; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_BAD_ARG); + + tree_args.max_nlines_per_leaf = 16; + CHK(sln_tree_create(sln, &tree_args, &tree) == RES_OK); + CHK(sln_tree_ref_put(tree) == RES_OK); +} + +/******************************************************************************* + * Test function + ******************************************************************************/ +int +main(int argc, char** argv) +{ + struct sln_device_create_args dev_args = SLN_DEVICE_CREATE_ARGS_DEFAULT; + struct sln_device* sln = NULL; + + struct shtr_create_args shtr_args = SHTR_CREATE_ARGS_DEFAULT; + struct shtr* shtr = NULL; + struct shtr_lines_list* lines_list = NULL; + struct shtr_isotope_metadata* metadata = NULL; + + FILE* fp_lines = NULL; + FILE* fp_mdata = NULL; + (void)argc, (void)argv; + + /* Generate the file of the isotope metadata */ + CHK(fp_mdata = tmpfile()); + fprintf(fp_mdata, "Molecule # Iso Abundance Q(296K) gj Molar Mass(g)\n"); + write_shtr_molecule(fp_mdata, &g_H2O); + write_shtr_molecule(fp_mdata, &g_CO2); + write_shtr_molecule(fp_mdata, &g_O3); + rewind(fp_mdata); + + /* Generate the file of lines */ + CHK(fp_lines = tmpfile()); + write_shtr_lines(fp_lines, g_lines, g_nlines); + rewind(fp_lines); + + /* Load the isotope metadata and the lines */ + shtr_args.verbose = 1; + CHK(shtr_create(&shtr_args, &shtr) == RES_OK); + CHK(shtr_isotope_metadata_load_stream(shtr, fp_mdata, NULL, &metadata) == RES_OK); + CHK(shtr_lines_list_load_stream(shtr, fp_lines, NULL, &lines_list) == RES_OK); + + CHK(fclose(fp_lines) == 0); + CHK(fclose(fp_mdata) == 0); + + dev_args.verbose = 1; + CHK(sln_device_create(&dev_args, &sln) == RES_OK); + + test_tree(sln, metadata, lines_list); + + CHK(sln_device_ref_put(sln) == RES_OK); + CHK(shtr_ref_put(shtr) == RES_OK); + CHK(shtr_lines_list_ref_put(lines_list) == RES_OK); + CHK(shtr_isotope_metadata_ref_put(metadata) == RES_OK); + CHK(mem_allocated_size() == 0); + return 0; +}