GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: dumux/dumux/material/components/tabulatedcomponent.hh
Date: 2025-04-12 19:19:20
Exec Total Coverage
Lines: 299 331 90.3%
Functions: 159 163 97.5%
Branches: 407 562 72.4%

Line Branch Exec Source
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 //
4 // SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5 // SPDX-License-Identifier: GPL-3.0-or-later
6 //
7 /*!
8 * \file
9 * \ingroup Components
10 * \brief Tabulates all thermodynamic properties of a given
11 * untabulated chemical species.
12 *
13 * At the moment, this class can only handle the sub-critical fluids
14 * since it tabulates along the vapor pressure curve.
15 */
16 #ifndef DUMUX_TABULATED_COMPONENT_HH
17 #define DUMUX_TABULATED_COMPONENT_HH
18
19 #include <cmath>
20 #include <limits>
21 #include <cassert>
22 #include <vector>
23 #include <iostream>
24 #include <iomanip>
25 #include <algorithm>
26
27 #include <dune/common/std/type_traits.hh>
28
29 #include <dumux/common/typetraits/typetraits.hh>
30 #include <dumux/common/exceptions.hh>
31 #include <dumux/parallel/multithreading.hh>
32 #include <dumux/parallel/parallel_for.hh>
33 #include <dumux/material/components/componenttraits.hh>
34
35 namespace Dumux::Components {
36 // forward declaration
37 template<class RawComponent, bool useVaporPressure>
38 class TabulatedComponent;
39 } // end namespace Dumux::Components
40
41 namespace Dumux {
42
43 //! component traits for tabulated component
44 template<class RawComponent, bool useVaporPressure>
45 struct ComponentTraits<Components::TabulatedComponent<RawComponent, useVaporPressure>>
46 {
47 using Scalar = typename RawComponent::Scalar;
48
49 //! if the component implements a solid state
50 static constexpr bool hasSolidState = std::is_base_of<Components::Solid<Scalar, RawComponent>, RawComponent>::value;
51
52 //! if the component implements a liquid state
53 static constexpr bool hasLiquidState = std::is_base_of<Components::Liquid<Scalar, RawComponent>, RawComponent>::value;
54
55 //! if the component implements a gaseous state
56 static constexpr bool hasGasState = std::is_base_of<Components::Gas<Scalar, RawComponent>, RawComponent>::value;
57 };
58 } // end namespace Dumux
59
60 namespace Dumux::Components::Detail {
61 struct DisableStaticAssert {};
62 } // end namespace Dumux::Components::Detail
63
64 namespace Dumux {
65 template<> struct AlwaysFalse<Components::Detail::DisableStaticAssert> : public std::true_type {};
66 }// end namespace Dumux
67
68 namespace Dumux::Components::Detail {
69
70 template<class C> using CompHasNoLiquidEnthalpy = decltype(C::template liquidEnthalpy<DisableStaticAssert>(0.0, 0.0));
71 template<class C> using CompHasNoLiquidDensity = decltype(C::template liquidDensity<DisableStaticAssert>(0.0, 0.0));
72 template<class C> using CompHasNoLiquidThermalCond = decltype(C::template liquidThermalConductivity<DisableStaticAssert>(0.0, 0.0));
73 template<class C> using CompHasNoLiquidHeatCapacity = decltype(C::template liquidHeatCapacity<DisableStaticAssert>(0.0, 0.0));
74 template<class C> using CompHasNoLiquidViscosity = decltype(C::template liquidViscosity<DisableStaticAssert>(0.0, 0.0));
75 template<class C> using CompHasLiquidPressure = decltype(C::liquidPressure(0.0, 0.0));
76
77 template<class C> using CompHasNoGasEnthalpy = decltype(C::template gasEnthalpy<DisableStaticAssert>(0.0, 0.0));
78 template<class C> using CompHasNoGasDensity = decltype(C::template gasDensity<DisableStaticAssert>(0.0, 0.0));
79 template<class C> using CompHasNoGasThermalCond = decltype(C::template gasThermalConductivity<DisableStaticAssert>(0.0, 0.0));
80 template<class C> using CompHasNoGasHeatCapacity = decltype(C::template gasHeatCapacity<DisableStaticAssert>(0.0, 0.0));
81 template<class C> using CompHasNoGasViscosity = decltype(C::template gasViscosity<DisableStaticAssert>(0.0, 0.0));
82 template<class C> using CompHasGasPressure = decltype(C::gasPressure(0.0, 0.0));
83
84 template<class C> constexpr inline bool hasLiquidEnthalpy()
85 { return !Dune::Std::is_detected<CompHasNoLiquidEnthalpy, C>::value && ComponentTraits<C>::hasLiquidState; }
86 template<class C> constexpr inline bool hasLiquidDensity()
87 { return !Dune::Std::is_detected<CompHasNoLiquidDensity, C>::value && ComponentTraits<C>::hasLiquidState; }
88 template<class C> constexpr inline bool hasLiquidThermalConductivity()
89 { return !Dune::Std::is_detected<CompHasNoLiquidThermalCond, C>::value && ComponentTraits<C>::hasLiquidState; }
90 template<class C> constexpr inline bool hasLiquidHeatCapacity()
91 { return !Dune::Std::is_detected<CompHasNoLiquidHeatCapacity, C>::value && ComponentTraits<C>::hasLiquidState; }
92 template<class C> constexpr inline bool hasLiquidViscosity()
93 { return !Dune::Std::is_detected<CompHasNoLiquidViscosity, C>::value && ComponentTraits<C>::hasLiquidState; }
94 template<class C> constexpr inline bool hasLiquidPressure()
95 { return Dune::Std::is_detected<CompHasLiquidPressure, C>::value && ComponentTraits<C>::hasLiquidState; }
96
97 template<class C> constexpr inline bool hasGasEnthalpy()
98 { return !Dune::Std::is_detected<CompHasNoGasEnthalpy, C>::value && ComponentTraits<C>::hasGasState; }
99 template<class C> constexpr inline bool hasGasDensity()
100 { return !Dune::Std::is_detected<CompHasNoGasDensity, C>::value && ComponentTraits<C>::hasGasState; }
101 template<class C> constexpr inline bool hasGasThermalConductivity()
102 { return !Dune::Std::is_detected<CompHasNoGasThermalCond, C>::value && ComponentTraits<C>::hasGasState; }
103 template<class C> constexpr inline bool hasGasHeatCapacity()
104 { return !Dune::Std::is_detected<CompHasNoGasHeatCapacity, C>::value && ComponentTraits<C>::hasGasState; }
105 template<class C> constexpr inline bool hasGasViscosity()
106 { return !Dune::Std::is_detected<CompHasNoGasViscosity, C>::value && ComponentTraits<C>::hasGasState; }
107 template<class C> constexpr inline bool hasGasPressure()
108 { return Dune::Std::is_detected<CompHasGasPressure, C>::value && ComponentTraits<C>::hasGasState; }
109
110 template<class RawComponent, bool useVaporPressure = true>
111 class TabulatedComponentTable
112 {
113 using Scalar = typename RawComponent::Scalar;
114 friend class TabulatedComponent<RawComponent, useVaporPressure>;
115
116 struct GasPolicy
117 {
118 173250 Scalar minP(std::size_t iT) const { return table.minGasPressure(iT); }
119 173250 Scalar maxP(std::size_t iT) const { return table.maxGasPressure(iT); }
120 const TabulatedComponentTable<RawComponent, useVaporPressure>& table;
121 };
122
123 struct LiquidPolicy
124 {
125 172632 Scalar minP(std::size_t iT) const { return table.minLiquidPressure(iT); }
126 172632 Scalar maxP(std::size_t iT) const { return table.maxLiquidPressure(iT); }
127 const TabulatedComponentTable<RawComponent, useVaporPressure>& table;
128 };
129 public:
130 209 void init(Scalar tempMin, Scalar tempMax, std::size_t nTemp,
131 Scalar pressMin, Scalar pressMax, std::size_t nPress)
132 {
133 209 tempMin_ = tempMin;
134 209 tempMax_ = tempMax;
135 209 nTemp_ = nTemp;
136 209 pressMin_ = pressMin;
137 209 pressMax_ = pressMax;
138 209 nPress_ = nPress;
139 209 nDensity_ = nPress;
140
141 // resize & initialize the arrays with NaN
142 assert(std::numeric_limits<Scalar>::has_quiet_NaN);
143 209 const auto NaN = std::numeric_limits<Scalar>::quiet_NaN();
144
145 // initialize vapor pressure array depending on useVaporPressure
146 209 vaporPressure_.resize(nTemp_, NaN);
147 191 tabularizeVaporPressure_();
148
149 if constexpr (ComponentTraits<RawComponent>::hasGasState)
150 {
151 202 minGasDensity_.resize(nTemp_, NaN);
152 202 maxGasDensity_.resize(nTemp_, NaN);
153 202 const std::size_t numEntriesTp = nTemp_*nPress_;
154 202 gasEnthalpy_.resize(numEntriesTp, NaN);
155 202 gasHeatCapacity_.resize(numEntriesTp, NaN);
156 202 gasDensity_.resize(numEntriesTp, NaN);
157 202 gasViscosity_.resize(numEntriesTp, NaN);
158 202 gasThermalConductivity_.resize(numEntriesTp, NaN);
159
160 if constexpr (RawComponent::gasIsCompressible())
161 202 gasPressure_.resize(numEntriesTp, NaN);
162
163 202 minMaxGasDensityInitialized_ = tabularizeMinMaxGasDensity_();
164 202 gasPressureInitialized_ = tabularizeGasPressure_();
165 202 gasEnthalpyInitialized_ = tabularizeGasEnthalpy_();
166 202 gasHeatCapacityInitialized_ = tabularizeGasHeatCapacity_();
167 202 gasDensityInitialized_ = tabularizeGasDensity_();
168 202 gasViscosityInitialized_ = tabularizeGasViscosity_();
169 202 gasThermalConductivityInitialized_ = tabularizeGasThermalConductivity_();
170 }
171
172 if constexpr (ComponentTraits<RawComponent>::hasLiquidState)
173 {
174 191 minLiquidDensity_.resize(nTemp_, NaN);
175 191 maxLiquidDensity_.resize(nTemp_, NaN);
176
177 191 const std::size_t numEntriesTp = nTemp_*nPress_;
178 191 liquidEnthalpy_.resize(numEntriesTp, NaN);
179 191 liquidHeatCapacity_.resize(numEntriesTp, NaN);
180 191 liquidDensity_.resize(numEntriesTp, NaN);
181 191 liquidViscosity_.resize(numEntriesTp, NaN);
182 191 liquidThermalConductivity_.resize(numEntriesTp, NaN);
183
184 if constexpr (RawComponent::liquidIsCompressible())
185 181 liquidPressure_.resize(numEntriesTp, NaN);
186
187 191 minMaxLiquidDensityInitialized_ = tabularizeMinMaxLiquidDensity_();
188 191 liquidPressureInitialized_ = tabularizeLiquidPressure_();
189 191 liquidEnthalpyInitialized_ = tabularizeLiquidEnthalpy_();
190 191 liquidHeatCapacityInitialized_ = tabularizeLiquidHeatCapacity_();
191 191 liquidDensityInitialized_ = tabularizeLiquidDensity_();
192 191 liquidViscosityInitialized_ = tabularizeLiquidViscosity_();
193 191 liquidThermalConductivityInitialized_ = tabularizeLiquidThermalConductivity_();
194 }
195 202 }
196
197 //! returns the index of an entry in a temperature field
198 1248144113 inline Scalar tempIdx(Scalar temperature) const
199 {
200 1248144113 return (nTemp_ - 1)*(temperature - tempMin_)/(tempMax_ - tempMin_);
201 }
202
203 //! returns the index of an entry in a pressure field
204 1788517418 inline Scalar pressLiquidIdx(Scalar pressure, std::size_t tempIdx) const
205 {
206 3570035104 const Scalar plMin = minLiquidPressure(tempIdx);
207 3570035104 const Scalar plMax = maxLiquidPressure(tempIdx);
208 1788517418 return (nPress_ - 1)*(pressure - plMin)/(plMax - plMin);
209 }
210
211 //! returns the index of an entry in a temperature field
212 600576382 inline Scalar pressGasIdx(Scalar pressure, std::size_t tempIdx) const
213 {
214 1192387912 const Scalar pgMin = minGasPressure(tempIdx);
215 1192387912 const Scalar pgMax = maxGasPressure(tempIdx);
216 600576382 return (nPress_ - 1)*(pressure - pgMin)/(pgMax - pgMin);
217 }
218
219 //! returns the index of an entry in a density field
220 1235383 inline Scalar densityLiquidIdx(Scalar density, std::size_t tempIdx) const
221 {
222 1235383 const Scalar densityMin = minLiquidDensity_[tempIdx];
223 1235383 const Scalar densityMax = maxLiquidDensity_[tempIdx];
224 1235383 return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
225 }
226
227 //! returns the index of an entry in a density field
228 205023 inline Scalar densityGasIdx(Scalar density, std::size_t tempIdx) const
229 {
230 205023 const Scalar densityMin = minGasDensity_[tempIdx];
231 205023 const Scalar densityMax = maxGasDensity_[tempIdx];
232 205023 return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
233 }
234
235 //! returns the minimum tabularized liquid pressure at a given temperature index
236 1788603749 inline Scalar minLiquidPressure(int tempIdx) const
237 {
238 using std::max;
239 if (!useVaporPressure)
240 6999744 return pressMin_;
241 else
242
2/2
✓ Branch 0 taken 559399003 times.
✓ Branch 1 taken 1222118683 times.
1781517686 return max(pressMin_, vaporPressure_[tempIdx] / 1.1);
243 }
244
245 //! returns the maximum tabularized liquid pressure at a given temperature index
246 1788603749 inline Scalar maxLiquidPressure(int tempIdx) const
247 {
248 using std::max;
249 if (!useVaporPressure)
250 6999744 return pressMax_;
251 else
252
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1781517686 times.
1781517686 return max(pressMax_, vaporPressure_[tempIdx] * 1.1);
253 }
254
255 //! returns the minimum tabularized gas pressure at a given temperature index
256 600663272 inline Scalar minGasPressure(int tempIdx) const
257 {
258 using std::min;
259 if (!useVaporPressure)
260 8765410 return pressMin_;
261 else
262
2/2
✓ Branch 0 taken 119024155 times.
✓ Branch 1 taken 472787375 times.
591845238 return min(pressMin_, vaporPressure_[tempIdx] / 1.1 );
263 }
264
265 //! returns the maximum tabularized gas pressure at a given temperature index
266 600663272 inline Scalar maxGasPressure(int tempIdx) const
267 {
268 using std::min;
269 if (!useVaporPressure)
270 8765410 return pressMax_;
271 else
272
2/2
✓ Branch 0 taken 591612874 times.
✓ Branch 1 taken 198656 times.
591845238 return min(pressMax_, vaporPressure_[tempIdx] * 1.1);
273 }
274
275 //! returns the minimum tabularized liquid density at a given temperature index
276 inline Scalar minLiquidDensity(int tempIdx) const
277 { return minLiquidDensity_[tempIdx]; }
278
279 //! returns the maximum tabularized liquid density at a given temperature index
280 inline Scalar maxLiquidDensity(int tempIdx) const
281 { return maxLiquidDensity_[tempIdx]; }
282
283 //! returns the minimum tabularized gas density at a given temperature index
284 inline Scalar minGasDensity(int tempIdx) const
285 { return minGasDensity_[tempIdx]; }
286
287 //! returns the maximum tabularized gas density at a given temperature index
288 inline Scalar maxGasDensity(int tempIdx) const
289 { return maxGasDensity_[tempIdx]; }
290
291 1216666941 inline std::size_t nTemp() const { return nTemp_; }
292 1215226522 inline std::size_t nPress() const { return nPress_; }
293 1440406 inline std::size_t nDensity() const { return nDensity_; }
294
295 13 inline Scalar tempMax() const { return tempMax_; }
296 13 inline Scalar tempMin() const { return tempMin_; }
297
298 private:
299
300 //! initializes vapor pressure if useVaporPressure = true
301 template< bool useVP = useVaporPressure, std::enable_if_t<useVP, int> = 0 >
302 182 void tabularizeVaporPressure_()
303 {
304 // fill the temperature-pressure arrays
305 3448 Dumux::parallelFor(nTemp_, [this](std::size_t iT)
306 {
307 14396 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
308 14387 vaporPressure_[iT] = RawComponent::vaporPressure(temperature);
309 });
310 }
311
312 //! if !useVaporPressure, do nothing here
313 template< bool useVP = useVaporPressure, std::enable_if_t<!useVP, int> = 0 >
314 void tabularizeVaporPressure_() {}
315
316 /*!
317 * \brief Initializes property values as function of temperature and pressure.
318 *
319 * \tparam PropFunc Function to evaluate the property prop(T, p)
320 * \tparam MinPFunc Function to evaluate the minimum pressure for a
321 * temperature index (depends on useVaporPressure)
322 * \tparam MaxPFunc Function to evaluate the maximum pressure for a
323 * temperature index (depends on useVaporPressure)
324 *
325 * \param f property function
326 * \param minP function to evaluate minimum pressure for temp idx
327 * \param maxP function to evaluate maximum pressure for temp idx
328 * \param values container to store property values
329 */
330 template<class PropFunc, class Policy>
331 1829 void tabularizeTPArray_(const PropFunc& f, Policy policy, std::vector<typename RawComponent::Scalar>& values) const
332 {
333 579547 Dumux::parallelFor(nTemp_, [this,f,policy,&values](std::size_t iT)
334 {
335 144337 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
336
337
46/72
✓ Branch 0 taken 43 times.
✓ Branch 1 taken 13535 times.
✓ Branch 2 taken 43 times.
✓ Branch 3 taken 13535 times.
✓ Branch 4 taken 43 times.
✓ Branch 5 taken 13535 times.
✓ Branch 6 taken 43 times.
✓ Branch 7 taken 13535 times.
✓ Branch 8 taken 46 times.
✓ Branch 9 taken 13532 times.
✓ Branch 10 taken 13535 times.
✓ Branch 11 taken 43 times.
✓ Branch 12 taken 13535 times.
✓ Branch 13 taken 43 times.
✓ Branch 14 taken 13535 times.
✓ Branch 15 taken 43 times.
✓ Branch 16 taken 13532 times.
✓ Branch 17 taken 46 times.
✓ Branch 18 taken 13535 times.
✓ Branch 19 taken 43 times.
✓ Branch 20 taken 3 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✓ Branch 23 taken 3 times.
✗ Branch 24 not taken.
✓ Branch 25 taken 3 times.
✗ Branch 26 not taken.
✓ Branch 27 taken 3 times.
✗ Branch 28 not taken.
✓ Branch 29 taken 3 times.
✗ Branch 30 not taken.
✓ Branch 31 taken 3 times.
✗ Branch 32 not taken.
✓ Branch 33 taken 3 times.
✗ Branch 34 not taken.
✓ Branch 35 taken 3 times.
✗ Branch 36 not taken.
✓ Branch 37 taken 3 times.
✓ Branch 38 taken 3 times.
✗ Branch 39 not taken.
✓ Branch 40 taken 3 times.
✗ Branch 41 not taken.
✓ Branch 42 taken 3 times.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✓ Branch 45 taken 3 times.
✓ Branch 46 taken 3 times.
✗ Branch 47 not taken.
✓ Branch 48 taken 3 times.
✗ Branch 49 not taken.
✓ Branch 50 taken 3 times.
✗ Branch 51 not taken.
✗ Branch 52 not taken.
✓ Branch 53 taken 803 times.
✗ Branch 54 not taken.
✓ Branch 55 taken 803 times.
✗ Branch 56 not taken.
✓ Branch 57 taken 803 times.
✗ Branch 58 not taken.
✓ Branch 59 taken 803 times.
✗ Branch 60 not taken.
✓ Branch 61 taken 803 times.
✓ Branch 62 taken 803 times.
✗ Branch 63 not taken.
✓ Branch 64 taken 803 times.
✗ Branch 65 not taken.
✓ Branch 66 taken 803 times.
✗ Branch 67 not taken.
✓ Branch 68 taken 803 times.
✗ Branch 69 not taken.
✓ Branch 70 taken 803 times.
✗ Branch 71 not taken.
144286 Scalar pMax = policy.maxP(iT);
338
56/72
✓ Branch 0 taken 12199 times.
✓ Branch 1 taken 1379 times.
✓ Branch 2 taken 12199 times.
✓ Branch 3 taken 1379 times.
✓ Branch 4 taken 12199 times.
✓ Branch 5 taken 1379 times.
✓ Branch 6 taken 12199 times.
✓ Branch 7 taken 1379 times.
✓ Branch 8 taken 12202 times.
✓ Branch 9 taken 1376 times.
✓ Branch 10 taken 1379 times.
✓ Branch 11 taken 12199 times.
✓ Branch 12 taken 1379 times.
✓ Branch 13 taken 12199 times.
✓ Branch 14 taken 1379 times.
✓ Branch 15 taken 12199 times.
✓ Branch 16 taken 1376 times.
✓ Branch 17 taken 12202 times.
✓ Branch 18 taken 1379 times.
✓ Branch 19 taken 12199 times.
✓ Branch 20 taken 3 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✓ Branch 23 taken 3 times.
✗ Branch 24 not taken.
✓ Branch 25 taken 3 times.
✗ Branch 26 not taken.
✓ Branch 27 taken 3 times.
✗ Branch 28 not taken.
✓ Branch 29 taken 3 times.
✗ Branch 30 not taken.
✓ Branch 31 taken 3 times.
✗ Branch 32 not taken.
✓ Branch 33 taken 3 times.
✗ Branch 34 not taken.
✓ Branch 35 taken 3 times.
✗ Branch 36 not taken.
✓ Branch 37 taken 3 times.
✓ Branch 38 taken 3 times.
✗ Branch 39 not taken.
✓ Branch 40 taken 3 times.
✗ Branch 41 not taken.
✓ Branch 42 taken 3 times.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✓ Branch 45 taken 3 times.
✓ Branch 46 taken 3 times.
✗ Branch 47 not taken.
✓ Branch 48 taken 3 times.
✗ Branch 49 not taken.
✓ Branch 50 taken 3 times.
✗ Branch 51 not taken.
✓ Branch 52 taken 692 times.
✓ Branch 53 taken 111 times.
✓ Branch 54 taken 692 times.
✓ Branch 55 taken 111 times.
✓ Branch 56 taken 692 times.
✓ Branch 57 taken 111 times.
✓ Branch 58 taken 692 times.
✓ Branch 59 taken 111 times.
✓ Branch 60 taken 692 times.
✓ Branch 61 taken 111 times.
✓ Branch 62 taken 111 times.
✓ Branch 63 taken 692 times.
✓ Branch 64 taken 111 times.
✓ Branch 65 taken 692 times.
✓ Branch 66 taken 111 times.
✓ Branch 67 taken 692 times.
✓ Branch 68 taken 111 times.
✓ Branch 69 taken 692 times.
✓ Branch 70 taken 111 times.
✓ Branch 71 taken 692 times.
144286 Scalar pMin = policy.minP(iT);
339
106/106
✓ Branch 0 taken 2685209 times.
✓ Branch 1 taken 13680 times.
✓ Branch 2 taken 2685209 times.
✓ Branch 3 taken 13680 times.
✓ Branch 4 taken 2685209 times.
✓ Branch 5 taken 13680 times.
✓ Branch 6 taken 2685209 times.
✓ Branch 7 taken 13680 times.
✓ Branch 8 taken 2665209 times.
✓ Branch 9 taken 13580 times.
✓ Branch 10 taken 2665209 times.
✓ Branch 11 taken 13580 times.
✓ Branch 12 taken 2665209 times.
✓ Branch 13 taken 13580 times.
✓ Branch 14 taken 2665209 times.
✓ Branch 15 taken 13580 times.
✓ Branch 16 taken 2665209 times.
✓ Branch 17 taken 13580 times.
✓ Branch 18 taken 2665209 times.
✓ Branch 19 taken 13580 times.
✓ Branch 20 taken 409 times.
✓ Branch 21 taken 5 times.
✓ Branch 22 taken 409 times.
✓ Branch 23 taken 5 times.
✓ Branch 24 taken 409 times.
✓ Branch 25 taken 5 times.
✓ Branch 26 taken 409 times.
✓ Branch 27 taken 5 times.
✓ Branch 28 taken 9 times.
✓ Branch 29 taken 3 times.
✓ Branch 30 taken 9 times.
✓ Branch 31 taken 3 times.
✓ Branch 32 taken 9 times.
✓ Branch 33 taken 3 times.
✓ Branch 34 taken 9 times.
✓ Branch 35 taken 3 times.
✓ Branch 36 taken 9 times.
✓ Branch 37 taken 3 times.
✓ Branch 38 taken 9 times.
✓ Branch 39 taken 3 times.
✓ Branch 40 taken 9 times.
✓ Branch 41 taken 3 times.
✓ Branch 42 taken 9 times.
✓ Branch 43 taken 3 times.
✓ Branch 44 taken 9 times.
✓ Branch 45 taken 3 times.
✓ Branch 46 taken 9 times.
✓ Branch 47 taken 3 times.
✓ Branch 48 taken 9 times.
✓ Branch 49 taken 3 times.
✓ Branch 50 taken 9 times.
✓ Branch 51 taken 3 times.
✓ Branch 52 taken 9 times.
✓ Branch 53 taken 3 times.
✓ Branch 54 taken 9 times.
✓ Branch 55 taken 3 times.
✓ Branch 56 taken 9 times.
✓ Branch 57 taken 3 times.
✓ Branch 58 taken 9 times.
✓ Branch 59 taken 3 times.
✓ Branch 60 taken 9 times.
✓ Branch 61 taken 3 times.
✓ Branch 62 taken 9 times.
✓ Branch 63 taken 3 times.
✓ Branch 64 taken 9 times.
✓ Branch 65 taken 3 times.
✓ Branch 66 taken 9 times.
✓ Branch 67 taken 3 times.
✓ Branch 68 taken 9 times.
✓ Branch 69 taken 3 times.
✓ Branch 70 taken 9 times.
✓ Branch 71 taken 3 times.
✓ Branch 72 taken 9 times.
✓ Branch 73 taken 3 times.
✓ Branch 74 taken 9 times.
✓ Branch 75 taken 3 times.
✓ Branch 76 taken 9 times.
✓ Branch 77 taken 3 times.
✓ Branch 78 taken 9 times.
✓ Branch 79 taken 3 times.
✓ Branch 80 taken 9 times.
✓ Branch 81 taken 3 times.
✓ Branch 82 taken 9 times.
✓ Branch 83 taken 3 times.
✓ Branch 84 taken 9 times.
✓ Branch 85 taken 3 times.
✓ Branch 86 taken 160009 times.
✓ Branch 87 taken 803 times.
✓ Branch 88 taken 160009 times.
✓ Branch 89 taken 803 times.
✓ Branch 90 taken 160009 times.
✓ Branch 91 taken 803 times.
✓ Branch 92 taken 160009 times.
✓ Branch 93 taken 803 times.
✓ Branch 94 taken 160009 times.
✓ Branch 95 taken 803 times.
✓ Branch 96 taken 160009 times.
✓ Branch 97 taken 803 times.
✓ Branch 98 taken 160009 times.
✓ Branch 99 taken 803 times.
✓ Branch 100 taken 160009 times.
✓ Branch 101 taken 803 times.
✓ Branch 102 taken 160009 times.
✓ Branch 103 taken 803 times.
✓ Branch 104 taken 160009 times.
✓ Branch 105 taken 803 times.
28478414 for (std::size_t iP = 0; iP < nPress_; ++ iP)
340 {
341 28334077 Scalar pressure = iP * (pMax - pMin)/(nPress_ - 1) + pMin;
342 28272625 values[iT + iP*nTemp_] = f(temperature, pressure);
343 }
344 });
345 }
346
347 /*!
348 * \brief Initializes the minimum/maximum densities on the temperature range.
349 *
350 * \tparam RhoFunc Function to evaluate the density rho(T, p)
351 * \tparam MinPFunc Function to evaluate the minimum pressure for a
352 * temperature index (depends on useVaporPressure)
353 * \tparam MaxPFunc Function to evaluate the maximum pressure for a
354 * temperature index (depends on useVaporPressure)
355 *
356 * \param rho density function
357 * \param minP function to evaluate minimum pressure for temp idx
358 * \param maxP function to evaluate maximum pressure for temp idx
359 * \param rhoMin container to store minimum density values
360 * \param rhoMax container to store maximum density values
361 */
362 template<class RhoFunc, class Policy>
363 368 void tabularizeMinMaxRhoArray_(const RhoFunc& rho, Policy policy,
364 std::vector<typename RawComponent::Scalar>& rhoMin,
365 std::vector<typename RawComponent::Scalar>& rhoMax) const
366 {
367 116018 Dumux::parallelFor(nTemp_, [this,rho,policy,&rhoMin,&rhoMax](std::size_t iT)
368 {
369 28901 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
370
371 using std::min;
372
22/28
✓ Branch 0 taken 12199 times.
✓ Branch 1 taken 1379 times.
✓ Branch 2 taken 176 times.
✓ Branch 3 taken 13402 times.
✓ Branch 4 taken 1379 times.
✓ Branch 5 taken 12199 times.
✓ Branch 6 taken 176 times.
✓ Branch 7 taken 13402 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 3 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 3 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 3 times.
✓ Branch 14 taken 3 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 3 times.
✓ Branch 18 taken 3 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 692 times.
✓ Branch 21 taken 111 times.
✓ Branch 22 taken 2 times.
✓ Branch 23 taken 801 times.
✓ Branch 24 taken 111 times.
✓ Branch 25 taken 692 times.
✓ Branch 26 taken 2 times.
✓ Branch 27 taken 801 times.
57646 rhoMin[iT] = rho(temperature, policy.minP(iT));
373
38/48
✓ Branch 0 taken 174 times.
✓ Branch 1 taken 13404 times.
✓ Branch 2 taken 45 times.
✓ Branch 3 taken 13533 times.
✓ Branch 4 taken 174 times.
✓ Branch 5 taken 13404 times.
✓ Branch 6 taken 13533 times.
✓ Branch 7 taken 45 times.
✓ Branch 8 taken 1 times.
✓ Branch 9 taken 2 times.
✓ Branch 10 taken 3 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 1 times.
✓ Branch 13 taken 2 times.
✓ Branch 14 taken 3 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 1 times.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✓ Branch 19 taken 3 times.
✓ Branch 20 taken 1 times.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✓ Branch 23 taken 3 times.
✓ Branch 24 taken 1 times.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✓ Branch 27 taken 3 times.
✓ Branch 28 taken 1 times.
✓ Branch 29 taken 2 times.
✓ Branch 30 taken 3 times.
✗ Branch 31 not taken.
✓ Branch 32 taken 1 times.
✓ Branch 33 taken 2 times.
✗ Branch 34 not taken.
✓ Branch 35 taken 3 times.
✓ Branch 36 taken 1 times.
✓ Branch 37 taken 2 times.
✓ Branch 38 taken 3 times.
✗ Branch 39 not taken.
✓ Branch 40 taken 2 times.
✓ Branch 41 taken 801 times.
✗ Branch 42 not taken.
✓ Branch 43 taken 803 times.
✓ Branch 44 taken 2 times.
✓ Branch 45 taken 801 times.
✓ Branch 46 taken 803 times.
✗ Branch 47 not taken.
29256 rhoMax[iT] = rho(temperature, policy.maxP(min(iT + 1, nTemp_ - 1)));
374 });
375 }
376
377 /*!
378 * \brief Initializes pressure arrays as function of temperature and density.
379 *
380 * \tparam PFunc Function to evaluate the pressure p(T, rho)
381 *
382 * \param pressure container to store pressure values
383 * \param p pressure function p(T, rho)
384 * \param rhoMin container with minimum density values
385 * \param rhoMax container with maximum density values
386 */
387 template<class PFunc>
388 359 void tabularizePressureArray_(std::vector<typename RawComponent::Scalar>& pressure,
389 const PFunc& p,
390 const std::vector<typename RawComponent::Scalar>& rhoMin,
391 const std::vector<typename RawComponent::Scalar>& rhoMax) const
392 {
393 116026 Dumux::parallelFor(nTemp_, [this,p,&pressure,&rhoMin,&rhoMax](std::size_t iT)
394 {
395 28877 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
396
397
14/14
✓ Branch 0 taken 2685209 times.
✓ Branch 1 taken 13680 times.
✓ Branch 2 taken 2665209 times.
✓ Branch 3 taken 13580 times.
✓ Branch 4 taken 409 times.
✓ Branch 5 taken 5 times.
✓ Branch 6 taken 9 times.
✓ Branch 7 taken 3 times.
✓ Branch 8 taken 9 times.
✓ Branch 9 taken 3 times.
✓ Branch 10 taken 160009 times.
✓ Branch 11 taken 803 times.
✓ Branch 12 taken 160009 times.
✓ Branch 13 taken 803 times.
5699740 for (std::size_t iRho = 0; iRho < nDensity_; ++ iRho)
398 {
399 5670863 Scalar density = Scalar(iRho)/(nDensity_ - 1)
400 5670863 * (rhoMax[iT] - rhoMin[iT])
401 5670863 + rhoMin[iT];
402 5670418 pressure[iT + iRho*nTemp_] = p(temperature, density);
403 }
404 });
405 }
406
407 template<class RC = RawComponent>
408 186 bool tabularizeGasEnthalpy_()
409 {
410 if constexpr (Detail::hasGasEnthalpy<RC>())
411 {
412 2845681 const auto gasEnth = [] (auto T, auto p) { return RC::gasEnthalpy(T, p); };
413 372 tabularizeTPArray_(gasEnth, GasPolicy{ *this }, gasEnthalpy_);
414 return true;
415 }
416
417 return false;
418 }
419
420 template<class RC = RawComponent>
421 180 bool tabularizeLiquidEnthalpy_()
422 {
423 if constexpr (Detail::hasLiquidEnthalpy<RC>())
424 {
425 2825245 const auto liqEnth = [] (auto T, auto p) { return RC::liquidEnthalpy(T, p); };
426 360 tabularizeTPArray_(liqEnth, LiquidPolicy{ *this }, liquidEnthalpy_);
427 return true;
428 }
429
430 return false;
431 }
432
433 template<class RC = RawComponent>
434 183 bool tabularizeGasHeatCapacity_()
435 {
436 if constexpr (Detail::hasGasHeatCapacity<RC>())
437 {
438 2845645 const auto gasHC = [] (auto T, auto p) { return RC::gasHeatCapacity(T, p); };
439 183 tabularizeTPArray_(gasHC, GasPolicy{ *this }, gasHeatCapacity_);
440 return true;
441 }
442
443 return false;
444 }
445
446 template<class RC = RawComponent>
447 180 bool tabularizeLiquidHeatCapacity_()
448 {
449 if constexpr (Detail::hasLiquidHeatCapacity<RC>())
450 {
451 2825227 const auto liqHC = [] (auto T, auto p) { return RC::liquidHeatCapacity(T, p); };
452 180 tabularizeTPArray_(liqHC, LiquidPolicy{ *this }, liquidHeatCapacity_);
453 return true;
454 }
455
456 return false;
457 }
458
459 template<class RC = RawComponent>
460 188 bool tabularizeMinMaxGasDensity_()
461 {
462 if constexpr (Detail::hasGasDensity<RC>())
463 {
464
14/14
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 1 times.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 2 times.
✓ Branch 8 taken 1 times.
✓ Branch 9 taken 2 times.
✓ Branch 10 taken 1 times.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 1 times.
✓ Branch 13 taken 2 times.
28879 const auto gasRho = [] (auto T, auto p) { return RC::gasDensity(T, p); };
465 188 tabularizeMinMaxRhoArray_(gasRho, GasPolicy{ *this }, minGasDensity_, maxGasDensity_);
466 return true;
467 }
468
469 return false;
470 }
471
472 template<class RC = RawComponent>
473 180 bool tabularizeMinMaxLiquidDensity_()
474 {
475 if constexpr (Detail::hasGasEnthalpy<RC>())
476 {
477
6/6
✓ Branch 1 taken 1 times.
✓ Branch 2 taken 2 times.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 2 times.
✓ Branch 8 taken 1 times.
✓ Branch 9 taken 2 times.
28773 const auto liqRho = [] (auto T, auto p) { return RC::liquidDensity(T, p); };
478 180 tabularizeMinMaxRhoArray_(liqRho, LiquidPolicy{ *this }, minLiquidDensity_, maxLiquidDensity_);
479 return true;
480 }
481
482 return false;
483 }
484
485 template<class RC = RawComponent>
486 183 bool tabularizeGasPressure_()
487 {
488 // pressure is only defined if the gas is compressible (this is usually the case)
489 if constexpr (Detail::hasGasPressure<RC>() && RC::gasIsCompressible())
490 {
491 2845654 const auto gasPFunc = [] (auto T, auto rho) { return RC::gasPressure(T, rho); };
492 183 tabularizePressureArray_(gasPressure_, gasPFunc, minGasDensity_, maxGasDensity_);
493 return true;
494 }
495
496 return false;
497 }
498
499 template<class RC = RawComponent>
500 176 bool tabularizeLiquidPressure_()
501 {
502 // pressure is only defined if the liquid is compressible (this is often not the case)
503 if constexpr (Detail::hasLiquidPressure<RC>() && RC::liquidIsCompressible())
504 {
505 2825209 const auto liqPFunc = [] (auto T, auto rho) { return RC::liquidPressure(T, rho); };
506 176 tabularizePressureArray_(liquidPressure_, liqPFunc, minLiquidDensity_, maxLiquidDensity_);
507 return true;
508 }
509
510 return false;
511 }
512
513 template<class RC = RawComponent>
514 188 bool tabularizeGasDensity_()
515 {
516 if constexpr (Detail::hasGasDensity<RC>())
517 {
518 2845699 const auto gasRho = [] (auto T, auto p) { return RC::gasDensity(T, p); };
519 188 tabularizeTPArray_(gasRho, GasPolicy{ *this }, gasDensity_);
520 return true;
521 }
522
523 return false;
524 }
525
526 template<class RC = RawComponent>
527 182 bool tabularizeLiquidDensity_()
528 {
529 if constexpr (Detail::hasLiquidDensity<RC>())
530 {
531 // TODO: we could get rid of the lambdas and pass the functor directly. But,
532 // currently Brine is a component (and not a fluid system) expecting a
533 // third argument with a default, which cannot be wrapped in a function pointer.
534 // For this reason we have to wrap this into a lambda here.
535 2825236 const auto liqRho = [] (auto T, auto p) { return RC::liquidDensity(T, p); };
536 182 tabularizeTPArray_(liqRho, LiquidPolicy{ *this }, liquidDensity_);
537 return true;
538 }
539
540 return false;
541 }
542
543 template<class RC = RawComponent>
544 186 bool tabularizeGasViscosity_()
545 {
546 if constexpr (Detail::hasGasViscosity<RC>())
547 {
548 2845672 const auto gasVisc = [] (auto T, auto p) { return RC::gasViscosity(T, p); };
549 186 tabularizeTPArray_(gasVisc, GasPolicy{ *this }, gasViscosity_);
550 return true;
551 }
552
553 return false;
554 }
555
556 template<class RC = RawComponent>
557 182 bool tabularizeLiquidViscosity_()
558 {
559 if constexpr (Detail::hasLiquidViscosity<RC>())
560 {
561 2825236 const auto liqVisc = [] (auto T, auto p) { return RC::liquidViscosity(T, p); };
562 182 tabularizeTPArray_(liqVisc, LiquidPolicy{ *this }, liquidViscosity_);
563 return true;
564 }
565
566 return false;
567 }
568
569 template<class RC = RawComponent>
570 182 bool tabularizeGasThermalConductivity_()
571 {
572 if constexpr (Detail::hasGasThermalConductivity<RC>())
573 {
574 2825227 const auto gasTC = [] (auto T, auto p) { return RC::gasThermalConductivity(T, p); };
575 182 tabularizeTPArray_(gasTC, GasPolicy{ *this }, gasThermalConductivity_);
576 return true;
577 }
578
579 return false;
580 }
581
582 template<class RC = RawComponent>
583 180 bool tabularizeLiquidThermalConductivity_()
584 {
585 if constexpr (Detail::hasLiquidThermalConductivity<RC>())
586 {
587 2825209 const auto liqTC = [] (auto T, auto p) { return RC::liquidThermalConductivity(T, p); };
588 180 tabularizeTPArray_(liqTC, LiquidPolicy{ *this }, liquidThermalConductivity_);
589 return true;
590 }
591
592 return false;
593 }
594
595 private:
596 // 1D fields with the temperature as degree of freedom
597 std::vector<Scalar> vaporPressure_;
598
599 std::vector<Scalar> minLiquidDensity_;
600 std::vector<Scalar> maxLiquidDensity_;
601 bool minMaxLiquidDensityInitialized_;
602
603 std::vector<Scalar> minGasDensity_;
604 std::vector<Scalar> maxGasDensity_;
605 bool minMaxGasDensityInitialized_;
606
607 // 2D fields with the temperature and pressure as degrees of freedom
608 std::vector<Scalar> gasEnthalpy_;
609 std::vector<Scalar> liquidEnthalpy_;
610 bool gasEnthalpyInitialized_;
611 bool liquidEnthalpyInitialized_;
612
613 std::vector<Scalar> gasHeatCapacity_;
614 std::vector<Scalar> liquidHeatCapacity_;
615 bool gasHeatCapacityInitialized_;
616 bool liquidHeatCapacityInitialized_;
617
618 std::vector<Scalar> gasDensity_;
619 std::vector<Scalar> liquidDensity_;
620 bool gasDensityInitialized_;
621 bool liquidDensityInitialized_;
622
623 std::vector<Scalar> gasViscosity_;
624 std::vector<Scalar> liquidViscosity_;
625 bool gasViscosityInitialized_;
626 bool liquidViscosityInitialized_;
627
628 std::vector<Scalar> gasThermalConductivity_;
629 std::vector<Scalar> liquidThermalConductivity_;
630 bool gasThermalConductivityInitialized_;
631 bool liquidThermalConductivityInitialized_;
632
633 // 2D fields with the temperature and density as degrees of freedom
634 std::vector<Scalar> gasPressure_;
635 std::vector<Scalar> liquidPressure_;
636 bool gasPressureInitialized_;
637 bool liquidPressureInitialized_;
638
639 // temperature, pressure and density ranges
640 Scalar tempMin_;
641 Scalar tempMax_;
642 std::size_t nTemp_;
643
644 Scalar pressMin_;
645 Scalar pressMax_;
646 std::size_t nPress_;
647
648 Scalar densityMin_;
649 Scalar densityMax_;
650 std::size_t nDensity_;
651 };
652
653 } // end namespace Dumux::Components::Detail
654
655 namespace Dumux::Components {
656
657 /*!
658 * \ingroup Components
659 * \brief Tabulates all thermodynamic properties of a given component
660 *
661 * At the moment, this class can only handle the sub-critical fluids
662 * since it tabulates along the vapor pressure curve.
663 *
664 * \tparam Scalar The type used for scalar values
665 * \tparam RawComponent The component which ought to be tabulated
666 * \tparam useVaporPressure If set to true, the min/max pressure
667 * values for gas&liquid phase will be set
668 * depending on the vapor pressure.
669 */
670 template <class RawComponent, bool useVaporPressure=true>
671 class TabulatedComponent
672 {
673 using ThisType = TabulatedComponent<RawComponent, useVaporPressure>;
674 using Table = Detail::TabulatedComponentTable<RawComponent, useVaporPressure>;
675
676 struct InterpolatePolicy
677 {
678 using Scalar = typename RawComponent::Scalar;
679
680 const Table& table;
681
682 1248144113 Scalar tempIdx(Scalar T) const { return table.tempIdx(T); }
683 13 Scalar tempMin() const { return table.tempMin(); }
684 13 Scalar tempMax() const { return table.tempMax(); }
685 1248144126 std::size_t nTemp() const { return table.nTemp(); }
686 1215226522 std::size_t nPress() const { return table.nPress(); }
687 1440406 std::size_t nDensity() const { return table.nDensity(); }
688 };
689
690 struct InterpolateGasPolicy : public InterpolatePolicy
691 {
692 using Scalar = typename RawComponent::Scalar;
693 304670617 Scalar pressIdx(Scalar p, std::size_t tempIdx) const { return this->table.pressGasIdx(p, tempIdx); }
694 205023 Scalar rhoIdx(Scalar rho, std::size_t tempIdx) const { return this->table.densityGasIdx(rho, tempIdx); }
695 2 Scalar minP(int tempIdx) const { return this->table.minGasPressure(tempIdx); }
696 2 Scalar maxP(int tempIdx) const { return this->table.maxGasPressure(tempIdx); }
697 Scalar minRho(int tempIdx) const { return this->table.minGasDensity(tempIdx); }
698 Scalar maxRho(int tempIdx) const { return this->table.maxGasDensity(tempIdx); }
699 };
700
701 struct InterpolateLiquidPolicy : public InterpolatePolicy
702 {
703 using Scalar = typename RawComponent::Scalar;
704 897758575 Scalar pressIdx(Scalar p, std::size_t tempIdx) const { return this->table.pressLiquidIdx(p, tempIdx); }
705 1235383 Scalar rhoIdx(Scalar rho, std::size_t tempIdx) const { return this->table.densityLiquidIdx(rho, tempIdx); }
706 24 Scalar minP(int tempIdx) const { return this->table.minLiquidPressure(tempIdx); }
707 24 Scalar maxP(int tempIdx) const { return this->table.maxLiquidPressure(tempIdx); }
708 Scalar minRho(int tempIdx) const { return this->table.minLiquidDensity(tempIdx); }
709 Scalar maxRho(int tempIdx) const { return this->table.maxLiquidDensity(tempIdx); }
710 };
711 public:
712 //! export scalar type
713 using Scalar = typename RawComponent::Scalar;
714
715 //! state that we are tabulated
716 static constexpr bool isTabulated = true;
717
718 /*!
719 * \brief Initialize the tables.
720 *
721 * \param tempMin The minimum of the temperature range in \f$\mathrm{[K]}\f$
722 * \param tempMax The maximum of the temperature range in \f$\mathrm{[K]}\f$
723 * \param nTemp The number of entries/steps within the temperature range
724 * \param pressMin The minimum of the pressure range in \f$\mathrm{[Pa]}\f$
725 * \param pressMax The maximum of the pressure range in \f$\mathrm{[Pa]}\f$
726 * \param nPress The number of entries/steps within the pressure range
727 */
728 216 static void init(Scalar tempMin, Scalar tempMax, std::size_t nTemp,
729 Scalar pressMin, Scalar pressMax, std::size_t nPress)
730 {
731 #ifndef NDEBUG
732 216 warningPrinted_ = false;
733 #endif
734 std::cout << "-------------------------------------------------------------------------\n"
735 << "Initializing tables for the " << RawComponent::name()
736
3/6
✓ Branch 3 taken 193 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 193 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 193 times.
✗ Branch 10 not taken.
428 << " fluid properties (" << nTemp*nPress << " entries).\n"
737
2/4
✓ Branch 1 taken 193 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 193 times.
✗ Branch 6 not taken.
216 << "Temperature -> min: " << std::scientific << std::setprecision(3)
738
4/8
✓ Branch 1 taken 193 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 193 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 193 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 193 times.
✗ Branch 11 not taken.
216 << tempMin << ", max: " << tempMax << ", n: " << nTemp << '\n'
739
3/6
✓ Branch 1 taken 193 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 193 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 193 times.
✗ Branch 9 not taken.
216 << "Pressure -> min: " << std::scientific << std::setprecision(3)
740
4/8
✓ Branch 1 taken 193 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 193 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 193 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 193 times.
✗ Branch 11 not taken.
216 << pressMin << ", max: " << pressMax << ", n: " << nPress << '\n'
741
3/6
✓ Branch 1 taken 193 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 193 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 193 times.
✗ Branch 8 not taken.
216 << "-------------------------------------------------------------------------" << std::endl;
742
743 216 table().init(tempMin, tempMax, nTemp, pressMin, pressMax, nPress);
744
745 216 initialized_ = true;
746 216 }
747
748 /*!
749 * \brief A human readable name for the component.
750 */
751 458184 static std::string name()
752 458176 { return RawComponent::name(); }
753
754 /*!
755 * \brief The molar mass in \f$\mathrm{[kg/mol]}\f$ of the component.
756 */
757 static constexpr Scalar molarMass()
758 { return RawComponent::molarMass(); }
759
760 /*!
761 * \brief Returns the critical temperature in \f$\mathrm{[K]}\f$ of the component.
762 */
763 static Scalar criticalTemperature()
764 { return RawComponent::criticalTemperature(); }
765
766 /*!
767 * \brief Returns the critical pressure in \f$\mathrm{[Pa]}\f$ of the component.
768 */
769 static Scalar criticalPressure()
770 { return RawComponent::criticalPressure(); }
771
772 /*!
773 * \brief Returns the temperature in \f$\mathrm{[K]}\f$ at the component's triple point.
774 */
775 static Scalar tripleTemperature()
776 { return RawComponent::tripleTemperature(); }
777
778 /*!
779 * \brief Returns the pressure in \f$\mathrm{[Pa]}\f$ at the component's triple point.
780 */
781 static Scalar triplePressure()
782 { return RawComponent::triplePressure(); }
783
784 /*!
785 * \brief The vapor pressure in \f$\mathrm{[Pa]}\f$ of the component at a given
786 * temperature.
787 *
788 * \param T temperature of component
789 */
790 31477185 static Scalar vaporPressure(Scalar T)
791 {
792 using std::isnan;
793 31477185 Scalar result = interpolateT_(table().vaporPressure_, T, InterpolatePolicy{{ table() }});
794
2/2
✓ Branch 0 taken 362047 times.
✓ Branch 1 taken 31115138 times.
31477185 if (isnan(result))
795 362047 return RawComponent::vaporPressure(T);
796 return result;
797 }
798
799 /*!
800 * \brief The vapor pressure in \f$\mathrm{[Pa]}\f$ of the component at a given
801 * temperature.
802 *
803 * The method is only called by the sequential flash, so tabulating is omitted.
804 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
805 */
806 101966 static Scalar vaporTemperature(Scalar pressure)
807 {
808 101966 return RawComponent::vaporTemperature(pressure);
809 }
810
811 /*!
812 * \brief Specific enthalpy of the gas \f$\mathrm{[J/kg]}\f$.
813 *
814 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
815 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
816 */
817 127555272 static const Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
818 {
819 127555272 Scalar result = interpolateTP_(table().gasEnthalpy_, temperature, pressure, InterpolateGasPolicy{{ table() }});
820 using std::isnan;
821
2/2
✓ Branch 0 taken 3148076 times.
✓ Branch 1 taken 124407196 times.
127555272 if (isnan(result))
822 {
823
1/2
✓ Branch 3 taken 3148076 times.
✗ Branch 4 not taken.
3148076 printWarningTP_("gasEnthalpy", temperature, pressure, InterpolateGasPolicy{{ table() }});
824 3148076 return RawComponent::gasEnthalpy(temperature, pressure);
825 }
826 return result;
827 }
828
829 /*!
830 * \brief Specific enthalpy of the liquid \f$\mathrm{[J/kg]}\f$.
831 *
832 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
833 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
834 */
835 148424949 static const Scalar liquidEnthalpy(Scalar temperature, Scalar pressure)
836 {
837 148424949 Scalar result = interpolateTP_(table().liquidEnthalpy_, temperature, pressure, InterpolateLiquidPolicy{{ table() }});
838 using std::isnan;
839
2/2
✓ Branch 0 taken 216246 times.
✓ Branch 1 taken 148208703 times.
148424949 if (isnan(result))
840 {
841
1/2
✓ Branch 3 taken 216246 times.
✗ Branch 4 not taken.
216246 printWarningTP_("liquidEnthalpy", temperature, pressure, InterpolateLiquidPolicy{{ table() }});
842 216246 return RawComponent::liquidEnthalpy(temperature, pressure);
843 }
844 return result;
845 }
846
847 /*!
848 * \brief Specific isobaric heat capacity of the gas \f$\mathrm{[J/(kg*K)]}\f$.
849 *
850 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
851 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
852 */
853 18189597 static const Scalar gasHeatCapacity(Scalar temperature, Scalar pressure)
854 {
855 18189597 Scalar result = interpolateTP_(table().gasHeatCapacity_, temperature, pressure, InterpolateGasPolicy{{ table() }});
856 using std::isnan;
857
2/2
✓ Branch 0 taken 301 times.
✓ Branch 1 taken 18189296 times.
18189597 if (isnan(result))
858 {
859
1/2
✓ Branch 3 taken 301 times.
✗ Branch 4 not taken.
301 printWarningTP_("gasHeatCapacity", temperature, pressure, InterpolateGasPolicy{{ table() }});
860 301 return RawComponent::gasHeatCapacity(temperature, pressure);
861 }
862 return result;
863 }
864
865 /*!
866 * \brief Specific isobaric heat capacity of the liquid \f$\mathrm{[J/(kg*K)]}\f$.
867 *
868 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
869 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
870 */
871 2178072 static const Scalar liquidHeatCapacity(Scalar temperature, Scalar pressure)
872 {
873 2178072 Scalar result = interpolateTP_(table().liquidHeatCapacity_, temperature, pressure, InterpolateLiquidPolicy{{ table() }});
874 using std::isnan;
875
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 2178072 times.
2178072 if (isnan(result))
876 {
877 printWarningTP_("liquidHeatCapacity", temperature, pressure, InterpolateLiquidPolicy{{ table() }});
878 return RawComponent::liquidHeatCapacity(temperature, pressure);
879 }
880 return result;
881 }
882
883 /*!
884 * \brief Specific internal energy of the gas \f$\mathrm{[J/kg]}\f$.
885 *
886 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
887 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
888 */
889 205023 static const Scalar gasInternalEnergy(Scalar temperature, Scalar pressure)
890 {
891 205023 return gasEnthalpy(temperature, pressure) - pressure/gasDensity(temperature, pressure);
892 }
893
894 /*!
895 * \brief Specific internal energy of the liquid \f$\mathrm{[J/kg]}\f$.
896 *
897 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
898 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
899 */
900 1235383 static const Scalar liquidInternalEnergy(Scalar temperature, Scalar pressure)
901 {
902 1235383 return liquidEnthalpy(temperature, pressure) - pressure/liquidDensity(temperature, pressure);
903 }
904
905 /*!
906 * \brief The pressure of gas in \f$\mathrm{[Pa]}\f$ at a given density and temperature.
907 *
908 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
909 * \param density density of component in \f$\mathrm{[kg/m^3]}\f$
910 */
911 205023 static Scalar gasPressure(Scalar temperature, Scalar density)
912 {
913 205023 Scalar result = interpolateTRho_(table().gasPressure_, temperature, density, InterpolateGasPolicy{{ table() }});
914 using std::isnan;
915
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 205023 times.
205023 if (isnan(result))
916 {
917 printWarningTRho_("gasPressure", temperature, density, InterpolateGasPolicy{{ table() }});
918 return RawComponent::gasPressure(temperature, density);
919 }
920 return result;
921 }
922
923 /*!
924 * \brief The pressure of liquid in \f$\mathrm{[Pa]}\f$ at a given density and temperature.
925 *
926 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
927 * \param density density of component in \f$\mathrm{[kg/m^3]}\f$
928 */
929 1235383 static Scalar liquidPressure(Scalar temperature, Scalar density)
930 {
931 1235383 Scalar result = interpolateTRho_(table().liquidPressure_, temperature, density, InterpolateLiquidPolicy{{ table() }});
932 using std::isnan;
933
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1235383 times.
1235383 if (isnan(result))
934 {
935 printWarningTRho_("liquidPressure", temperature, density, InterpolateLiquidPolicy{{ table() }});
936 return RawComponent::liquidPressure(temperature, density);
937 }
938 return result;
939 }
940
941 /*!
942 * \brief Returns true if the gas phase is assumed to be compressible
943 */
944 static constexpr bool gasIsCompressible()
945 { return RawComponent::gasIsCompressible(); }
946
947 /*!
948 * \brief Returns true if the liquid phase is assumed to be compressible
949 */
950 static constexpr bool liquidIsCompressible()
951 { return RawComponent::liquidIsCompressible(); }
952
953 /*!
954 * \brief Returns true if the gas phase is assumed to be ideal
955 */
956 static constexpr bool gasIsIdeal()
957 { return RawComponent::gasIsIdeal(); }
958
959
960 /*!
961 * \brief The density of gas at a given pressure and temperature
962 * \f$\mathrm{[kg/m^3]}\f$.
963 *
964 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
965 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
966 */
967 147404032 static Scalar gasDensity(Scalar temperature, Scalar pressure)
968 {
969 147404032 Scalar result = interpolateTP_(table().gasDensity_, temperature, pressure, InterpolateGasPolicy{{ table() }});
970 using std::isnan;
971
2/2
✓ Branch 0 taken 7862587 times.
✓ Branch 1 taken 139541445 times.
147404032 if (isnan(result))
972 {
973
1/2
✓ Branch 3 taken 7862587 times.
✗ Branch 4 not taken.
7862587 printWarningTP_("gasDensity", temperature, pressure, InterpolateGasPolicy{{ table() }});
974 7862587 return RawComponent::gasDensity(temperature, pressure);
975 }
976 return result;
977 }
978
979 /*!
980 * \brief The molar density of gas in \f$\mathrm{[mol/m^3]}\f$
981 * at a given pressure and temperature.
982 *
983 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
984 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
985 *
986 */
987 68179966 static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
988 68179966 { return gasDensity(temperature, pressure)/molarMass(); }
989
990 /*!
991 * \brief The density of liquid at a given pressure and
992 * temperature \f$\mathrm{[kg/m^3]}\f$.
993 *
994 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
995 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
996 */
997 583160543 static Scalar liquidDensity(Scalar temperature, Scalar pressure)
998 {
999 583160543 Scalar result = interpolateTP_(table().liquidDensity_, temperature, pressure, InterpolateLiquidPolicy{{ table() }});
1000 using std::isnan;
1001
2/2
✓ Branch 0 taken 868085 times.
✓ Branch 1 taken 582292458 times.
583160543 if (isnan(result))
1002 {
1003
1/2
✓ Branch 3 taken 868085 times.
✗ Branch 4 not taken.
868085 printWarningTP_("liquidDensity", temperature, pressure, InterpolateLiquidPolicy{{ table() }});
1004 868085 return RawComponent::liquidDensity(temperature, pressure);
1005 }
1006
1007 return result;
1008 }
1009
1010 /*!
1011 * \brief The molar density of liquid in \f$\mathrm{[mol/m^3]}\f$
1012 * at a given pressure and temperature.
1013 *
1014 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
1015 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
1016 *
1017 */
1018 103560161 static Scalar liquidMolarDensity(Scalar temperature, Scalar pressure)
1019 103560161 { return liquidDensity(temperature, pressure)/molarMass(); }
1020
1021 /*!
1022 * \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of gas.
1023 *
1024 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
1025 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
1026 */
1027 5814033 static Scalar gasViscosity(Scalar temperature, Scalar pressure)
1028 {
1029 5814033 Scalar result = interpolateTP_(table().gasViscosity_, temperature, pressure, InterpolateGasPolicy{{ table() }});
1030 using std::isnan;
1031
2/2
✓ Branch 0 taken 105 times.
✓ Branch 1 taken 5813928 times.
5814033 if (isnan(result))
1032 {
1033
1/2
✓ Branch 3 taken 105 times.
✗ Branch 4 not taken.
105 printWarningTP_("gasViscosity", temperature, pressure, InterpolateGasPolicy{{ table() }});
1034 105 return RawComponent::gasViscosity(temperature, pressure);
1035 }
1036 return result;
1037 }
1038
1039 /*!
1040 * \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of liquid.
1041 *
1042 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
1043 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
1044 */
1045 70313493 static Scalar liquidViscosity(Scalar temperature, Scalar pressure)
1046 {
1047 70313493 Scalar result = interpolateTP_(table().liquidViscosity_, temperature, pressure, InterpolateLiquidPolicy{{ table() }});
1048 using std::isnan;
1049
2/2
✓ Branch 0 taken 352712 times.
✓ Branch 1 taken 69960781 times.
70313493 if (isnan(result))
1050 {
1051
1/2
✓ Branch 3 taken 352712 times.
✗ Branch 4 not taken.
352712 printWarningTP_("liquidViscosity",temperature, pressure, InterpolateLiquidPolicy{{ table() }});
1052 352712 return RawComponent::liquidViscosity(temperature, pressure);
1053 }
1054 return result;
1055 }
1056
1057 /*!
1058 * \brief The thermal conductivity of gaseous water \f$\mathrm{[W/(m*K)]}\f$.
1059 *
1060 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
1061 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
1062 */
1063 16858156 static Scalar gasThermalConductivity(Scalar temperature, Scalar pressure)
1064 {
1065 16858156 Scalar result = interpolateTP_(table().gasThermalConductivity_, temperature, pressure, InterpolateGasPolicy{{ table() }});
1066 using std::isnan;
1067
2/2
✓ Branch 0 taken 139404 times.
✓ Branch 1 taken 16718752 times.
16858156 if (isnan(result))
1068 {
1069
1/2
✓ Branch 3 taken 139404 times.
✗ Branch 4 not taken.
139404 printWarningTP_("gasThermalConductivity", temperature, pressure, InterpolateGasPolicy{{ table() }});
1070 139404 return RawComponent::gasThermalConductivity(temperature, pressure);
1071 }
1072 return result;
1073 }
1074
1075 /*!
1076 * \brief The thermal conductivity of liquid water \f$\mathrm{[W/(m*K)]}\f$.
1077 *
1078 * \param temperature temperature of component in \f$\mathrm{[K]}\f$
1079 * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
1080 */
1081 95328375 static Scalar liquidThermalConductivity(Scalar temperature, Scalar pressure)
1082 {
1083 95328375 Scalar result = interpolateTP_(table().liquidThermalConductivity_, temperature, pressure, InterpolateLiquidPolicy{{ table() }});
1084 using std::isnan;
1085
2/2
✓ Branch 0 taken 209814 times.
✓ Branch 1 taken 95118561 times.
95328375 if (isnan(result))
1086 {
1087
1/2
✓ Branch 3 taken 209814 times.
✗ Branch 4 not taken.
209814 printWarningTP_("liquidThermalConductivity", temperature, pressure, InterpolateLiquidPolicy{{ table() }});
1088 209814 return RawComponent::liquidThermalConductivity(temperature, pressure);
1089 }
1090 return result;
1091 }
1092
1093
1094 private:
1095 //! prints a warning if the result is not in range or the table has not been initialized
1096 template<class Policy>
1097 25594658 static void printWarningTP_(const std::string& quantity, Scalar T, Scalar p, Policy policy)
1098 {
1099 #ifndef NDEBUG
1100
2/2
✓ Branch 0 taken 31 times.
✓ Branch 1 taken 12797299 times.
25594658 if (warningPrinted_)
1101 return;
1102
1103
2/2
✓ Branch 0 taken 18 times.
✓ Branch 1 taken 13 times.
61 if (!initialized_)
1104 std::cerr << "Warning: tabulated component '" << name()
1105 << "' has not been initialized. "
1106
2/4
✓ Branch 2 taken 18 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 18 times.
✗ Branch 6 not taken.
105 << "Call FluidSystem::init() to use the tabulation in order to reduce runtime. \n";
1107 else
1108 26 std::cerr << "Warning: "<<quantity<<"(T="<<T<<", p="<<p<<") of component '"<<name()
1109
7/12
✓ Branch 2 taken 13 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 13 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 13 times.
✗ Branch 9 not taken.
✓ Branch 11 taken 13 times.
✗ Branch 12 not taken.
✓ Branch 14 taken 13 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 5 times.
✓ Branch 17 taken 8 times.
78 << "' is outside tabulation range: ("<< policy.tempMin() <<"<=T<="<< policy.tempMax() <<"), ("
1110
6/12
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 13 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 12 times.
✓ Branch 9 taken 13 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 13 times.
✗ Branch 13 not taken.
✗ Branch 8 not taken.
✗ Branch 11 not taken.
28 << policy.minP(0) <<"<=p<=" << policy.maxP(policy.nTemp()-1) <<"). "
1111
3/6
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 13 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 13 times.
✗ Branch 8 not taken.
26 << "Forwarded to FluidSystem for direct evaluation of "<<quantity<<". \n";
1112 61 warningPrinted_ = true;
1113 #endif
1114 }
1115
1116 //! prints a warning if the result is not in range or the table has not been initialized
1117 template<class Policy>
1118 static void printWarningTRho_(const std::string& quantity, Scalar T, Scalar rho, Policy policy)
1119 {
1120 #ifndef NDEBUG
1121 if (warningPrinted_)
1122 return;
1123
1124 if (!initialized_)
1125 std::cerr << "Warning: tabulated component '" << name()
1126 << "' has not been initialized. "
1127 << "Call FluidSystem::init() to use the tabulation in order to reduce runtime. \n";
1128 else
1129 {
1130 const auto [densityMin, densityMax] = [&]
1131 {
1132 const Scalar alphaT = policy.tempIdx(T);
1133 using std::clamp;
1134 const auto iT = clamp<int>(static_cast<int>(alphaT), 0, policy.nTemp() - 2);
1135 return std::make_pair( policy.minRho(iT), policy.maxRho(iT) );
1136 }();
1137
1138 std::cerr << "Warning: "<<quantity<<"(T="<<T<<", density="<<rho<<") of component '"<<name()
1139 << "' is outside tabulation range: ("<< policy.tempMin() <<"<=T<="<< policy.tempMax() <<"), ("
1140 << densityMin <<"<=density<=" << densityMax <<"). "
1141 << "Forwarded to FluidSystem for direct evaluation of "<<quantity<<". \n";
1142 }
1143 warningPrinted_ = true;
1144 #endif
1145 }
1146
1147 //! returns an interpolated value depending on temperature
1148 template<class Policy>
1149 31477185 static Scalar interpolateT_(const std::vector<Scalar>& values, Scalar T, Policy policy)
1150 {
1151 31477185 Scalar alphaT = policy.tempIdx(T);
1152 31477185 const auto nTemp = policy.nTemp();
1153
4/4
✓ Branch 0 taken 31156789 times.
✓ Branch 1 taken 320396 times.
✓ Branch 2 taken 31115138 times.
✓ Branch 3 taken 41651 times.
31477185 if (alphaT < 0 - 1e-7*nTemp || alphaT >= nTemp - 1 + 1e-7*nTemp)
1154 return std::numeric_limits<Scalar>::quiet_NaN();
1155
1156 using std::clamp;
1157
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 31115138 times.
31115138 const auto iT = clamp<int>(static_cast<int>(alphaT), 0, nTemp - 2);
1158 31115138 alphaT -= iT;
1159
1160 31115138 return values[iT ]*(1 - alphaT) +
1161 31115138 values[iT + 1]*( alphaT);
1162 }
1163
1164 //! returns an interpolated value depending on temperature and pressure
1165 template<class Policy>
1166 2254229975 static Scalar interpolateTP_(const std::vector<Scalar>& values, const Scalar T, const Scalar p, Policy policy)
1167 {
1168 2254229975 const auto nTemp = policy.nTemp();
1169 2254229975 const auto nPress = policy.nPress();
1170
1171 2254229975 Scalar alphaT = policy.tempIdx(T);
1172
4/4
✓ Branch 0 taken 1213071101 times.
✓ Branch 1 taken 2155421 times.
✓ Branch 2 taken 1202429192 times.
✓ Branch 3 taken 10641909 times.
2254229975 if (alphaT < 0 - 1e-7*nTemp || alphaT >= nTemp - 1 + 1e-7*nTemp)
1173 return std::numeric_limits<Scalar>::quiet_NaN();
1174
1175 using std::clamp;
1176
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1202429192 times.
2228635317 const auto iT = clamp<int>(static_cast<int>(alphaT), 0, nTemp - 2);
1177 2228635317 alphaT -= iT;
1178
1179 2228635317 Scalar alphaP1 = policy.pressIdx(p, iT);
1180 2228635317 Scalar alphaP2 = policy.pressIdx(p, iT + 1);
1181
1182
2/2
✓ Branch 0 taken 2676832 times.
✓ Branch 1 taken 1199752360 times.
2228635317 const auto iP1 = clamp<int>(static_cast<int>(alphaP1), 0, nPress - 2);
1183
2/2
✓ Branch 0 taken 7868629 times.
✓ Branch 1 taken 1194560563 times.
2228635317 const auto iP2 = clamp<int>(static_cast<int>(alphaP2), 0, nPress - 2);
1184 2228635317 alphaP1 -= iP1;
1185 2228635317 alphaP2 -= iP2;
1186
1187 #if 0 && !defined NDEBUG
1188 const auto tempMin = policy.tempMin();
1189 const auto tempMin = policy.tempMax();
1190 if(!(0 <= alphaT && alphaT <= 1.0))
1191 DUNE_THROW(NumericalProblem, "Temperature out of range: "
1192 << "T=" << T << " range: [" << tempMin_ << ", " << tempMax_ << "]");
1193 if(!(0 <= alphaP1 && alphaP1 <= 1.0))
1194 DUNE_THROW(NumericalProblem, "First pressure out of range: "
1195 << "p=" << p << " range: [" << policy.minP(policy.tempIdx(T)) << ", " << policy.maxP(policy.tempIdx(T)) << "]");
1196 if(!(0 <= alphaP2 && alphaP2 <= 1.0))
1197 DUNE_THROW(NumericalProblem, "Second pressure out of range: "
1198 << "p=" << p << " range: [" << policy.minP(policy.tempIdx(T) + 1) << ", " << policy.maxP(policy.tempIdx(T) + 1) << "]");
1199 #endif
1200 2228635317 return values[(iT ) + (iP1 )*nTemp]*(1 - alphaT)*(1 - alphaP1) +
1201 2228635317 values[(iT ) + (iP1 + 1)*nTemp]*(1 - alphaT)*( alphaP1) +
1202 2228635317 values[(iT + 1) + (iP2 )*nTemp]*( alphaT)*(1 - alphaP2) +
1203 2228635317 values[(iT + 1) + (iP2 + 1)*nTemp]*( alphaT)*( alphaP2);
1204 }
1205
1206 //! returns an interpolated value for gas depending on temperature and density
1207 template<class Policy>
1208 2880812 static Scalar interpolateTRho_(const std::vector<Scalar>& values, const Scalar T, const Scalar rho, Policy policy)
1209 {
1210 2880812 const auto nTemp = policy.nTemp();
1211 2880812 const auto nDensity = policy.nDensity();
1212
1213 using std::clamp;
1214 2880812 Scalar alphaT = policy.tempIdx(T);
1215
2/4
✓ Branch 0 taken 1440406 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1440406 times.
✗ Branch 3 not taken.
2880812 if (alphaT < 0 - 1e-7*nTemp || alphaT >= nTemp - 1 + 1e-7*nTemp)
1216 return std::numeric_limits<Scalar>::quiet_NaN();
1217
1218
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1440406 times.
2880812 const auto iT = clamp<int>(static_cast<int>(alphaT), 0, nTemp - 2);
1219 2880812 alphaT -= iT;
1220
1221
2/2
✓ Branch 0 taken 6862 times.
✓ Branch 1 taken 1433544 times.
2880812 Scalar alphaP1 = policy.rhoIdx(rho, iT);
1222 2880812 Scalar alphaP2 = policy.rhoIdx(rho, iT + 1);
1223
2/2
✓ Branch 0 taken 6862 times.
✓ Branch 1 taken 1433544 times.
2880812 const auto iP1 = clamp<int>(static_cast<int>(alphaP1), 0, nDensity - 2);
1224
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1440406 times.
2880812 const auto iP2 = clamp<int>(static_cast<int>(alphaP2), 0, nDensity - 2);
1225 2880812 alphaP1 -= iP1;
1226 2880812 alphaP2 -= iP2;
1227
1228 2880812 return values[(iT ) + (iP1 )*nTemp]*(1 - alphaT)*(1 - alphaP1) +
1229 2880812 values[(iT ) + (iP1 + 1)*nTemp]*(1 - alphaT)*( alphaP1) +
1230 2880812 values[(iT + 1) + (iP2 )*nTemp]*( alphaT)*(1 - alphaP2) +
1231 2880812 values[(iT + 1) + (iP2 + 1)*nTemp]*( alphaT)*( alphaP2);
1232 }
1233
1234 // specifies whether the table was initialized
1235 static bool initialized_;
1236
1237 #ifndef NDEBUG
1238 // specifies whether some warning was printed
1239 static bool warningPrinted_;
1240 #endif
1241
1242 2509085751 static Table& table()
1243 {
1244
3/4
✓ Branch 0 taken 159 times.
✓ Branch 1 taken 2463796971 times.
✓ Branch 3 taken 159 times.
✗ Branch 4 not taken.
2509085751 static Table t;
1245 2509085751 return t;
1246 }
1247 };
1248
1249 template <class RawComponent, bool useVaporPressure>
1250 bool TabulatedComponent<RawComponent, useVaporPressure>::initialized_ = false;
1251
1252 #ifndef NDEBUG
1253 template <class RawComponent, bool useVaporPressure>
1254 bool TabulatedComponent<RawComponent, useVaporPressure>::warningPrinted_ = false;
1255 #endif
1256
1257 // forward declaration
1258 template <class Component>
1259 struct IsAqueous;
1260
1261 // we are aqueous if the raw compont is so
1262 template <class RawComponent, bool useVaporPressure>
1263 struct IsAqueous<TabulatedComponent<RawComponent, useVaporPressure>> : public IsAqueous<RawComponent> {};
1264
1265 } // end namespace Dumux::Components
1266
1267 #endif
1268