All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
TabulatedComponent.hpp
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  This file is part of the Open Porous Media project (OPM).
5 
6  OPM 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 2 of the License, or
9  (at your option) any later version.
10 
11  OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
28 #ifndef OPM_TABULATED_COMPONENT_HPP
29 #define OPM_TABULATED_COMPONENT_HPP
30 
31 #include <cmath>
32 #include <limits>
33 #include <cassert>
34 #include <iostream>
35 
36 #include <opm/common/Exceptions.hpp>
37 #include <opm/common/ErrorMacros.hpp>
38 
40 
41 namespace Opm {
57 template <class ScalarT, class RawComponent, bool useVaporPressure=true>
59 {
60 public:
61  typedef ScalarT Scalar;
62 
63  static const bool isTabulated = true;
64 
75  static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
76  Scalar pressMin, Scalar pressMax, unsigned nPress)
77  {
78  tempMin_ = tempMin;
79  tempMax_ = tempMax;
80  nTemp_ = nTemp;
81  pressMin_ = pressMin;
82  pressMax_ = pressMax;
83  nPress_ = nPress;
84  nDensity_ = nPress_;
85 
86  // allocate the arrays
87  vaporPressure_ = new Scalar[nTemp_];
88  minGasDensity__ = new Scalar[nTemp_];
89  maxGasDensity__ = new Scalar[nTemp_];
90  minLiquidDensity__ = new Scalar[nTemp_];
91  maxLiquidDensity__ = new Scalar[nTemp_];
92 
93  gasEnthalpy_ = new Scalar[nTemp_*nPress_];
94  liquidEnthalpy_ = new Scalar[nTemp_*nPress_];
95  gasHeatCapacity_ = new Scalar[nTemp_*nPress_];
96  liquidHeatCapacity_ = new Scalar[nTemp_*nPress_];
97  gasDensity_ = new Scalar[nTemp_*nPress_];
98  liquidDensity_ = new Scalar[nTemp_*nPress_];
99  gasViscosity_ = new Scalar[nTemp_*nPress_];
100  liquidViscosity_ = new Scalar[nTemp_*nPress_];
101  gasThermalConductivity_ = new Scalar[nTemp_*nPress_];
102  liquidThermalConductivity_ = new Scalar[nTemp_*nPress_];
103  gasPressure_ = new Scalar[nTemp_*nDensity_];
104  liquidPressure_ = new Scalar[nTemp_*nDensity_];
105 
106  assert(std::numeric_limits<Scalar>::has_quiet_NaN);
107  Scalar NaN = std::numeric_limits<Scalar>::quiet_NaN();
108 
109  // fill the temperature-pressure arrays
110  for (unsigned iT = 0; iT < nTemp_; ++ iT) {
111  Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
112 
113  try { vaporPressure_[iT] = RawComponent::vaporPressure(temperature); }
114  catch (const std::exception&) { vaporPressure_[iT] = NaN; }
115 
116  Scalar pgMax = maxGasPressure_(iT);
117  Scalar pgMin = minGasPressure_(iT);
118 
119  // fill the temperature, pressure gas arrays
120  for (unsigned iP = 0; iP < nPress_; ++ iP) {
121  Scalar pressure = iP * (pgMax - pgMin)/(nPress_ - 1) + pgMin;
122 
123  unsigned i = iT + iP*nTemp_;
124 
125  try { gasEnthalpy_[i] = RawComponent::gasEnthalpy(temperature, pressure); }
126  catch (const std::exception&) { gasEnthalpy_[i] = NaN; }
127 
128  try { gasHeatCapacity_[i] = RawComponent::gasHeatCapacity(temperature, pressure); }
129  catch (const std::exception&) { gasHeatCapacity_[i] = NaN; }
130 
131  try { gasDensity_[i] = RawComponent::gasDensity(temperature, pressure); }
132  catch (const std::exception&) { gasDensity_[i] = NaN; }
133 
134  try { gasViscosity_[i] = RawComponent::gasViscosity(temperature, pressure); }
135  catch (const std::exception&) { gasViscosity_[i] = NaN; }
136 
137  try { gasThermalConductivity_[i] = RawComponent::gasThermalConductivity(temperature, pressure); }
138  catch (const std::exception&) { gasThermalConductivity_[i] = NaN; }
139  };
140 
141  Scalar plMin = minLiquidPressure_(iT);
142  Scalar plMax = maxLiquidPressure_(iT);
143  for (unsigned iP = 0; iP < nPress_; ++ iP) {
144  Scalar pressure = iP * (plMax - plMin)/(nPress_ - 1) + plMin;
145 
146  unsigned i = iT + iP*nTemp_;
147 
148  try { liquidEnthalpy_[i] = RawComponent::liquidEnthalpy(temperature, pressure); }
149  catch (const std::exception&) { liquidEnthalpy_[i] = NaN; }
150 
151  try { liquidHeatCapacity_[i] = RawComponent::liquidHeatCapacity(temperature, pressure); }
152  catch (const std::exception&) { liquidHeatCapacity_[i] = NaN; }
153 
154  try { liquidDensity_[i] = RawComponent::liquidDensity(temperature, pressure); }
155  catch (const std::exception&) { liquidDensity_[i] = NaN; }
156 
157  try { liquidViscosity_[i] = RawComponent::liquidViscosity(temperature, pressure); }
158  catch (const std::exception&) { liquidViscosity_[i] = NaN; }
159 
160  try { liquidThermalConductivity_[i] = RawComponent::liquidThermalConductivity(temperature, pressure); }
161  catch (const std::exception&) { liquidThermalConductivity_[i] = NaN; }
162  }
163  }
164 
165  // fill the temperature-density arrays
166  for (unsigned iT = 0; iT < nTemp_; ++ iT) {
167  Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
168 
169  // calculate the minimum and maximum values for the gas
170  // densities
171  minGasDensity__[iT] = RawComponent::gasDensity(temperature, minGasPressure_(iT));
172  if (iT < nTemp_ - 1)
173  maxGasDensity__[iT] = RawComponent::gasDensity(temperature, maxGasPressure_(iT + 1));
174  else
175  maxGasDensity__[iT] = RawComponent::gasDensity(temperature, maxGasPressure_(iT));
176 
177  // fill the temperature, density gas arrays
178  for (unsigned iRho = 0; iRho < nDensity_; ++ iRho) {
179  Scalar density =
180  Scalar(iRho)/(nDensity_ - 1) *
181  (maxGasDensity__[iT] - minGasDensity__[iT])
182  +
183  minGasDensity__[iT];
184 
185  unsigned i = iT + iRho*nTemp_;
186 
187  try { gasPressure_[i] = RawComponent::gasPressure(temperature, density); }
188  catch (const std::exception&) { gasPressure_[i] = NaN; };
189  };
190 
191  // calculate the minimum and maximum values for the liquid
192  // densities
193  minLiquidDensity__[iT] = RawComponent::liquidDensity(temperature, minLiquidPressure_(iT));
194  if (iT < nTemp_ - 1)
195  maxLiquidDensity__[iT] = RawComponent::liquidDensity(temperature, maxLiquidPressure_(iT + 1));
196  else
197  maxLiquidDensity__[iT] = RawComponent::liquidDensity(temperature, maxLiquidPressure_(iT));
198 
199  // fill the temperature, density liquid arrays
200  for (unsigned iRho = 0; iRho < nDensity_; ++ iRho) {
201  Scalar density =
202  Scalar(iRho)/(nDensity_ - 1) *
203  (maxLiquidDensity__[iT] - minLiquidDensity__[iT])
204  +
205  minLiquidDensity__[iT];
206 
207  unsigned i = iT + iRho*nTemp_;
208 
209  try { liquidPressure_[i] = RawComponent::liquidPressure(temperature, density); }
210  catch (const std::exception&) { liquidPressure_[i] = NaN; };
211  };
212  }
213  }
214 
218  static const char* name()
219  { return RawComponent::name(); }
220 
224  static Scalar molarMass()
225  { return RawComponent::molarMass(); }
226 
230  static Scalar criticalTemperature()
231  { return RawComponent::criticalTemperature(); }
232 
236  static Scalar criticalPressure()
237  { return RawComponent::criticalPressure(); }
238 
242  static Scalar tripleTemperature()
243  { return RawComponent::tripleTemperature(); }
244 
248  static Scalar triplePressure()
249  { return RawComponent::triplePressure(); }
250 
257  template <class Evaluation>
258  static Evaluation vaporPressure(const Evaluation& temperature)
259  {
260  const Evaluation& result = interpolateT_(vaporPressure_, temperature);
261  if (std::isnan(Opm::scalarValue(result)))
262  return RawComponent::vaporPressure(temperature);
263  return result;
264  }
265 
272  template <class Evaluation>
273  static Evaluation gasEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
274  {
275  const Evaluation& result = interpolateGasTP_(gasEnthalpy_,
276  temperature,
277  pressure);
278  if (std::isnan(Opm::scalarValue(result)))
279  return RawComponent::gasEnthalpy(temperature, pressure);
280  return result;
281  }
282 
289  template <class Evaluation>
290  static Evaluation liquidEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
291  {
292  const Evaluation& result = interpolateLiquidTP_(liquidEnthalpy_,
293  temperature,
294  pressure);
295  if (std::isnan(Opm::scalarValue(result)))
296  return RawComponent::liquidEnthalpy(temperature, pressure);
297  return result;
298  }
299 
306  template <class Evaluation>
307  static Evaluation gasHeatCapacity(const Evaluation& temperature, const Evaluation& pressure)
308  {
309  const Evaluation& result = interpolateGasTP_(gasHeatCapacity_,
310  temperature,
311  pressure);
312  if (std::isnan(Opm::scalarValue(result)))
313  return RawComponent::gasHeatCapacity(temperature, pressure);
314  return result;
315  }
316 
323  template <class Evaluation>
324  static Evaluation liquidHeatCapacity(const Evaluation& temperature, const Evaluation& pressure)
325  {
326  const Evaluation& result = interpolateLiquidTP_(liquidHeatCapacity_,
327  temperature,
328  pressure);
329  if (std::isnan(Opm::scalarValue(result)))
330  return RawComponent::liquidHeatCapacity(temperature, pressure);
331  return result;
332  }
333 
340  template <class Evaluation>
341  static Evaluation gasInternalEnergy(const Evaluation& temperature, const Evaluation& pressure)
342  { return gasEnthalpy(temperature, pressure) - pressure/gasDensity(temperature, pressure); }
343 
350  template <class Evaluation>
351  static Evaluation liquidInternalEnergy(const Evaluation& temperature, const Evaluation& pressure)
352  { return liquidEnthalpy(temperature, pressure) - pressure/liquidDensity(temperature, pressure); }
353 
360  template <class Evaluation>
361  static Evaluation gasPressure(const Evaluation& temperature, Scalar density)
362  {
363  const Evaluation& result = interpolateGasTRho_(gasPressure_,
364  temperature,
365  density);
366  if (std::isnan(Opm::scalarValue(result)))
367  return RawComponent::gasPressure(temperature,
368  density);
369  return result;
370  }
371 
378  template <class Evaluation>
379  static Evaluation liquidPressure(const Evaluation& temperature, Scalar density)
380  {
381  const Evaluation& result = interpolateLiquidTRho_(liquidPressure_,
382  temperature,
383  density);
384  if (std::isnan(Opm::scalarValue(result)))
385  return RawComponent::liquidPressure(temperature,
386  density);
387  return result;
388  }
389 
393  static bool gasIsCompressible()
394  { return RawComponent::gasIsCompressible(); }
395 
399  static bool liquidIsCompressible()
400  { return RawComponent::liquidIsCompressible(); }
401 
405  static bool gasIsIdeal()
406  { return RawComponent::gasIsIdeal(); }
407 
408 
416  template <class Evaluation>
417  static Evaluation gasDensity(const Evaluation& temperature, const Evaluation& pressure)
418  {
419  const Evaluation& result = interpolateGasTP_(gasDensity_,
420  temperature,
421  pressure);
422  if (std::isnan(Opm::scalarValue(result)))
423  return RawComponent::gasDensity(temperature, pressure);
424  return result;
425  }
426 
434  template <class Evaluation>
435  static Evaluation liquidDensity(const Evaluation& temperature, const Evaluation& pressure)
436  {
437  const Evaluation& result = interpolateLiquidTP_(liquidDensity_,
438  temperature,
439  pressure);
440  if (std::isnan(Opm::scalarValue(result)))
441  return RawComponent::liquidDensity(temperature, pressure);
442  return result;
443  }
444 
451  template <class Evaluation>
452  static Evaluation gasViscosity(const Evaluation& temperature, const Evaluation& pressure)
453  {
454  const Evaluation& result = interpolateGasTP_(gasViscosity_,
455  temperature,
456  pressure);
457  if (std::isnan(Opm::scalarValue(result)))
458  return RawComponent::gasViscosity(temperature, pressure);
459  return result;
460  }
461 
468  template <class Evaluation>
469  static Evaluation liquidViscosity(const Evaluation& temperature, const Evaluation& pressure)
470  {
471  const Evaluation& result = interpolateLiquidTP_(liquidViscosity_,
472  temperature,
473  pressure);
474  if (std::isnan(Opm::scalarValue(result)))
475  return RawComponent::liquidViscosity(temperature, pressure);
476  return result;
477  }
478 
485  template <class Evaluation>
486  static Evaluation gasThermalConductivity(const Evaluation& temperature, const Evaluation& pressure)
487  {
488  const Evaluation& result = interpolateGasTP_(gasThermalConductivity_,
489  temperature,
490  pressure);
491  if (std::isnan(Opm::scalarValue(result)))
492  return RawComponent::gasThermalConductivity(temperature, pressure);
493  return result;
494  }
495 
502  template <class Evaluation>
503  static Evaluation liquidThermalConductivity(const Evaluation& temperature, const Evaluation& pressure)
504  {
505  const Evaluation& result = interpolateLiquidTP_(liquidThermalConductivity_,
506  temperature,
507  pressure);
508  if (std::isnan(Opm::scalarValue(result)))
509  return RawComponent::liquidThermalConductivity(temperature, pressure);
510  return result;
511  }
512 
513 private:
514  // returns an interpolated value depending on temperature
515  template <class Evaluation>
516  static Evaluation interpolateT_(const Scalar* values, const Evaluation& T)
517  {
518  Evaluation alphaT = tempIdx_(T);
519  if (alphaT < 0 || alphaT >= nTemp_ - 1)
520  return std::numeric_limits<Scalar>::quiet_NaN();
521 
522  size_t iT = static_cast<size_t>(Opm::scalarValue(alphaT));
523  alphaT -= iT;
524 
525  return
526  values[iT ]*(1 - alphaT) +
527  values[iT + 1]*( alphaT);
528  }
529 
530  // returns an interpolated value for liquid depending on
531  // temperature and pressure
532  template <class Evaluation>
533  static Evaluation interpolateLiquidTP_(const Scalar* values, const Evaluation& T, const Evaluation& p)
534  {
535  Evaluation alphaT = tempIdx_(T);
536  if (alphaT < 0 || alphaT >= nTemp_ - 1)
537  return std::numeric_limits<Scalar>::quiet_NaN();
538 
539  size_t iT = static_cast<size_t>(Opm::scalarValue(alphaT));
540  alphaT -= iT;
541 
542  Evaluation alphaP1 = pressLiquidIdx_(p, iT);
543  Evaluation alphaP2 = pressLiquidIdx_(p, iT + 1);
544 
545  size_t iP1 =
546  static_cast<size_t>(
547  std::max<int>(0, std::min(static_cast<int>(nPress_) - 2,
548  static_cast<int>(Opm::scalarValue(alphaP1)))));
549  size_t iP2 =
550  static_cast<size_t>(
551  std::max(0, std::min(static_cast<int>(nPress_) - 2,
552  static_cast<int>(Opm::scalarValue(alphaP2)))));
553  alphaP1 -= iP1;
554  alphaP2 -= iP2;
555 
556  return
557  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
558  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
559  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
560  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
561  }
562 
563  // returns an interpolated value for gas depending on
564  // temperature and pressure
565  template <class Evaluation>
566  static Evaluation interpolateGasTP_(const Scalar* values, const Evaluation& T, const Evaluation& p)
567  {
568  Evaluation alphaT = tempIdx_(T);
569  if (alphaT < 0 || alphaT >= nTemp_ - 1)
570  return std::numeric_limits<Scalar>::quiet_NaN();
571 
572  size_t iT =
573  static_cast<size_t>(
574  std::max(0, std::min(static_cast<int>(nTemp_) - 2,
575  static_cast<int>(Opm::scalarValue(alphaT)))));
576  alphaT -= iT;
577 
578  Evaluation alphaP1 = pressGasIdx_(p, iT);
579  Evaluation alphaP2 = pressGasIdx_(p, iT + 1);
580  size_t iP1 =
581  static_cast<size_t>(
582  std::max(0, std::min(static_cast<int>(nPress_) - 2,
583  static_cast<int>(Opm::scalarValue(alphaP1)))));
584  size_t iP2 =
585  static_cast<size_t>(
586  std::max(0, std::min(static_cast<int>(nPress_) - 2,
587  static_cast<int>(Opm::scalarValue(alphaP2)))));
588  alphaP1 -= iP1;
589  alphaP2 -= iP2;
590 
591  return
592  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
593  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
594  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
595  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
596  }
597 
598  // returns an interpolated value for gas depending on
599  // temperature and density
600  template <class Evaluation>
601  static Evaluation interpolateGasTRho_(const Scalar* values, const Evaluation& T, const Evaluation& rho)
602  {
603  Evaluation alphaT = tempIdx_(T);
604  unsigned iT = std::max(0,
605  std::min(static_cast<int>(nTemp_ - 2),
606  static_cast<int>(alphaT)));
607  alphaT -= iT;
608 
609  Evaluation alphaP1 = densityGasIdx_(rho, iT);
610  Evaluation alphaP2 = densityGasIdx_(rho, iT + 1);
611  unsigned iP1 =
612  std::max(0,
613  std::min(static_cast<int>(nDensity_ - 2),
614  static_cast<int>(alphaP1)));
615  unsigned iP2 =
616  std::max(0,
617  std::min(static_cast<int>(nDensity_ - 2),
618  static_cast<int>(alphaP2)));
619  alphaP1 -= iP1;
620  alphaP2 -= iP2;
621 
622  return
623  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
624  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
625  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
626  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
627  }
628 
629  // returns an interpolated value for liquid depending on
630  // temperature and density
631  template <class Evaluation>
632  static Evaluation interpolateLiquidTRho_(const Scalar* values, const Evaluation& T, const Evaluation& rho)
633  {
634  Evaluation alphaT = tempIdx_(T);
635  unsigned iT = std::max<int>(0, std::min<int>(nTemp_ - 2, static_cast<int>(alphaT)));
636  alphaT -= iT;
637 
638  Evaluation alphaP1 = densityLiquidIdx_(rho, iT);
639  Evaluation alphaP2 = densityLiquidIdx_(rho, iT + 1);
640  unsigned iP1 = std::max<int>(0, std::min<int>(nDensity_ - 2, static_cast<int>(alphaP1)));
641  unsigned iP2 = std::max<int>(0, std::min<int>(nDensity_ - 2, static_cast<int>(alphaP2)));
642  alphaP1 -= iP1;
643  alphaP2 -= iP2;
644 
645  return
646  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
647  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
648  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
649  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
650  }
651 
652 
653  // returns the index of an entry in a temperature field
654  template <class Evaluation>
655  static Evaluation tempIdx_(const Evaluation& temperature)
656  {
657  return (nTemp_ - 1)*(temperature - tempMin_)/(tempMax_ - tempMin_);
658  }
659 
660  // returns the index of an entry in a pressure field
661  template <class Evaluation>
662  static Evaluation pressLiquidIdx_(const Evaluation& pressure, size_t tempIdx)
663  {
664  Scalar plMin = minLiquidPressure_(tempIdx);
665  Scalar plMax = maxLiquidPressure_(tempIdx);
666 
667  return (nPress_ - 1)*(pressure - plMin)/(plMax - plMin);
668  }
669 
670  // returns the index of an entry in a temperature field
671  template <class Evaluation>
672  static Evaluation pressGasIdx_(const Evaluation& pressure, size_t tempIdx)
673  {
674  Scalar pgMin = minGasPressure_(tempIdx);
675  Scalar pgMax = maxGasPressure_(tempIdx);
676 
677  return (nPress_ - 1)*(pressure - pgMin)/(pgMax - pgMin);
678  }
679 
680  // returns the index of an entry in a density field
681  template <class Evaluation>
682  static Evaluation densityLiquidIdx_(const Evaluation& density, size_t tempIdx)
683  {
684  Scalar densityMin = minLiquidDensity_(tempIdx);
685  Scalar densityMax = maxLiquidDensity_(tempIdx);
686  return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
687  }
688 
689  // returns the index of an entry in a density field
690  template <class Evaluation>
691  static Evaluation densityGasIdx_(const Evaluation& density, size_t tempIdx)
692  {
693  Scalar densityMin = minGasDensity_(tempIdx);
694  Scalar densityMax = maxGasDensity_(tempIdx);
695  return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
696  }
697 
698  // returns the minimum tabulized liquid pressure at a given
699  // temperature index
700  static Scalar minLiquidPressure_(size_t tempIdx)
701  {
702  if (!useVaporPressure)
703  return pressMin_;
704  else
705  return std::max<Scalar>(pressMin_, vaporPressure_[tempIdx] / 1.1);
706  }
707 
708  // returns the maximum tabulized liquid pressure at a given
709  // temperature index
710  static Scalar maxLiquidPressure_(size_t tempIdx)
711  {
712  if (!useVaporPressure)
713  return pressMax_;
714  else
715  return std::max<Scalar>(pressMax_, vaporPressure_[tempIdx] * 1.1);
716  }
717 
718  // returns the minumum tabulized gas pressure at a given
719  // temperature index
720  static Scalar minGasPressure_(size_t tempIdx)
721  {
722  if (!useVaporPressure)
723  return pressMin_;
724  else
725  return std::min<Scalar>(pressMin_, vaporPressure_[tempIdx] / 1.1 );
726  }
727 
728  // returns the maximum tabulized gas pressure at a given
729  // temperature index
730  static Scalar maxGasPressure_(size_t tempIdx)
731  {
732  if (!useVaporPressure)
733  return pressMax_;
734  else
735  return std::min<Scalar>(pressMax_, vaporPressure_[tempIdx] * 1.1);
736  }
737 
738 
739  // returns the minimum tabulized liquid density at a given
740  // temperature index
741  static Scalar minLiquidDensity_(size_t tempIdx)
742  { return minLiquidDensity__[tempIdx]; }
743 
744  // returns the maximum tabulized liquid density at a given
745  // temperature index
746  static Scalar maxLiquidDensity_(size_t tempIdx)
747  { return maxLiquidDensity__[tempIdx]; }
748 
749  // returns the minumum tabulized gas density at a given
750  // temperature index
751  static Scalar minGasDensity_(size_t tempIdx)
752  { return minGasDensity__[tempIdx]; }
753 
754  // returns the maximum tabulized gas density at a given
755  // temperature index
756  static Scalar maxGasDensity_(size_t tempIdx)
757  { return maxGasDensity__[tempIdx]; }
758 
759  // 1D fields with the temperature as degree of freedom
760  static Scalar* vaporPressure_;
761 
762  static Scalar* minLiquidDensity__;
763  static Scalar* maxLiquidDensity__;
764 
765  static Scalar* minGasDensity__;
766  static Scalar* maxGasDensity__;
767 
768  // 2D fields with the temperature and pressure as degrees of
769  // freedom
770  static Scalar* gasEnthalpy_;
771  static Scalar* liquidEnthalpy_;
772 
773  static Scalar* gasHeatCapacity_;
774  static Scalar* liquidHeatCapacity_;
775 
776  static Scalar* gasDensity_;
777  static Scalar* liquidDensity_;
778 
779  static Scalar* gasViscosity_;
780  static Scalar* liquidViscosity_;
781 
782  static Scalar* gasThermalConductivity_;
783  static Scalar* liquidThermalConductivity_;
784 
785  // 2D fields with the temperature and density as degrees of
786  // freedom
787  static Scalar* gasPressure_;
788  static Scalar* liquidPressure_;
789 
790  // temperature, pressure and density ranges
791  static Scalar tempMin_;
792  static Scalar tempMax_;
793  static unsigned nTemp_;
794 
795  static Scalar pressMin_;
796  static Scalar pressMax_;
797  static unsigned nPress_;
798 
799  static Scalar densityMin_;
800  static Scalar densityMax_;
801  static unsigned nDensity_;
802 };
803 
804 template <class Scalar, class RawComponent, bool useVaporPressure>
805 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::vaporPressure_;
806 template <class Scalar, class RawComponent, bool useVaporPressure>
807 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::minLiquidDensity__;
808 template <class Scalar, class RawComponent, bool useVaporPressure>
809 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::maxLiquidDensity__;
810 template <class Scalar, class RawComponent, bool useVaporPressure>
811 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::minGasDensity__;
812 template <class Scalar, class RawComponent, bool useVaporPressure>
813 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::maxGasDensity__;
814 template <class Scalar, class RawComponent, bool useVaporPressure>
815 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasEnthalpy_;
816 template <class Scalar, class RawComponent, bool useVaporPressure>
817 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidEnthalpy_;
818 template <class Scalar, class RawComponent, bool useVaporPressure>
819 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasHeatCapacity_;
820 template <class Scalar, class RawComponent, bool useVaporPressure>
821 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidHeatCapacity_;
822 template <class Scalar, class RawComponent, bool useVaporPressure>
823 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasDensity_;
824 template <class Scalar, class RawComponent, bool useVaporPressure>
825 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidDensity_;
826 template <class Scalar, class RawComponent, bool useVaporPressure>
827 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasViscosity_;
828 template <class Scalar, class RawComponent, bool useVaporPressure>
829 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidViscosity_;
830 template <class Scalar, class RawComponent, bool useVaporPressure>
831 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasThermalConductivity_;
832 template <class Scalar, class RawComponent, bool useVaporPressure>
833 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidThermalConductivity_;
834 template <class Scalar, class RawComponent, bool useVaporPressure>
835 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasPressure_;
836 template <class Scalar, class RawComponent, bool useVaporPressure>
837 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidPressure_;
838 template <class Scalar, class RawComponent, bool useVaporPressure>
839 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::tempMin_;
840 template <class Scalar, class RawComponent, bool useVaporPressure>
841 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::tempMax_;
842 template <class Scalar, class RawComponent, bool useVaporPressure>
843 unsigned TabulatedComponent<Scalar, RawComponent, useVaporPressure>::nTemp_;
844 template <class Scalar, class RawComponent, bool useVaporPressure>
845 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::pressMin_;
846 template <class Scalar, class RawComponent, bool useVaporPressure>
847 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::pressMax_;
848 template <class Scalar, class RawComponent, bool useVaporPressure>
849 unsigned TabulatedComponent<Scalar, RawComponent, useVaporPressure>::nPress_;
850 template <class Scalar, class RawComponent, bool useVaporPressure>
851 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::densityMin_;
852 template <class Scalar, class RawComponent, bool useVaporPressure>
853 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::densityMax_;
854 template <class Scalar, class RawComponent, bool useVaporPressure>
855 unsigned TabulatedComponent<Scalar, RawComponent, useVaporPressure>::nDensity_;
856 
857 
858 } // namespace Opm
859 
860 #endif
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of liquid at a given pressure and temperature .
Definition: TabulatedComponent.hpp:435
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the liquid .
Definition: TabulatedComponent.hpp:290
static Evaluation liquidHeatCapacity(const Evaluation &temperature, const Evaluation &pressure)
Specific isobaric heat capacity of the liquid .
Definition: TabulatedComponent.hpp:324
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
static Scalar molarMass()
The molar mass in of the component.
Definition: TabulatedComponent.hpp:224
static Scalar triplePressure()
Returns the pressure in at the component&#39;s triple point.
Definition: TabulatedComponent.hpp:248
static Evaluation liquidInternalEnergy(const Evaluation &temperature, const Evaluation &pressure)
Specific internal energy of the liquid .
Definition: TabulatedComponent.hpp:351
static Evaluation gasViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of gas.
Definition: TabulatedComponent.hpp:452
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of gas at a given pressure and temperature .
Definition: TabulatedComponent.hpp:417
static Evaluation liquidThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
The thermal conductivity of liquid water .
Definition: TabulatedComponent.hpp:503
static const char * name()
A human readable name for the component.
Definition: TabulatedComponent.hpp:218
static Evaluation vaporPressure(const Evaluation &temperature)
The vapor pressure in of the component at a given temperature.
Definition: TabulatedComponent.hpp:258
static Evaluation gasInternalEnergy(const Evaluation &temperature, const Evaluation &pressure)
Specific internal energy of the gas .
Definition: TabulatedComponent.hpp:341
static Scalar criticalPressure()
Returns the critical pressure in of the component.
Definition: TabulatedComponent.hpp:236
static Evaluation gasHeatCapacity(const Evaluation &temperature, const Evaluation &pressure)
Specific isobaric heat capacity of the gas .
Definition: TabulatedComponent.hpp:307
static Scalar criticalTemperature()
Returns the critical temperature in of the component.
Definition: TabulatedComponent.hpp:230
static Evaluation liquidPressure(const Evaluation &temperature, Scalar density)
The pressure of liquid in at a given density and temperature.
Definition: TabulatedComponent.hpp:379
A generic class which tabulates all thermodynamic properties of a given component.
Definition: TabulatedComponent.hpp:58
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the tables.
Definition: TabulatedComponent.hpp:75
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the gas .
Definition: TabulatedComponent.hpp:273
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: TabulatedComponent.hpp:405
static Evaluation liquidViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of liquid.
Definition: TabulatedComponent.hpp:469
static Evaluation gasThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
The thermal conductivity of gaseous water .
Definition: TabulatedComponent.hpp:486
static bool gasIsCompressible()
Returns true iff the gas phase is assumed to be compressible.
Definition: TabulatedComponent.hpp:393
static Evaluation gasPressure(const Evaluation &temperature, Scalar density)
The pressure of gas in at a given density and temperature.
Definition: TabulatedComponent.hpp:361
static Scalar tripleTemperature()
Returns the temperature in at the component&#39;s triple point.
Definition: TabulatedComponent.hpp:242
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: TabulatedComponent.hpp:399