1.3.49
 
Loading...
Searching...
No Matches
selfTest.cpp
1#include "EnergyBalanceModel.h"
2
3#define DOCTEST_CONFIG_IMPLEMENT
4#include <doctest.h>
5#include "doctest_utils.h"
6
7using namespace helios;
8
9float err_tol = 1e-3;
10
11DOCTEST_TEST_CASE("EnergyBalanceModel Equilibrium Test") {
12 Context context_test;
13 std::vector<uint> UUIDs;
14 UUIDs.push_back(context_test.addPatch(make_vec3(1, 2, 3), make_vec2(3, 2)));
15 UUIDs.push_back(context_test.addTriangle(make_vec3(4, 5, 6), make_vec3(5, 5, 6), make_vec3(5, 6, 6)));
16
17 float Tref = 350;
18 context_test.setPrimitiveData(UUIDs, "radiation_flux_LW", float(2.f * 5.67e-8 * pow(Tref, 4)));
19 context_test.setPrimitiveData(UUIDs, "air_temperature", Tref);
20
21 EnergyBalanceModel energymodeltest(&context_test);
22 energymodeltest.disableMessages();
23 energymodeltest.addRadiationBand("LW");
24 DOCTEST_CHECK_NOTHROW(energymodeltest.run());
25
26 for (int p = 0; p < UUIDs.size(); p++) {
27 float T;
28 DOCTEST_CHECK_NOTHROW(context_test.getPrimitiveData(UUIDs.at(p), "temperature", T));
29 DOCTEST_CHECK(T == doctest::Approx(Tref).epsilon(err_tol));
30 }
31}
32
33DOCTEST_TEST_CASE("EnergyBalanceModel Energy Budget Closure") {
34 Context context_2;
35 std::vector<uint> UUIDs_2;
36 UUIDs_2.push_back(context_2.addPatch(make_vec3(1, 2, 3), make_vec2(3, 2)));
37 UUIDs_2.push_back(context_2.addTriangle(make_vec3(4, 5, 6), make_vec3(5, 5, 6), make_vec3(5, 6, 6)));
38
39 float T = 300;
40 context_2.setPrimitiveData(UUIDs_2, "radiation_flux_LW", float(2.f * 5.67e-8 * pow(T, 4)));
41 context_2.setPrimitiveData(UUIDs_2, "radiation_flux_SW", float(300.f));
42 context_2.setPrimitiveData(UUIDs_2, "air_temperature", T);
43
44 EnergyBalanceModel energymodel_2(&context_2);
45 energymodel_2.disableMessages();
46 energymodel_2.addRadiationBand("LW");
47 energymodel_2.addRadiationBand("SW");
48 DOCTEST_CHECK_NOTHROW(energymodel_2.run());
49
50 for (int p = 0; p < UUIDs_2.size(); p++) {
51 float sensible_flux, latent_flux, R, temperature;
52 float Rin = 0;
53
54 DOCTEST_CHECK_NOTHROW(context_2.getPrimitiveData(UUIDs_2.at(p), "sensible_flux", sensible_flux));
55 DOCTEST_CHECK_NOTHROW(context_2.getPrimitiveData(UUIDs_2.at(p), "latent_flux", latent_flux));
56 DOCTEST_CHECK_NOTHROW(context_2.getPrimitiveData(UUIDs_2.at(p), "radiation_flux_LW", R));
57 Rin += R;
58 DOCTEST_CHECK_NOTHROW(context_2.getPrimitiveData(UUIDs_2.at(p), "radiation_flux_SW", R));
59 Rin += R;
60 DOCTEST_CHECK_NOTHROW(context_2.getPrimitiveData(UUIDs_2.at(p), "temperature", temperature));
61
62 float Rout = 2.f * 5.67e-8 * pow(temperature, 4);
63 float resid = Rin - Rout - sensible_flux - latent_flux;
64 DOCTEST_CHECK(resid == doctest::Approx(0.0f).epsilon(err_tol));
65 }
66}
67
68DOCTEST_TEST_CASE("EnergyBalanceModel Temperature Solution Check 1") {
69 Context context_3;
70 uint UUID_3 = context_3.addPatch(make_vec3(1, 2, 3), make_vec2(3, 2));
71
72 float T = 312.f;
73 context_3.setPrimitiveData(UUID_3, "radiation_flux_LW", float(5.67e-8 * pow(T, 4)));
74 context_3.setPrimitiveData(UUID_3, "radiation_flux_SW", float(350.f));
75 context_3.setPrimitiveData(UUID_3, "wind_speed", float(1.244f));
76 context_3.setPrimitiveData(UUID_3, "moisture_conductance", float(0.05f));
77 context_3.setPrimitiveData(UUID_3, "air_humidity", float(0.4f));
78 context_3.setPrimitiveData(UUID_3, "air_pressure", float(956789));
79 context_3.setPrimitiveData(UUID_3, "other_surface_flux", float(150.f));
80 context_3.setPrimitiveData(UUID_3, "air_temperature", T);
81 context_3.setPrimitiveData(UUID_3, "twosided_flag", uint(0));
82
83 EnergyBalanceModel energymodel_3(&context_3);
84 energymodel_3.disableMessages();
85 energymodel_3.addRadiationBand("LW");
86 energymodel_3.addRadiationBand("SW");
87 DOCTEST_CHECK_NOTHROW(energymodel_3.run());
88
89 float sensible_flux, latent_flux, temperature;
90 float sensible_flux_exact = 48.7017;
91 float latent_flux_exact = 21.6094;
92 float temperature_exact = 329.307;
93
94 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "sensible_flux", sensible_flux));
95 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "latent_flux", latent_flux));
96 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "temperature", temperature));
97
98 DOCTEST_CHECK(sensible_flux == doctest::Approx(sensible_flux_exact).epsilon(err_tol));
99 DOCTEST_CHECK(latent_flux == doctest::Approx(latent_flux_exact).epsilon(err_tol));
100 DOCTEST_CHECK(temperature == doctest::Approx(temperature_exact).epsilon(err_tol));
101}
102
103DOCTEST_TEST_CASE("EnergyBalanceModel Temperature Solution Check 2 - Object Length") {
104 Context context_3;
105 uint UUID_3 = context_3.addPatch(make_vec3(1, 2, 3), make_vec2(3, 2));
106
107 float T = 312.f;
108 context_3.setPrimitiveData(UUID_3, "radiation_flux_LW", float(5.67e-8 * pow(T, 4)));
109 context_3.setPrimitiveData(UUID_3, "radiation_flux_SW", float(350.f));
110 context_3.setPrimitiveData(UUID_3, "wind_speed", float(1.244f));
111 context_3.setPrimitiveData(UUID_3, "moisture_conductance", float(0.05f));
112 context_3.setPrimitiveData(UUID_3, "air_humidity", float(0.4f));
113 context_3.setPrimitiveData(UUID_3, "air_pressure", float(956789));
114 context_3.setPrimitiveData(UUID_3, "other_surface_flux", float(150.f));
115 context_3.setPrimitiveData(UUID_3, "air_temperature", T);
116 context_3.setPrimitiveData(UUID_3, "twosided_flag", uint(0));
117
118 EnergyBalanceModel energymodel_3(&context_3);
119 energymodel_3.disableMessages();
120 energymodel_3.addRadiationBand("LW");
121 energymodel_3.addRadiationBand("SW");
122
123 // Use object length instead of sqrt(area)
124 context_3.setPrimitiveData(UUID_3, "object_length", float(0.374f));
125 DOCTEST_CHECK_NOTHROW(energymodel_3.run());
126
127 float sensible_flux, latent_flux, temperature;
128 float sensible_flux_exact = 89.2024;
129 float latent_flux_exact = 20.0723;
130 float temperature_exact = 324.386;
131
132 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "sensible_flux", sensible_flux));
133 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "latent_flux", latent_flux));
134 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "temperature", temperature));
135
136 DOCTEST_CHECK(sensible_flux == doctest::Approx(sensible_flux_exact).epsilon(err_tol));
137 DOCTEST_CHECK(latent_flux == doctest::Approx(latent_flux_exact).epsilon(err_tol));
138 DOCTEST_CHECK(temperature == doctest::Approx(temperature_exact).epsilon(err_tol));
139}
140
141DOCTEST_TEST_CASE("EnergyBalanceModel Temperature Solution Check 3 - Manual Boundary Layer Conductance") {
142 Context context_3;
143 uint UUID_3 = context_3.addPatch(make_vec3(1, 2, 3), make_vec2(3, 2));
144
145 float T = 312.f;
146 context_3.setPrimitiveData(UUID_3, "radiation_flux_LW", float(5.67e-8 * pow(T, 4)));
147 context_3.setPrimitiveData(UUID_3, "radiation_flux_SW", float(350.f));
148 context_3.setPrimitiveData(UUID_3, "wind_speed", float(1.244f));
149 context_3.setPrimitiveData(UUID_3, "moisture_conductance", float(0.05f));
150 context_3.setPrimitiveData(UUID_3, "air_humidity", float(0.4f));
151 context_3.setPrimitiveData(UUID_3, "air_pressure", float(956789));
152 context_3.setPrimitiveData(UUID_3, "other_surface_flux", float(150.f));
153 context_3.setPrimitiveData(UUID_3, "air_temperature", T);
154 context_3.setPrimitiveData(UUID_3, "twosided_flag", uint(0));
155
156 EnergyBalanceModel energymodel_3(&context_3);
157 energymodel_3.disableMessages();
158 energymodel_3.addRadiationBand("LW");
159 energymodel_3.addRadiationBand("SW");
160
161 // Manually set boundary-layer conductance
162 context_3.setPrimitiveData(UUID_3, "boundarylayer_conductance", float(0.134f));
163 DOCTEST_CHECK_NOTHROW(energymodel_3.run());
164
165 float sensible_flux, latent_flux, temperature;
166 float sensible_flux_exact = 61.5411f;
167 float latent_flux_exact = 21.6718f;
168 float temperature_exact = 327.701f;
169
170 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "sensible_flux", sensible_flux));
171 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "latent_flux", latent_flux));
172 DOCTEST_CHECK_NOTHROW(context_3.getPrimitiveData(UUID_3, "temperature", temperature));
173
174 DOCTEST_CHECK(sensible_flux == doctest::Approx(sensible_flux_exact).epsilon(err_tol));
175 DOCTEST_CHECK(latent_flux == doctest::Approx(latent_flux_exact).epsilon(err_tol));
176 DOCTEST_CHECK(temperature == doctest::Approx(temperature_exact).epsilon(err_tol));
177}
178
179DOCTEST_TEST_CASE("EnergyBalanceModel Optional Primitive Data Output Check") {
180 Context context_4;
181 uint UUID_4 = context_4.addPatch(make_vec3(1, 2, 3), make_vec2(3, 2));
182
183 EnergyBalanceModel energymodel_4(&context_4);
184 energymodel_4.disableMessages();
185 energymodel_4.addRadiationBand("LW");
186 DOCTEST_CHECK_NOTHROW(energymodel_4.optionalOutputPrimitiveData("boundarylayer_conductance_out"));
187 DOCTEST_CHECK_NOTHROW(energymodel_4.optionalOutputPrimitiveData("vapor_pressure_deficit"));
188
189 context_4.setPrimitiveData(UUID_4, "radiation_flux_LW", 0.f);
190 DOCTEST_CHECK_NOTHROW(energymodel_4.run());
191
192 DOCTEST_CHECK(context_4.doesPrimitiveDataExist(UUID_4, "vapor_pressure_deficit"));
193 DOCTEST_CHECK(context_4.doesPrimitiveDataExist(UUID_4, "boundarylayer_conductance_out"));
194}
195
196DOCTEST_TEST_CASE("EnergyBalanceModel Dynamic Model Check") {
197 Context context_5;
198 float dt_5 = 1.f, T_5 = 3600, To_5 = 300.f, cp_5 = 2000;
199 float Rlow = 50.f, Rhigh = 500.f;
200
201 uint UUID_5 = context_5.addPatch(make_vec3(0, 0, 0), make_vec2(1, 1));
202 context_5.setPrimitiveData(UUID_5, "radiation_flux_SW", Rlow);
203 context_5.setPrimitiveData(UUID_5, "temperature", To_5);
204 context_5.setPrimitiveData(UUID_5, "heat_capacity", cp_5);
205 context_5.setPrimitiveData(UUID_5, "twosided_flag", uint(0));
206 context_5.setPrimitiveData(UUID_5, "emissivity_SW", 0.f);
207
208 EnergyBalanceModel energybalance_5(&context_5);
209 energybalance_5.disableMessages();
210 energybalance_5.addRadiationBand("SW");
211 DOCTEST_CHECK_NOTHROW(energybalance_5.optionalOutputPrimitiveData("boundarylayer_conductance_out"));
212
213 std::vector<float> temperature_dyn;
214 int N = round(T_5 / dt_5);
215 for (int t = 0; t < N; t++) {
216 if (t > 0.5f * N) {
217 context_5.setPrimitiveData(UUID_5, "radiation_flux_SW", Rhigh);
218 }
219 DOCTEST_CHECK_NOTHROW(energybalance_5.run(dt_5));
220
221 float temp;
222 DOCTEST_CHECK_NOTHROW(context_5.getPrimitiveData(UUID_5, "temperature", temp));
223 temperature_dyn.push_back(temp);
224 }
225
226 float gH_5;
227 DOCTEST_CHECK_NOTHROW(context_5.getPrimitiveData(UUID_5, "boundarylayer_conductance_out", gH_5));
228 float tau_5 = cp_5 / gH_5 / 29.25f;
229
230 context_5.setPrimitiveData(UUID_5, "radiation_flux_SW", Rlow);
231 DOCTEST_CHECK_NOTHROW(energybalance_5.run());
232 float Tlow;
233 DOCTEST_CHECK_NOTHROW(context_5.getPrimitiveData(UUID_5, "temperature", Tlow));
234
235 context_5.setPrimitiveData(UUID_5, "radiation_flux_SW", Rhigh);
236 DOCTEST_CHECK_NOTHROW(energybalance_5.run());
237 float Thigh;
238 DOCTEST_CHECK_NOTHROW(context_5.getPrimitiveData(UUID_5, "temperature", Thigh));
239
240 float err = 0;
241 for (int t = round(0.5f * N); t < N; t++) {
242 float time = dt_5 * (t - round(0.5f * N));
243 float temperature_ref = Tlow + (Thigh - Tlow) * (1.f - exp(-time / tau_5));
244 err += pow(temperature_ref - temperature_dyn.at(t), 2);
245 }
246
247 err = sqrt(err / float(N));
248 DOCTEST_CHECK(err == doctest::Approx(0.0f).epsilon(0.2f));
249}
250
251DOCTEST_TEST_CASE("EnergyBalanceModel Enable/Disable Messages") {
252 Context context_enable;
253 EnergyBalanceModel testModel(&context_enable);
254
255 DOCTEST_CHECK_NOTHROW(testModel.enableMessages());
256 DOCTEST_CHECK_NOTHROW(testModel.disableMessages());
257}
258
259DOCTEST_TEST_CASE("EnergyBalanceModel Radiation Band Management") {
260 Context context_radiation;
261 EnergyBalanceModel testModel(&context_radiation);
262
263 DOCTEST_CHECK_NOTHROW(testModel.addRadiationBand("LW"));
264 DOCTEST_CHECK_NOTHROW(testModel.addRadiationBand("LW")); // Should not duplicate
265 DOCTEST_CHECK_NOTHROW(testModel.addRadiationBand("SW"));
266}
267
268DOCTEST_TEST_CASE("EnergyBalanceModel Optional Output Primitive Data") {
269 Context context_output;
270 EnergyBalanceModel testModel(&context_output);
271
272 DOCTEST_CHECK_NOTHROW(testModel.optionalOutputPrimitiveData("boundarylayer_conductance_out"));
273
274 capture_cerr cerr_buffer;
275 DOCTEST_CHECK_NOTHROW(testModel.optionalOutputPrimitiveData("invalid_label")); // Should print warning
276 DOCTEST_CHECK(cerr_buffer.has_output());
277}
278
279DOCTEST_TEST_CASE("EnergyBalanceModel Print Default Value Report") {
280 Context context_print;
281 EnergyBalanceModel testModel(&context_print);
282
283 std::stringstream buffer;
284 std::streambuf *old = std::cout.rdbuf(buffer.rdbuf()); // Redirect std::cout
285
286 DOCTEST_CHECK_NOTHROW(testModel.printDefaultValueReport());
287
288 std::cout.rdbuf(old); // Restore std::cout
289
290 std::string output = buffer.str();
291 std::vector<std::string> required_keywords = {"surface temperature", "air pressure", "air temperature", "air humidity", "boundary-layer conductance", "moisture conductance", "surface humidity", "two-sided flag", "evaporating faces"};
292
293 for (const auto &keyword: required_keywords) {
294 DOCTEST_CHECK(output.find(keyword) != std::string::npos);
295 }
296}
297
298DOCTEST_TEST_CASE("EnergyBalanceModel Print Default Value Report with UUIDs") {
299 Context context_print;
300 EnergyBalanceModel testModel(&context_print);
301 std::vector<uint> testUUIDs;
302 testUUIDs.push_back(context_print.addPatch(make_vec3(1, 2, 3), make_vec2(3, 2)));
303
304 std::stringstream buffer;
305 std::streambuf *old = std::cout.rdbuf(buffer.rdbuf()); // Redirect std::cout
306
307 DOCTEST_CHECK_NOTHROW(testModel.printDefaultValueReport(testUUIDs));
308
309 std::cout.rdbuf(old); // Restore std::cout
310
311 std::string output = buffer.str();
312 std::vector<std::string> required_keywords = {"surface temperature", "air pressure", "air temperature", "air humidity", "boundary-layer conductance", "moisture conductance", "surface humidity", "two-sided flag", "evaporating faces"};
313
314 for (const auto &keyword: required_keywords) {
315 DOCTEST_CHECK(output.find(keyword) != std::string::npos);
316 }
317}
318
319DOCTEST_TEST_CASE("EnergyBalanceModel Additional Dynamic Model Check") {
320 Context context_dyn;
321 float dt = 1.f, Tfinal = 3600, To = 300.f, cp = 2000;
322 float Rlow = 50.f, Rhigh = 500.f;
323
324 uint UUID_dyn = context_dyn.addPatch(make_vec3(0, 0, 0), make_vec2(1, 1));
325 context_dyn.setPrimitiveData(UUID_dyn, "radiation_flux_SW", Rlow);
326 context_dyn.setPrimitiveData(UUID_dyn, "temperature", To);
327 context_dyn.setPrimitiveData(UUID_dyn, "heat_capacity", cp);
328 context_dyn.setPrimitiveData(UUID_dyn, "twosided_flag", uint(0));
329 context_dyn.setPrimitiveData(UUID_dyn, "emissivity_SW", 0.f);
330
331 EnergyBalanceModel energybalance_dyn(&context_dyn);
332 energybalance_dyn.disableMessages();
333 energybalance_dyn.addRadiationBand("SW");
334
335 std::vector<float> temperature_dyn;
336 int N = round(Tfinal / dt);
337 for (int t = 0; t < N; t++) {
338 if (t > 0.5f * N) {
339 context_dyn.setPrimitiveData(UUID_dyn, "radiation_flux_SW", Rhigh);
340 }
341 DOCTEST_CHECK_NOTHROW(energybalance_dyn.run(dt));
342
343 float temp;
344 DOCTEST_CHECK_NOTHROW(context_dyn.getPrimitiveData(UUID_dyn, "temperature", temp));
345 temperature_dyn.push_back(temp);
346 }
347
348 DOCTEST_CHECK(!temperature_dyn.empty());
349}
350
351int EnergyBalanceModel::selfTest(int argc, char **argv) {
352 return helios::runDoctestWithValidation(argc, argv);
353}