3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
tabulatedcomponent.hh
Go to the documentation of this file.
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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
28#ifndef DUMUX_TABULATED_COMPONENT_HH
29#define DUMUX_TABULATED_COMPONENT_HH
30
31#include <cmath>
32#include <limits>
33#include <cassert>
34#include <vector>
35#include <iostream>
36#include <iomanip>
37
40
41namespace Dumux {
42namespace Components {
43// forward declaration
44template<class RawComponent, bool useVaporPressure>
45class TabulatedComponent;
46} // end namespace Components
47
49template<class RawComponent, bool useVaporPressure>
50struct ComponentTraits<Components::TabulatedComponent<RawComponent, useVaporPressure>>
51{
52 using Scalar = typename RawComponent::Scalar;
53
55 static constexpr bool hasSolidState = std::is_base_of<Components::Solid<Scalar, RawComponent>, RawComponent>::value;
56
58 static constexpr bool hasLiquidState = std::is_base_of<Components::Liquid<Scalar, RawComponent>, RawComponent>::value;
59
61 static constexpr bool hasGasState = std::is_base_of<Components::Gas<Scalar, RawComponent>, RawComponent>::value;
62};
63
64namespace Components {
65
80template <class RawComponent, bool useVaporPressure=true>
82{
83public:
85 using Scalar = typename RawComponent::Scalar;
86
88 static constexpr bool isTabulated = true;
89
100 static void init(Scalar tempMin, Scalar tempMax, std::size_t nTemp,
101 Scalar pressMin, Scalar pressMax, std::size_t nPress)
102 {
103 tempMin_ = tempMin;
104 tempMax_ = tempMax;
105 nTemp_ = nTemp;
106 pressMin_ = pressMin;
107 pressMax_ = pressMax;
108 nPress_ = nPress;
109 nDensity_ = nPress_;
110
111 std::cout << "-------------------------------------------------------------------------\n"
112 << "Initializing tables for the " << RawComponent::name()
113 << " fluid properties (" << nTemp*nPress << " entries).\n"
114 << "Temperature -> min: " << std::scientific << std::setprecision(3)
115 << tempMin << ", max: " << tempMax << ", n: " << nTemp << '\n'
116 << "Pressure -> min: " << std::scientific << std::setprecision(3)
117 << pressMin << ", max: " << pressMax << ", n: " << nPress << '\n'
118 << "-------------------------------------------------------------------------" << std::endl;
119
120 // resize & initilialize the arrays with NaN
121 assert(std::numeric_limits<Scalar>::has_quiet_NaN);
122 const auto NaN = std::numeric_limits<Scalar>::quiet_NaN();
123
124 vaporPressure_.resize(nTemp_, NaN);
125 minGasDensity_.resize(nTemp_, NaN);
126 maxGasDensity_.resize(nTemp_, NaN);
127 minLiquidDensity_.resize(nTemp_, NaN);
128 maxLiquidDensity_.resize(nTemp_, NaN);
129
130 const std::size_t numEntriesTp = nTemp_*nPress_; // = nTemp_*nDensity_
131 gasEnthalpy_.resize(numEntriesTp, NaN);
132 liquidEnthalpy_.resize(numEntriesTp, NaN);
133 gasHeatCapacity_.resize(numEntriesTp, NaN);
134 liquidHeatCapacity_.resize(numEntriesTp, NaN);
135 gasDensity_.resize(numEntriesTp, NaN);
136 liquidDensity_.resize(numEntriesTp, NaN);
137 gasViscosity_.resize(numEntriesTp, NaN);
138 liquidViscosity_.resize(numEntriesTp, NaN);
139 gasThermalConductivity_.resize(numEntriesTp, NaN);
140 liquidThermalConductivity_.resize(numEntriesTp, NaN);
141 gasPressure_.resize(numEntriesTp, NaN);
142 liquidPressure_.resize(numEntriesTp, NaN);
143
144 // reset all flags
145 minMaxLiquidDensityInitialized_ = false;
146 minMaxGasDensityInitialized_ = false;
147 gasEnthalpyInitialized_ = false;
148 liquidEnthalpyInitialized_ = false;
149 gasHeatCapacityInitialized_ = false;
150 liquidHeatCapacityInitialized_ = false;
151 gasDensityInitialized_ = false;
152 liquidDensityInitialized_ = false;
153 gasViscosityInitialized_ = false;
154 liquidViscosityInitialized_ = false;
155 gasThermalConductivityInitialized_ = false;
156 liquidThermalConductivityInitialized_ = false;
157 gasPressureInitialized_ = false;
158 liquidPressureInitialized_ = false;
159
161 initVaporPressure_();
162
163#ifndef NDEBUG
164 initialized_ = true;
165 warningPrinted_ = false;
166#endif
167 }
168
172 static std::string name()
173 { return RawComponent::name(); }
174
178 static constexpr Scalar molarMass()
179 { return RawComponent::molarMass(); }
180
185 { return RawComponent::criticalTemperature(); }
186
191 { return RawComponent::criticalPressure(); }
192
197 { return RawComponent::tripleTemperature(); }
198
203 { return RawComponent::triplePressure(); }
204
212 {
213 using std::isnan;
214 Scalar result = interpolateT_(vaporPressure_, T);
215 if (isnan(result))
216 return RawComponent::vaporPressure(T);
217 return result;
218 }
219
228 {
229 return RawComponent::vaporTemperature(pressure);
230 }
231
239 {
240 Scalar result = interpolateTP_(gasEnthalpy_, temperature, pressure,
241 pressGasIdx_, minGasPressure_, maxGasPressure_, "gas");
242 using std::isnan;
243 if (isnan(result))
244 {
245 if (!gasEnthalpyInitialized_)
246 {
247 auto gasEnth = [] (auto T, auto p) { return RawComponent::gasEnthalpy(T, p); };
248 initTPArray_(gasEnth, minGasPressure_, maxGasPressure_, gasEnthalpy_);
249 gasEnthalpyInitialized_ = true;
251 }
252
253 printWarning_("gasEnthalpy", temperature, pressure);
254 return RawComponent::gasEnthalpy(temperature, pressure);
255 }
256 return result;
257 }
258
266 {
267 Scalar result = interpolateTP_(liquidEnthalpy_, temperature, pressure,
268 pressLiquidIdx_, minLiquidPressure_, maxLiquidPressure_, "liquid");
269 using std::isnan;
270 if (isnan(result))
271 {
272 if (!liquidEnthalpyInitialized_)
273 {
274 auto liqEnth = [] (auto T, auto p) { return RawComponent::liquidEnthalpy(T, p); };
275 initTPArray_(liqEnth, minLiquidPressure_, maxLiquidPressure_, liquidEnthalpy_);
276 liquidEnthalpyInitialized_ = true;
278 }
279
280 printWarning_("liquidEnthalpy", temperature, pressure);
281 return RawComponent::liquidEnthalpy(temperature, pressure);
282 }
283 return result;
284 }
285
293 {
294 Scalar result = interpolateTP_(gasHeatCapacity_, temperature, pressure,
295 pressGasIdx_, minGasPressure_, maxGasPressure_, "gas");
296 using std::isnan;
297 if (isnan(result))
298 {
299 if (!gasHeatCapacityInitialized_)
300 {
301 auto gasHC = [] (auto T, auto p) { return RawComponent::gasHeatCapacity(T, p); };
302 initTPArray_(gasHC, minGasPressure_, maxGasPressure_, gasHeatCapacity_);
303 gasHeatCapacityInitialized_ = true;
305 }
306
307 printWarning_("gasHeatCapacity", temperature, pressure);
308 return RawComponent::gasHeatCapacity(temperature, pressure);
309 }
310 return result;
311 }
312
320 {
321 Scalar result = interpolateTP_(liquidHeatCapacity_, temperature, pressure,
322 pressLiquidIdx_, minLiquidPressure_, maxLiquidPressure_, "liquid");
323 using std::isnan;
324 if (isnan(result))
325 {
326 if (!liquidHeatCapacityInitialized_)
327 {
328 auto liqHC = [] (auto T, auto p) { return RawComponent::liquidHeatCapacity(T, p); };
329 initTPArray_(liqHC, minLiquidPressure_, maxLiquidPressure_, liquidHeatCapacity_);
330 liquidHeatCapacityInitialized_ = true;
332 }
333
334 printWarning_("liquidHeatCapacity", temperature, pressure);
335 return RawComponent::liquidHeatCapacity(temperature, pressure);
336 }
337 return result;
338 }
339
347 {
349 }
350
358 {
360 }
361
369 {
371 if (!minMaxGasDensityInitialized_)
372 {
373 auto gasRho = [] (auto T, auto p) { return RawComponent::gasDensity(T, p); };
374 initMinMaxRhoArray_(gasRho, minGasPressure_, maxGasPressure_, minGasDensity_, maxGasDensity_);
375 minMaxGasDensityInitialized_ = true;
376 }
377
378 Scalar result = interpolateTRho_(gasPressure_, temperature, density, densityGasIdx_);
379 using std::isnan;
380 if (isnan(result))
381 {
382 if (!gasPressureInitialized_)
383 {
384 auto gasPFunc = [] (auto T, auto rho) { return RawComponent::gasPressure(T, rho); };
385 initPressureArray_(gasPressure_, gasPFunc, minGasDensity_, maxGasDensity_);
386 gasPressureInitialized_ = true;
388 }
389
390 printWarning_("gasPressure", temperature, density);
391 return RawComponent::gasPressure(temperature, density);
392 }
393 return result;
394 }
395
403 {
405 if (!minMaxLiquidDensityInitialized_)
406 {
407 auto liqRho = [] (auto T, auto p) { return RawComponent::liquidDensity(T, p); };
408 initMinMaxRhoArray_(liqRho, minLiquidPressure_, maxLiquidPressure_, minLiquidDensity_, maxLiquidDensity_);
409 minMaxLiquidDensityInitialized_ = true;
410 }
411
412 Scalar result = interpolateTRho_(liquidPressure_, temperature, density, densityLiquidIdx_);
413 using std::isnan;
414 if (isnan(result))
415 {
416 if (!liquidPressureInitialized_)
417 {
418 auto liqPFunc = [] (auto T, auto rho) { return RawComponent::liquidPressure(T, rho); };
419 initPressureArray_(liquidPressure_, liqPFunc, minLiquidDensity_, maxLiquidDensity_);
420 liquidPressureInitialized_ = true;
422 }
423
424 printWarning_("liquidPressure", temperature, density);
425 return RawComponent::liquidPressure(temperature, density);
426 }
427 return result;
428 }
429
433 static constexpr bool gasIsCompressible()
434 { return RawComponent::gasIsCompressible(); }
435
439 static constexpr bool liquidIsCompressible()
440 { return RawComponent::liquidIsCompressible(); }
441
445 static constexpr bool gasIsIdeal()
446 { return RawComponent::gasIsIdeal(); }
447
448
457 {
458 Scalar result = interpolateTP_(gasDensity_, temperature, pressure,
459 pressGasIdx_, minGasPressure_, maxGasPressure_, "gas");
460 using std::isnan;
461 if (isnan(result))
462 {
463 if (!gasDensityInitialized_)
464 {
465 auto gasRho = [] (auto T, auto p) { return RawComponent::gasDensity(T, p); };
466 initTPArray_(gasRho, minGasPressure_, maxGasPressure_, gasDensity_);
467 gasDensityInitialized_ = true;
469 }
470
471 printWarning_("gasDensity", temperature, pressure);
472 return RawComponent::gasDensity(temperature, pressure);
473 }
474 return result;
475 }
476
487
496 {
497 Scalar result = interpolateTP_(liquidDensity_, temperature, pressure,
498 pressLiquidIdx_, minLiquidPressure_, maxLiquidPressure_, "liquid");
499 using std::isnan;
500 if (isnan(result))
501 {
502 if (!liquidDensityInitialized_)
503 {
504 // TODO: we could get rid of the lambdas and pass the functor irectly. But,
505 // currently Brine is a component (and not a fluid system) expecting a
506 // third argument with a default, which cannot be wrapped in a function pointer.
507 // For this reason we have to wrap this into a lambda here.
508 auto liqRho = [] (auto T, auto p) { return RawComponent::liquidDensity(T, p); };
509 initTPArray_(liqRho, minLiquidPressure_, maxLiquidPressure_, liquidDensity_);
510 liquidDensityInitialized_ = true;
512 }
513
514 printWarning_("liquidDensity", temperature, pressure);
515 return RawComponent::liquidDensity(temperature, pressure);
516 }
517
518 return result;
519 }
520
531
539 {
540 Scalar result = interpolateTP_(gasViscosity_, temperature, pressure,
541 pressGasIdx_, minGasPressure_, maxGasPressure_, "gas");
542 using std::isnan;
543 if (isnan(result))
544 {
545 if (!gasViscosityInitialized_)
546 {
547 auto gasVisc = [] (auto T, auto p) { return RawComponent::gasViscosity(T, p); };
548 initTPArray_(gasVisc, minGasPressure_, maxGasPressure_, gasViscosity_);
549 gasViscosityInitialized_ = true;
551 }
552
553 printWarning_("gasViscosity", temperature, pressure);
554 return RawComponent::gasViscosity(temperature, pressure);
555 }
556 return result;
557 }
558
566 {
567 Scalar result = interpolateTP_(liquidViscosity_, temperature, pressure,
568 pressLiquidIdx_, minLiquidPressure_, maxLiquidPressure_, "liquid");
569 using std::isnan;
570 if (isnan(result))
571 {
572 if (!liquidViscosityInitialized_)
573 {
574 auto liqVisc = [] (auto T, auto p) { return RawComponent::liquidViscosity(T, p); };
575 initTPArray_(liqVisc, minLiquidPressure_, maxLiquidPressure_, liquidViscosity_);
576 liquidViscosityInitialized_ = true;
578 }
579
580 printWarning_("liquidViscosity",temperature, pressure);
581 return RawComponent::liquidViscosity(temperature, pressure);
582 }
583 return result;
584 }
585
593 {
594 Scalar result = interpolateTP_(gasThermalConductivity_, temperature, pressure,
595 pressGasIdx_, minGasPressure_, maxGasPressure_, "gas");
596 using std::isnan;
597 if (isnan(result))
598 {
599 if (!gasThermalConductivityInitialized_)
600 {
601 auto gasTC = [] (auto T, auto p) { return RawComponent::gasThermalConductivity(T, p); };
602 initTPArray_(gasTC, minGasPressure_, maxGasPressure_, gasThermalConductivity_);
603 gasThermalConductivityInitialized_ = true;
605 }
606
607 printWarning_("gasThermalConductivity", temperature, pressure);
608 return RawComponent::gasThermalConductivity(temperature, pressure);
609 }
610 return result;
611 }
612
620 {
621 Scalar result = interpolateTP_(liquidThermalConductivity_, temperature, pressure,
622 pressLiquidIdx_, minLiquidPressure_, maxLiquidPressure_, "liquid");
623 using std::isnan;
624 if (isnan(result))
625 {
626 if (!liquidThermalConductivityInitialized_)
627 {
628 auto liqTC = [] (auto T, auto p) { return RawComponent::liquidThermalConductivity(T, p); };
629 initTPArray_(liqTC, minLiquidPressure_, maxLiquidPressure_, liquidThermalConductivity_);
630 liquidThermalConductivityInitialized_ = true;
632 }
633
634 printWarning_("liquidThermalConductivity", temperature, pressure);
635 return RawComponent::liquidThermalConductivity(temperature, pressure);
636 }
637 return result;
638 }
639
640
641private:
643 static void printWarning_(const std::string& quantity, Scalar arg1, Scalar arg2)
644 {
645#ifndef NDEBUG
646 if (warningPrinted_)
647 return;
648
649 if (!initialized_)
650 std::cerr << "Warning: tabulated component '" << name()
651 << "' has not been initialized. "
652 << "Call FluidSystem::init() to use the tabulation in order to reduce runtime. \n";
653 else
654 std::cerr << "Warning: "<<quantity<<"(T="<<arg1<<", p="<<arg2<<") of component '"<<name()
655 << "' is outside tabulation range: ("<< tempMin_<<"<=T<="<<tempMax_<<"), ("
656 << pressMin_<<"<=p<=" <<pressMax_<<"). "
657 << "Forwarded to FluidSystem for direct evaluation of "<<quantity<<". \n";
658 warningPrinted_ = true;
659#endif
660 }
661
663 template< bool useVP = useVaporPressure, std::enable_if_t<useVP, int> = 0 >
664 static void initVaporPressure_()
665 {
666 // fill the temperature-pressure arrays
667 for (unsigned iT = 0; iT < nTemp_; ++ iT)
668 {
669 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
670 vaporPressure_[iT] = RawComponent::vaporPressure(temperature);
671 }
672 }
673
675 template< bool useVP = useVaporPressure, std::enable_if_t<!useVP, int> = 0 >
676 static void initVaporPressure_() {}
677
692 template<class PropFunc, class MinPFunc, class MaxPFunc>
693 static void initTPArray_(PropFunc&& f, MinPFunc&& minP, MaxPFunc&& maxP, std::vector<typename RawComponent::Scalar>& values)
694 {
695 for (unsigned iT = 0; iT < nTemp_; ++ iT)
696 {
697 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
698
699 Scalar pMax = maxP(iT);
700 Scalar pMin = minP(iT);
701 for (unsigned iP = 0; iP < nPress_; ++ iP)
702 {
703 Scalar pressure = iP * (pMax - pMin)/(nPress_ - 1) + pMin;
704 values[iT + iP*nTemp_] = f(temperature, pressure);
705 }
706 }
707 }
708
724 template<class RhoFunc, class MinPFunc, class MaxPFunc>
725 static void initMinMaxRhoArray_(RhoFunc&& rho,
726 MinPFunc&& minP,
727 MaxPFunc&& maxP,
728 std::vector<typename RawComponent::Scalar>& rhoMin,
729 std::vector<typename RawComponent::Scalar>& rhoMax)
730 {
731 for (unsigned iT = 0; iT < nTemp_; ++ iT)
732 {
733 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
734
735 rhoMin[iT] = rho(temperature, minP(iT));
736 if (iT < nTemp_ - 1)
737 rhoMax[iT] = rho(temperature, maxP(iT + 1));
738 else
739 rhoMax[iT] = rho(temperature, maxP(iT));
740 }
741 }
742
753 template<class PFunc>
754 static void initPressureArray_(std::vector<typename RawComponent::Scalar>& pressure, PFunc&& p,
755 const std::vector<typename RawComponent::Scalar>& rhoMin,
756 const std::vector<typename RawComponent::Scalar>& rhoMax)
757 {
758 for (unsigned iT = 0; iT < nTemp_; ++ iT)
759 {
760 Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
761
762 for (unsigned iRho = 0; iRho < nDensity_; ++ iRho)
763 {
764 Scalar density = Scalar(iRho)/(nDensity_ - 1)
765 * (rhoMax[iT] - rhoMin[iT])
766 + rhoMin[iT];
767 pressure[iT + iRho*nTemp_] = p(temperature, density);
768 }
769 }
770 }
771
773 static Scalar interpolateT_(const std::vector<typename RawComponent::Scalar>& values, Scalar T)
774 {
775 Scalar alphaT = tempIdx_(T);
776 if (alphaT < 0 || alphaT >= nTemp_ - 1)
777 return std::numeric_limits<Scalar>::quiet_NaN();
778
779 unsigned iT = (unsigned) alphaT;
780 alphaT -= iT;
781
782 return values[iT ]*(1 - alphaT) +
783 values[iT + 1]*( alphaT);
784 }
785
787 template<class GetPIdx, class MinPFunc, class MaxPFunc>
788 static Scalar interpolateTP_(const std::vector<typename RawComponent::Scalar>& values, Scalar T, Scalar p,
789 GetPIdx&& getPIdx, MinPFunc&& minP, MaxPFunc&& maxP,
790 const std::string& phaseName)
791 {
792 Scalar alphaT = tempIdx_(T);
793 if (alphaT < 0 || alphaT >= nTemp_ - 1) {
794 return std::numeric_limits<Scalar>::quiet_NaN();
795 }
796 using std::min;
797 using std::max;
798 unsigned iT = max<int>(0, min<int>(nTemp_ - 2, (int) alphaT));
799 alphaT -= iT;
800
801 Scalar alphaP1 = getPIdx(p, iT);
802 Scalar alphaP2 = getPIdx(p, iT + 1);
803
804 unsigned iP1 = max<int>(0, min<int>(nPress_ - 2, (int) alphaP1));
805 unsigned iP2 = max<int>(0, min<int>(nPress_ - 2, (int) alphaP2));
806 alphaP1 -= iP1;
807 alphaP2 -= iP2;
808
809#if 0 && !defined NDEBUG
810 if(!(0 <= alphaT && alphaT <= 1.0))
811 DUNE_THROW(NumericalProblem, "Temperature out of range: "
812 << "T=" << T << " range: [" << tempMin_ << ", " << tempMax_ << "]");
813 if(!(0 <= alphaP1 && alphaP1 <= 1.0))
814 DUNE_THROW(NumericalProblem, "First " << phaseName " pressure out of range: "
815 << "p=" << p << " range: [" << minP(tempIdx_(T)) << ", " << maxP(tempIdx_(T)) << "]");
816 if(!(0 <= alphaP2 && alphaP2 <= 1.0))
817 DUNE_THROW(NumericalProblem, "Second " << phaseName " pressure out of range: "
818 << "p=" << p << " range: [" << minP(tempIdx_(T) + 1) << ", " << maxP(tempIdx_(T) + 1) << "]");
819#endif
820
821 return values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
822 values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
823 values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
824 values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
825 }
826
828 template<class GetRhoIdx>
829 static Scalar interpolateTRho_(const std::vector<typename RawComponent::Scalar>& values, Scalar T, Scalar rho, GetRhoIdx&& rhoIdx)
830 {
831 using std::min;
832 using std::max;
833 Scalar alphaT = tempIdx_(T);
834 unsigned iT = max<int>(0, min<int>(nTemp_ - 2, (int) alphaT));
835 alphaT -= iT;
836
837 Scalar alphaP1 = rhoIdx(rho, iT);
838 Scalar alphaP2 = rhoIdx(rho, iT + 1);
839 unsigned iP1 = max<int>(0, min<int>(nDensity_ - 2, (int) alphaP1));
840 unsigned iP2 = max<int>(0, min<int>(nDensity_ - 2, (int) alphaP2));
841 alphaP1 -= iP1;
842 alphaP2 -= iP2;
843
844 return values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
845 values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
846 values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
847 values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
848 }
849
851 static Scalar tempIdx_(Scalar temperature)
852 {
853 return (nTemp_ - 1)*(temperature - tempMin_)/(tempMax_ - tempMin_);
854 }
855
857 static Scalar pressLiquidIdx_(Scalar pressure, unsigned tempIdx)
858 {
859 Scalar plMin = minLiquidPressure_(tempIdx);
860 Scalar plMax = maxLiquidPressure_(tempIdx);
861
862 return (nPress_ - 1)*(pressure - plMin)/(plMax - plMin);
863 }
864
866 static Scalar pressGasIdx_(Scalar pressure, unsigned tempIdx)
867 {
868 Scalar pgMin = minGasPressure_(tempIdx);
869 Scalar pgMax = maxGasPressure_(tempIdx);
870
871 return (nPress_ - 1)*(pressure - pgMin)/(pgMax - pgMin);
872 }
873
875 static Scalar densityLiquidIdx_(Scalar density, unsigned tempIdx)
876 {
877 Scalar densityMin = minLiquidDensity_[tempIdx];
878 Scalar densityMax = maxLiquidDensity_[tempIdx];
879 return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
880 }
881
883 static Scalar densityGasIdx_(Scalar density, unsigned tempIdx)
884 {
885 Scalar densityMin = minGasDensity_[tempIdx];
886 Scalar densityMax = maxGasDensity_[tempIdx];
887 return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
888 }
889
891 static Scalar minLiquidPressure_(int tempIdx)
892 {
893 using std::max;
894 if (!useVaporPressure)
895 return pressMin_;
896 else
897 return max(pressMin_, vaporPressure_[tempIdx] / 1.1);
898 }
899
901 static Scalar maxLiquidPressure_(int tempIdx)
902 {
903 using std::max;
904 if (!useVaporPressure)
905 return pressMax_;
906 else
907 return max(pressMax_, vaporPressure_[tempIdx] * 1.1);
908 }
909
911 static Scalar minGasPressure_(int tempIdx)
912 {
913 using std::min;
914 if (!useVaporPressure)
915 return pressMin_;
916 else
917 return min(pressMin_, vaporPressure_[tempIdx] / 1.1 );
918 }
919
921 static Scalar maxGasPressure_(int tempIdx)
922 {
923 using std::min;
924 if (!useVaporPressure)
925 return pressMax_;
926 else
927 return min(pressMax_, vaporPressure_[tempIdx] * 1.1);
928 }
929
930#ifndef NDEBUG
931 // specifies whether the table was initialized
932 static bool initialized_;
933 // specifies whether some warning was printed
934 static bool warningPrinted_;
935#endif
936
937 // 1D fields with the temperature as degree of freedom
938 static std::vector<typename RawComponent::Scalar> vaporPressure_;
939
940 static std::vector<typename RawComponent::Scalar> minLiquidDensity_;
941 static std::vector<typename RawComponent::Scalar> maxLiquidDensity_;
942 static bool minMaxLiquidDensityInitialized_;
943
944 static std::vector<typename RawComponent::Scalar> minGasDensity_;
945 static std::vector<typename RawComponent::Scalar> maxGasDensity_;
946 static bool minMaxGasDensityInitialized_;
947
948 // 2D fields with the temperature and pressure as degrees of freedom
949 static std::vector<typename RawComponent::Scalar> gasEnthalpy_;
950 static std::vector<typename RawComponent::Scalar> liquidEnthalpy_;
951 static bool gasEnthalpyInitialized_;
952 static bool liquidEnthalpyInitialized_;
953
954 static std::vector<typename RawComponent::Scalar> gasHeatCapacity_;
955 static std::vector<typename RawComponent::Scalar> liquidHeatCapacity_;
956 static bool gasHeatCapacityInitialized_;
957 static bool liquidHeatCapacityInitialized_;
958
959 static std::vector<typename RawComponent::Scalar> gasDensity_;
960 static std::vector<typename RawComponent::Scalar> liquidDensity_;
961 static bool gasDensityInitialized_;
962 static bool liquidDensityInitialized_;
963
964 static std::vector<typename RawComponent::Scalar> gasViscosity_;
965 static std::vector<typename RawComponent::Scalar> liquidViscosity_;
966 static bool gasViscosityInitialized_;
967 static bool liquidViscosityInitialized_;
968
969 static std::vector<typename RawComponent::Scalar> gasThermalConductivity_;
970 static std::vector<typename RawComponent::Scalar> liquidThermalConductivity_;
971 static bool gasThermalConductivityInitialized_;
972 static bool liquidThermalConductivityInitialized_;
973
974 // 2D fields with the temperature and density as degrees of freedom
975 static std::vector<typename RawComponent::Scalar> gasPressure_;
976 static std::vector<typename RawComponent::Scalar> liquidPressure_;
977 static bool gasPressureInitialized_;
978 static bool liquidPressureInitialized_;
979
980 // temperature, pressure and density ranges
981 static Scalar tempMin_;
982 static Scalar tempMax_;
983 static unsigned nTemp_;
984
985 static Scalar pressMin_;
986 static Scalar pressMax_;
987 static unsigned nPress_;
988
989 static Scalar densityMin_;
990 static Scalar densityMax_;
991 static unsigned nDensity_;
992};
993
994#ifndef NDEBUG
995template <class RawComponent, bool useVaporPressure>
996bool TabulatedComponent<RawComponent, useVaporPressure>::initialized_ = false;
997
998template <class RawComponent, bool useVaporPressure>
999bool TabulatedComponent<RawComponent, useVaporPressure>::warningPrinted_ = false;
1000#endif
1001
1002template <class RawComponent, bool useVaporPressure>
1003bool TabulatedComponent<RawComponent, useVaporPressure>::minMaxLiquidDensityInitialized_ = false;
1004template <class RawComponent, bool useVaporPressure>
1005bool TabulatedComponent<RawComponent, useVaporPressure>::minMaxGasDensityInitialized_ = false;
1006template <class RawComponent, bool useVaporPressure>
1007bool TabulatedComponent<RawComponent, useVaporPressure>::gasEnthalpyInitialized_ = false;
1008template <class RawComponent, bool useVaporPressure>
1009bool TabulatedComponent<RawComponent, useVaporPressure>::liquidEnthalpyInitialized_ = false;
1010template <class RawComponent, bool useVaporPressure>
1011bool TabulatedComponent<RawComponent, useVaporPressure>::gasHeatCapacityInitialized_ = false;
1012template <class RawComponent, bool useVaporPressure>
1013bool TabulatedComponent<RawComponent, useVaporPressure>::liquidHeatCapacityInitialized_ = false;
1014template <class RawComponent, bool useVaporPressure>
1015bool TabulatedComponent<RawComponent, useVaporPressure>::gasDensityInitialized_ = false;
1016template <class RawComponent, bool useVaporPressure>
1017bool TabulatedComponent<RawComponent, useVaporPressure>::liquidDensityInitialized_ = false;
1018template <class RawComponent, bool useVaporPressure>
1019bool TabulatedComponent<RawComponent, useVaporPressure>::gasViscosityInitialized_ = false;
1020template <class RawComponent, bool useVaporPressure>
1021bool TabulatedComponent<RawComponent, useVaporPressure>::liquidViscosityInitialized_ = false;
1022template <class RawComponent, bool useVaporPressure>
1023bool TabulatedComponent<RawComponent, useVaporPressure>::gasThermalConductivityInitialized_ = false;
1024template <class RawComponent, bool useVaporPressure>
1025bool TabulatedComponent<RawComponent, useVaporPressure>::liquidThermalConductivityInitialized_ = false;
1026template <class RawComponent, bool useVaporPressure>
1027bool TabulatedComponent<RawComponent, useVaporPressure>::gasPressureInitialized_ = false;
1028template <class RawComponent, bool useVaporPressure>
1029bool TabulatedComponent<RawComponent, useVaporPressure>::liquidPressureInitialized_ = false;
1030
1031template <class RawComponent, bool useVaporPressure>
1032std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::vaporPressure_;
1033template <class RawComponent, bool useVaporPressure>
1034std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::minLiquidDensity_;
1035template <class RawComponent, bool useVaporPressure>
1036std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::maxLiquidDensity_;
1037template <class RawComponent, bool useVaporPressure>
1038std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::minGasDensity_;
1039template <class RawComponent, bool useVaporPressure>
1040std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::maxGasDensity_;
1041template <class RawComponent, bool useVaporPressure>
1042std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::gasEnthalpy_;
1043template <class RawComponent, bool useVaporPressure>
1044std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::liquidEnthalpy_;
1045template <class RawComponent, bool useVaporPressure>
1046std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::gasHeatCapacity_;
1047template <class RawComponent, bool useVaporPressure>
1048std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::liquidHeatCapacity_;
1049template <class RawComponent, bool useVaporPressure>
1050std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::gasDensity_;
1051template <class RawComponent, bool useVaporPressure>
1052std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::liquidDensity_;
1053template <class RawComponent, bool useVaporPressure>
1054std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::gasViscosity_;
1055template <class RawComponent, bool useVaporPressure>
1056std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::liquidViscosity_;
1057template <class RawComponent, bool useVaporPressure>
1058std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::gasThermalConductivity_;
1059template <class RawComponent, bool useVaporPressure>
1060std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::liquidThermalConductivity_;
1061template <class RawComponent, bool useVaporPressure>
1062std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::gasPressure_;
1063template <class RawComponent, bool useVaporPressure>
1064std::vector<typename RawComponent::Scalar> TabulatedComponent<RawComponent, useVaporPressure>::liquidPressure_;
1065template <class RawComponent, bool useVaporPressure>
1066typename RawComponent::Scalar TabulatedComponent<RawComponent, useVaporPressure>::tempMin_;
1067template <class RawComponent, bool useVaporPressure>
1068typename RawComponent::Scalar TabulatedComponent<RawComponent, useVaporPressure>::tempMax_;
1069template <class RawComponent, bool useVaporPressure>
1070unsigned TabulatedComponent<RawComponent, useVaporPressure>::nTemp_;
1071template <class RawComponent, bool useVaporPressure>
1072typename RawComponent::Scalar TabulatedComponent<RawComponent, useVaporPressure>::pressMin_;
1073template <class RawComponent, bool useVaporPressure>
1074typename RawComponent::Scalar TabulatedComponent<RawComponent, useVaporPressure>::pressMax_;
1075template <class RawComponent, bool useVaporPressure>
1076unsigned TabulatedComponent<RawComponent, useVaporPressure>::nPress_;
1077template <class RawComponent, bool useVaporPressure>
1078typename RawComponent::Scalar TabulatedComponent<RawComponent, useVaporPressure>::densityMin_;
1079template <class RawComponent, bool useVaporPressure>
1080typename RawComponent::Scalar TabulatedComponent<RawComponent, useVaporPressure>::densityMax_;
1081template <class RawComponent, bool useVaporPressure>
1082unsigned TabulatedComponent<RawComponent, useVaporPressure>::nDensity_;
1083
1084// forward declaration
1085template <class Component>
1086struct IsAqueous;
1087
1088// we are aqueous if the raw compont is so
1089template <class RawComponent, bool useVaporPressure>
1090struct IsAqueous<TabulatedComponent<RawComponent, useVaporPressure>> : public IsAqueous<RawComponent> {};
1091
1092} // end namespace Components
1093
1094} // end namespace Dumux
1095
1096#endif
Some exceptions thrown in DuMux
Component traits, i.e. information extracted from components.
Definition: adapt.hh:29
std::string temperature() noexcept
I/O name of temperature for equilibrium models.
Definition: name.hh:51
std::string pressure(int phaseIdx) noexcept
I/O name of pressure for multiphase systems.
Definition: name.hh:34
std::string density(int phaseIdx) noexcept
I/O name of density for multiphase systems.
Definition: name.hh:65
IsAqueous struct.
Definition: components/base.hh:47
Component traits, i.e. information extracted from components.
Definition: componenttraits.hh:43
static constexpr bool hasLiquidState
if the component implements a liquid state
Definition: componenttraits.hh:50
static constexpr bool hasSolidState
if the component implements a solid state
Definition: componenttraits.hh:47
static constexpr bool hasGasState
if the component implements a gaseous state
Definition: componenttraits.hh:53
Tabulates all thermodynamic properties of a given untabulated chemical species.
Definition: tabulatedcomponent.hh:82
static const Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of the gas .
Definition: tabulatedcomponent.hh:238
static Scalar gasPressure(Scalar temperature, Scalar density)
The pressure of gas in at a given density and temperature.
Definition: tabulatedcomponent.hh:368
static Scalar criticalTemperature()
Returns the critical temperature in of the component.
Definition: tabulatedcomponent.hh:184
static Scalar liquidPressure(Scalar temperature, Scalar density)
The pressure of liquid in at a given density and temperature.
Definition: tabulatedcomponent.hh:402
static const Scalar gasHeatCapacity(Scalar temperature, Scalar pressure)
Specific isobaric heat capacity of the gas .
Definition: tabulatedcomponent.hh:292
static std::string name()
A human readable name for the component.
Definition: tabulatedcomponent.hh:172
static Scalar vaporTemperature(Scalar pressure)
The vapor pressure in of the component at a given temperature.
Definition: tabulatedcomponent.hh:227
static Scalar liquidThermalConductivity(Scalar temperature, Scalar pressure)
The thermal conductivity of liquid water .
Definition: tabulatedcomponent.hh:619
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of gas at a given pressure and temperature .
Definition: tabulatedcomponent.hh:456
static const Scalar liquidEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of the liquid .
Definition: tabulatedcomponent.hh:265
static Scalar triplePressure()
Returns the pressure in at the component's triple point.
Definition: tabulatedcomponent.hh:202
static Scalar criticalPressure()
Returns the critical pressure in of the component.
Definition: tabulatedcomponent.hh:190
static const Scalar gasInternalEnergy(Scalar temperature, Scalar pressure)
Specific internal energy of the gas .
Definition: tabulatedcomponent.hh:346
static Scalar liquidMolarDensity(Scalar temperature, Scalar pressure)
The molar density of liquid in at a given pressure and temperature.
Definition: tabulatedcomponent.hh:529
static constexpr Scalar molarMass()
The molar mass in of the component.
Definition: tabulatedcomponent.hh:178
static Scalar liquidViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of liquid.
Definition: tabulatedcomponent.hh:565
static Scalar gasViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of gas.
Definition: tabulatedcomponent.hh:538
static Scalar tripleTemperature()
Returns the temperature in at the component's triple point.
Definition: tabulatedcomponent.hh:196
static constexpr bool liquidIsCompressible()
Returns true if the liquid phase is assumed to be compressible.
Definition: tabulatedcomponent.hh:439
static Scalar liquidDensity(Scalar temperature, Scalar pressure)
The density of liquid at a given pressure and temperature .
Definition: tabulatedcomponent.hh:495
static Scalar gasThermalConductivity(Scalar temperature, Scalar pressure)
The thermal conductivity of gaseous water .
Definition: tabulatedcomponent.hh:592
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar density of gas in at a given pressure and temperature.
Definition: tabulatedcomponent.hh:485
static constexpr bool isTabulated
state that we are tabulated
Definition: tabulatedcomponent.hh:88
static const Scalar liquidInternalEnergy(Scalar temperature, Scalar pressure)
Specific internal energy of the liquid .
Definition: tabulatedcomponent.hh:357
static Scalar vaporPressure(Scalar T)
The vapor pressure in of the component at a given temperature.
Definition: tabulatedcomponent.hh:211
static constexpr bool gasIsCompressible()
Returns true if the gas phase is assumed to be compressible.
Definition: tabulatedcomponent.hh:433
static void init(Scalar tempMin, Scalar tempMax, std::size_t nTemp, Scalar pressMin, Scalar pressMax, std::size_t nPress)
Initialize the tables.
Definition: tabulatedcomponent.hh:100
static constexpr bool gasIsIdeal()
Returns true if the gas phase is assumed to be ideal.
Definition: tabulatedcomponent.hh:445
static const Scalar liquidHeatCapacity(Scalar temperature, Scalar pressure)
Specific isobaric heat capacity of the liquid .
Definition: tabulatedcomponent.hh:319
typename RawComponent::Scalar Scalar
export scalar type
Definition: tabulatedcomponent.hh:85
typename RawComponent::Scalar Scalar
Definition: tabulatedcomponent.hh:52