My Project
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 
37 
38 namespace Opm {
54 template <class ScalarT, class RawComponent, bool useVaporPressure=true>
56 {
57 public:
58  typedef ScalarT Scalar;
59 
60  static const bool isTabulated = true;
61 
72  static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
73  Scalar pressMin, Scalar pressMax, unsigned nPress)
74  {
75  tempMin_ = tempMin;
76  tempMax_ = tempMax;
77  nTemp_ = nTemp;
78  pressMin_ = pressMin;
79  pressMax_ = pressMax;
80  nPress_ = nPress;
81  nDensity_ = nPress_;
82 
83  // allocate the arrays
84  vaporPressure_ = new Scalar[nTemp_];
85  minGasDensity__ = new Scalar[nTemp_];
86  maxGasDensity__ = new Scalar[nTemp_];
87  minLiquidDensity__ = new Scalar[nTemp_];
88  maxLiquidDensity__ = new Scalar[nTemp_];
89 
90  gasEnthalpy_ = new Scalar[nTemp_*nPress_];
91  liquidEnthalpy_ = new Scalar[nTemp_*nPress_];
92  gasHeatCapacity_ = new Scalar[nTemp_*nPress_];
93  liquidHeatCapacity_ = new Scalar[nTemp_*nPress_];
94  gasDensity_ = new Scalar[nTemp_*nPress_];
95  liquidDensity_ = new Scalar[nTemp_*nPress_];
96  gasViscosity_ = new Scalar[nTemp_*nPress_];
97  liquidViscosity_ = new Scalar[nTemp_*nPress_];
98  gasThermalConductivity_ = new Scalar[nTemp_*nPress_];
99  liquidThermalConductivity_ = new Scalar[nTemp_*nPress_];
100  gasPressure_ = new Scalar[nTemp_*nDensity_];
101  liquidPressure_ = new Scalar[nTemp_*nDensity_];
102 
103  assert(std::numeric_limits<Scalar>::has_quiet_NaN);
104  Scalar NaN = std::numeric_limits<Scalar>::quiet_NaN();
105 
106  // fill the temperature-pressure arrays
107  for (unsigned iT = 0; iT < nTemp_; ++ iT) {
108  Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
109 
110  try { vaporPressure_[iT] = RawComponent::vaporPressure(temperature); }
111  catch (const std::exception&) { vaporPressure_[iT] = NaN; }
112 
113  Scalar pgMax = maxGasPressure_(iT);
114  Scalar pgMin = minGasPressure_(iT);
115 
116  // fill the temperature, pressure gas arrays
117  for (unsigned iP = 0; iP < nPress_; ++ iP) {
118  Scalar pressure = iP * (pgMax - pgMin)/(nPress_ - 1) + pgMin;
119 
120  unsigned i = iT + iP*nTemp_;
121 
122  try { gasEnthalpy_[i] = RawComponent::gasEnthalpy(temperature, pressure); }
123  catch (const std::exception&) { gasEnthalpy_[i] = NaN; }
124 
125  try { gasHeatCapacity_[i] = RawComponent::gasHeatCapacity(temperature, pressure); }
126  catch (const std::exception&) { gasHeatCapacity_[i] = NaN; }
127 
128  try { gasDensity_[i] = RawComponent::gasDensity(temperature, pressure); }
129  catch (const std::exception&) { gasDensity_[i] = NaN; }
130 
131  try { gasViscosity_[i] = RawComponent::gasViscosity(temperature, pressure); }
132  catch (const std::exception&) { gasViscosity_[i] = NaN; }
133 
134  try { gasThermalConductivity_[i] = RawComponent::gasThermalConductivity(temperature, pressure); }
135  catch (const std::exception&) { gasThermalConductivity_[i] = NaN; }
136  };
137 
138  Scalar plMin = minLiquidPressure_(iT);
139  Scalar plMax = maxLiquidPressure_(iT);
140  for (unsigned iP = 0; iP < nPress_; ++ iP) {
141  Scalar pressure = iP * (plMax - plMin)/(nPress_ - 1) + plMin;
142 
143  unsigned i = iT + iP*nTemp_;
144 
145  try { liquidEnthalpy_[i] = RawComponent::liquidEnthalpy(temperature, pressure); }
146  catch (const std::exception&) { liquidEnthalpy_[i] = NaN; }
147 
148  try { liquidHeatCapacity_[i] = RawComponent::liquidHeatCapacity(temperature, pressure); }
149  catch (const std::exception&) { liquidHeatCapacity_[i] = NaN; }
150 
151  try { liquidDensity_[i] = RawComponent::liquidDensity(temperature, pressure); }
152  catch (const std::exception&) { liquidDensity_[i] = NaN; }
153 
154  try { liquidViscosity_[i] = RawComponent::liquidViscosity(temperature, pressure); }
155  catch (const std::exception&) { liquidViscosity_[i] = NaN; }
156 
157  try { liquidThermalConductivity_[i] = RawComponent::liquidThermalConductivity(temperature, pressure); }
158  catch (const std::exception&) { liquidThermalConductivity_[i] = NaN; }
159  }
160  }
161 
162  // fill the temperature-density arrays
163  for (unsigned iT = 0; iT < nTemp_; ++ iT) {
164  Scalar temperature = iT * (tempMax_ - tempMin_)/(nTemp_ - 1) + tempMin_;
165 
166  // calculate the minimum and maximum values for the gas
167  // densities
168  minGasDensity__[iT] = RawComponent::gasDensity(temperature, minGasPressure_(iT));
169  if (iT < nTemp_ - 1)
170  maxGasDensity__[iT] = RawComponent::gasDensity(temperature, maxGasPressure_(iT + 1));
171  else
172  maxGasDensity__[iT] = RawComponent::gasDensity(temperature, maxGasPressure_(iT));
173 
174  // fill the temperature, density gas arrays
175  for (unsigned iRho = 0; iRho < nDensity_; ++ iRho) {
176  Scalar density =
177  Scalar(iRho)/(nDensity_ - 1) *
178  (maxGasDensity__[iT] - minGasDensity__[iT])
179  +
180  minGasDensity__[iT];
181 
182  unsigned i = iT + iRho*nTemp_;
183 
184  try { gasPressure_[i] = RawComponent::gasPressure(temperature, density); }
185  catch (const std::exception&) { gasPressure_[i] = NaN; };
186  };
187 
188  // calculate the minimum and maximum values for the liquid
189  // densities
190  minLiquidDensity__[iT] = RawComponent::liquidDensity(temperature, minLiquidPressure_(iT));
191  if (iT < nTemp_ - 1)
192  maxLiquidDensity__[iT] = RawComponent::liquidDensity(temperature, maxLiquidPressure_(iT + 1));
193  else
194  maxLiquidDensity__[iT] = RawComponent::liquidDensity(temperature, maxLiquidPressure_(iT));
195 
196  // fill the temperature, density liquid arrays
197  for (unsigned iRho = 0; iRho < nDensity_; ++ iRho) {
198  Scalar density =
199  Scalar(iRho)/(nDensity_ - 1) *
200  (maxLiquidDensity__[iT] - minLiquidDensity__[iT])
201  +
202  minLiquidDensity__[iT];
203 
204  unsigned i = iT + iRho*nTemp_;
205 
206  try { liquidPressure_[i] = RawComponent::liquidPressure(temperature, density); }
207  catch (const std::exception&) { liquidPressure_[i] = NaN; };
208  };
209  }
210  }
211 
215  static const char* name()
216  { return RawComponent::name(); }
217 
221  static Scalar molarMass()
222  { return RawComponent::molarMass(); }
223 
227  static Scalar criticalTemperature()
228  { return RawComponent::criticalTemperature(); }
229 
233  static Scalar criticalPressure()
234  { return RawComponent::criticalPressure(); }
235 
239  static Scalar acentricFactor()
240  { throw std::runtime_error("Not implemented: acentricFactor of the component"); }
241 
245  static Scalar criticalVolume()
246  { throw std::runtime_error("Not implemented: criticalVolume of the compoenent"); }
247 
251  static Scalar tripleTemperature()
252  { return RawComponent::tripleTemperature(); }
253 
257  static Scalar triplePressure()
258  { return RawComponent::triplePressure(); }
259 
266  template <class Evaluation>
267  static Evaluation vaporPressure(const Evaluation& temperature)
268  {
269  const Evaluation& result = interpolateT_(vaporPressure_, temperature);
270  if (std::isnan(scalarValue(result)))
271  return RawComponent::vaporPressure(temperature);
272  return result;
273  }
274 
281  template <class Evaluation>
282  static Evaluation gasEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
283  {
284  const Evaluation& result = interpolateGasTP_(gasEnthalpy_,
285  temperature,
286  pressure);
287  if (std::isnan(scalarValue(result)))
288  return RawComponent::gasEnthalpy(temperature, pressure);
289  return result;
290  }
291 
298  template <class Evaluation>
299  static Evaluation liquidEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
300  {
301  const Evaluation& result = interpolateLiquidTP_(liquidEnthalpy_,
302  temperature,
303  pressure);
304  if (std::isnan(scalarValue(result)))
305  return RawComponent::liquidEnthalpy(temperature, pressure);
306  return result;
307  }
308 
315  template <class Evaluation>
316  static Evaluation gasHeatCapacity(const Evaluation& temperature, const Evaluation& pressure)
317  {
318  const Evaluation& result = interpolateGasTP_(gasHeatCapacity_,
319  temperature,
320  pressure);
321  if (std::isnan(scalarValue(result)))
322  return RawComponent::gasHeatCapacity(temperature, pressure);
323  return result;
324  }
325 
332  template <class Evaluation>
333  static Evaluation liquidHeatCapacity(const Evaluation& temperature, const Evaluation& pressure)
334  {
335  const Evaluation& result = interpolateLiquidTP_(liquidHeatCapacity_,
336  temperature,
337  pressure);
338  if (std::isnan(scalarValue(result)))
339  return RawComponent::liquidHeatCapacity(temperature, pressure);
340  return result;
341  }
342 
349  template <class Evaluation>
350  static Evaluation gasInternalEnergy(const Evaluation& temperature, const Evaluation& pressure)
351  { return gasEnthalpy(temperature, pressure) - pressure/gasDensity(temperature, pressure); }
352 
359  template <class Evaluation>
360  static Evaluation liquidInternalEnergy(const Evaluation& temperature, const Evaluation& pressure)
361  { return liquidEnthalpy(temperature, pressure) - pressure/liquidDensity(temperature, pressure); }
362 
369  template <class Evaluation>
370  static Evaluation gasPressure(const Evaluation& temperature, Scalar density)
371  {
372  const Evaluation& result = interpolateGasTRho_(gasPressure_,
373  temperature,
374  density);
375  if (std::isnan(scalarValue(result)))
376  return RawComponent::gasPressure(temperature,
377  density);
378  return result;
379  }
380 
387  template <class Evaluation>
388  static Evaluation liquidPressure(const Evaluation& temperature, Scalar density)
389  {
390  const Evaluation& result = interpolateLiquidTRho_(liquidPressure_,
391  temperature,
392  density);
393  if (std::isnan(scalarValue(result)))
394  return RawComponent::liquidPressure(temperature,
395  density);
396  return result;
397  }
398 
402  static bool gasIsCompressible()
403  { return RawComponent::gasIsCompressible(); }
404 
408  static bool liquidIsCompressible()
409  { return RawComponent::liquidIsCompressible(); }
410 
414  static bool gasIsIdeal()
415  { return RawComponent::gasIsIdeal(); }
416 
417 
425  template <class Evaluation>
426  static Evaluation gasDensity(const Evaluation& temperature, const Evaluation& pressure)
427  {
428  const Evaluation& result = interpolateGasTP_(gasDensity_,
429  temperature,
430  pressure);
431  if (std::isnan(scalarValue(result)))
432  return RawComponent::gasDensity(temperature, pressure);
433  return result;
434  }
435 
443  template <class Evaluation>
444  static Evaluation liquidDensity(const Evaluation& temperature, const Evaluation& pressure)
445  {
446  const Evaluation& result = interpolateLiquidTP_(liquidDensity_,
447  temperature,
448  pressure);
449  if (std::isnan(scalarValue(result)))
450  return RawComponent::liquidDensity(temperature, pressure);
451  return result;
452  }
453 
460  template <class Evaluation>
461  static Evaluation gasViscosity(const Evaluation& temperature, const Evaluation& pressure)
462  {
463  const Evaluation& result = interpolateGasTP_(gasViscosity_,
464  temperature,
465  pressure);
466  if (std::isnan(scalarValue(result)))
467  return RawComponent::gasViscosity(temperature, pressure);
468  return result;
469  }
470 
477  template <class Evaluation>
478  static Evaluation liquidViscosity(const Evaluation& temperature, const Evaluation& pressure)
479  {
480  const Evaluation& result = interpolateLiquidTP_(liquidViscosity_,
481  temperature,
482  pressure);
483  if (std::isnan(scalarValue(result)))
484  return RawComponent::liquidViscosity(temperature, pressure);
485  return result;
486  }
487 
494  template <class Evaluation>
495  static Evaluation gasThermalConductivity(const Evaluation& temperature, const Evaluation& pressure)
496  {
497  const Evaluation& result = interpolateGasTP_(gasThermalConductivity_,
498  temperature,
499  pressure);
500  if (std::isnan(scalarValue(result)))
501  return RawComponent::gasThermalConductivity(temperature, pressure);
502  return result;
503  }
504 
511  template <class Evaluation>
512  static Evaluation liquidThermalConductivity(const Evaluation& temperature, const Evaluation& pressure)
513  {
514  const Evaluation& result = interpolateLiquidTP_(liquidThermalConductivity_,
515  temperature,
516  pressure);
517  if (std::isnan(scalarValue(result)))
518  return RawComponent::liquidThermalConductivity(temperature, pressure);
519  return result;
520  }
521 
522 private:
523  // returns an interpolated value depending on temperature
524  template <class Evaluation>
525  static Evaluation interpolateT_(const Scalar* values, const Evaluation& T)
526  {
527  Evaluation alphaT = tempIdx_(T);
528  if (alphaT < 0 || alphaT >= nTemp_ - 1)
529  return std::numeric_limits<Scalar>::quiet_NaN();
530 
531  size_t iT = static_cast<size_t>(scalarValue(alphaT));
532  alphaT -= iT;
533 
534  return
535  values[iT ]*(1 - alphaT) +
536  values[iT + 1]*( alphaT);
537  }
538 
539  // returns an interpolated value for liquid depending on
540  // temperature and pressure
541  template <class Evaluation>
542  static Evaluation interpolateLiquidTP_(const Scalar* values, const Evaluation& T, const Evaluation& p)
543  {
544  Evaluation alphaT = tempIdx_(T);
545  if (alphaT < 0 || alphaT >= nTemp_ - 1)
546  return std::numeric_limits<Scalar>::quiet_NaN();
547 
548  size_t iT = static_cast<size_t>(scalarValue(alphaT));
549  alphaT -= iT;
550 
551  Evaluation alphaP1 = pressLiquidIdx_(p, iT);
552  Evaluation alphaP2 = pressLiquidIdx_(p, iT + 1);
553 
554  size_t iP1 =
555  static_cast<size_t>(
556  std::max<int>(0, std::min(static_cast<int>(nPress_) - 2,
557  static_cast<int>(scalarValue(alphaP1)))));
558  size_t iP2 =
559  static_cast<size_t>(
560  std::max(0, std::min(static_cast<int>(nPress_) - 2,
561  static_cast<int>(scalarValue(alphaP2)))));
562  alphaP1 -= iP1;
563  alphaP2 -= iP2;
564 
565  return
566  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
567  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
568  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
569  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
570  }
571 
572  // returns an interpolated value for gas depending on
573  // temperature and pressure
574  template <class Evaluation>
575  static Evaluation interpolateGasTP_(const Scalar* values, const Evaluation& T, const Evaluation& p)
576  {
577  Evaluation alphaT = tempIdx_(T);
578  if (alphaT < 0 || alphaT >= nTemp_ - 1)
579  return std::numeric_limits<Scalar>::quiet_NaN();
580 
581  size_t iT =
582  static_cast<size_t>(
583  std::max(0, std::min(static_cast<int>(nTemp_) - 2,
584  static_cast<int>(scalarValue(alphaT)))));
585  alphaT -= iT;
586 
587  Evaluation alphaP1 = pressGasIdx_(p, iT);
588  Evaluation alphaP2 = pressGasIdx_(p, iT + 1);
589  size_t iP1 =
590  static_cast<size_t>(
591  std::max(0, std::min(static_cast<int>(nPress_) - 2,
592  static_cast<int>(scalarValue(alphaP1)))));
593  size_t iP2 =
594  static_cast<size_t>(
595  std::max(0, std::min(static_cast<int>(nPress_) - 2,
596  static_cast<int>(scalarValue(alphaP2)))));
597  alphaP1 -= iP1;
598  alphaP2 -= iP2;
599 
600  return
601  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
602  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
603  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
604  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
605  }
606 
607  // returns an interpolated value for gas depending on
608  // temperature and density
609  template <class Evaluation>
610  static Evaluation interpolateGasTRho_(const Scalar* values, const Evaluation& T, const Evaluation& rho)
611  {
612  Evaluation alphaT = tempIdx_(T);
613  unsigned iT = std::max(0,
614  std::min(static_cast<int>(nTemp_ - 2),
615  static_cast<int>(alphaT)));
616  alphaT -= iT;
617 
618  Evaluation alphaP1 = densityGasIdx_(rho, iT);
619  Evaluation alphaP2 = densityGasIdx_(rho, iT + 1);
620  unsigned iP1 =
621  std::max(0,
622  std::min(static_cast<int>(nDensity_ - 2),
623  static_cast<int>(alphaP1)));
624  unsigned iP2 =
625  std::max(0,
626  std::min(static_cast<int>(nDensity_ - 2),
627  static_cast<int>(alphaP2)));
628  alphaP1 -= iP1;
629  alphaP2 -= iP2;
630 
631  return
632  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
633  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
634  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
635  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
636  }
637 
638  // returns an interpolated value for liquid depending on
639  // temperature and density
640  template <class Evaluation>
641  static Evaluation interpolateLiquidTRho_(const Scalar* values, const Evaluation& T, const Evaluation& rho)
642  {
643  Evaluation alphaT = tempIdx_(T);
644  unsigned iT = std::max<int>(0, std::min<int>(nTemp_ - 2, static_cast<int>(alphaT)));
645  alphaT -= iT;
646 
647  Evaluation alphaP1 = densityLiquidIdx_(rho, iT);
648  Evaluation alphaP2 = densityLiquidIdx_(rho, iT + 1);
649  unsigned iP1 = std::max<int>(0, std::min<int>(nDensity_ - 2, static_cast<int>(alphaP1)));
650  unsigned iP2 = std::max<int>(0, std::min<int>(nDensity_ - 2, static_cast<int>(alphaP2)));
651  alphaP1 -= iP1;
652  alphaP2 -= iP2;
653 
654  return
655  values[(iT ) + (iP1 )*nTemp_]*(1 - alphaT)*(1 - alphaP1) +
656  values[(iT ) + (iP1 + 1)*nTemp_]*(1 - alphaT)*( alphaP1) +
657  values[(iT + 1) + (iP2 )*nTemp_]*( alphaT)*(1 - alphaP2) +
658  values[(iT + 1) + (iP2 + 1)*nTemp_]*( alphaT)*( alphaP2);
659  }
660 
661 
662  // returns the index of an entry in a temperature field
663  template <class Evaluation>
664  static Evaluation tempIdx_(const Evaluation& temperature)
665  {
666  return (nTemp_ - 1)*(temperature - tempMin_)/(tempMax_ - tempMin_);
667  }
668 
669  // returns the index of an entry in a pressure field
670  template <class Evaluation>
671  static Evaluation pressLiquidIdx_(const Evaluation& pressure, size_t tempIdx)
672  {
673  Scalar plMin = minLiquidPressure_(tempIdx);
674  Scalar plMax = maxLiquidPressure_(tempIdx);
675 
676  return (nPress_ - 1)*(pressure - plMin)/(plMax - plMin);
677  }
678 
679  // returns the index of an entry in a temperature field
680  template <class Evaluation>
681  static Evaluation pressGasIdx_(const Evaluation& pressure, size_t tempIdx)
682  {
683  Scalar pgMin = minGasPressure_(tempIdx);
684  Scalar pgMax = maxGasPressure_(tempIdx);
685 
686  return (nPress_ - 1)*(pressure - pgMin)/(pgMax - pgMin);
687  }
688 
689  // returns the index of an entry in a density field
690  template <class Evaluation>
691  static Evaluation densityLiquidIdx_(const Evaluation& density, size_t tempIdx)
692  {
693  Scalar densityMin = minLiquidDensity_(tempIdx);
694  Scalar densityMax = maxLiquidDensity_(tempIdx);
695  return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
696  }
697 
698  // returns the index of an entry in a density field
699  template <class Evaluation>
700  static Evaluation densityGasIdx_(const Evaluation& density, size_t tempIdx)
701  {
702  Scalar densityMin = minGasDensity_(tempIdx);
703  Scalar densityMax = maxGasDensity_(tempIdx);
704  return (nDensity_ - 1) * (density - densityMin)/(densityMax - densityMin);
705  }
706 
707  // returns the minimum tabulized liquid pressure at a given
708  // temperature index
709  static Scalar minLiquidPressure_(size_t tempIdx)
710  {
711  if (!useVaporPressure)
712  return pressMin_;
713  else
714  return std::max<Scalar>(pressMin_, vaporPressure_[tempIdx] / 1.1);
715  }
716 
717  // returns the maximum tabulized liquid pressure at a given
718  // temperature index
719  static Scalar maxLiquidPressure_(size_t tempIdx)
720  {
721  if (!useVaporPressure)
722  return pressMax_;
723  else
724  return std::max<Scalar>(pressMax_, vaporPressure_[tempIdx] * 1.1);
725  }
726 
727  // returns the minumum tabulized gas pressure at a given
728  // temperature index
729  static Scalar minGasPressure_(size_t tempIdx)
730  {
731  if (!useVaporPressure)
732  return pressMin_;
733  else
734  return std::min<Scalar>(pressMin_, vaporPressure_[tempIdx] / 1.1 );
735  }
736 
737  // returns the maximum tabulized gas pressure at a given
738  // temperature index
739  static Scalar maxGasPressure_(size_t tempIdx)
740  {
741  if (!useVaporPressure)
742  return pressMax_;
743  else
744  return std::min<Scalar>(pressMax_, vaporPressure_[tempIdx] * 1.1);
745  }
746 
747 
748  // returns the minimum tabulized liquid density at a given
749  // temperature index
750  static Scalar minLiquidDensity_(size_t tempIdx)
751  { return minLiquidDensity__[tempIdx]; }
752 
753  // returns the maximum tabulized liquid density at a given
754  // temperature index
755  static Scalar maxLiquidDensity_(size_t tempIdx)
756  { return maxLiquidDensity__[tempIdx]; }
757 
758  // returns the minumum tabulized gas density at a given
759  // temperature index
760  static Scalar minGasDensity_(size_t tempIdx)
761  { return minGasDensity__[tempIdx]; }
762 
763  // returns the maximum tabulized gas density at a given
764  // temperature index
765  static Scalar maxGasDensity_(size_t tempIdx)
766  { return maxGasDensity__[tempIdx]; }
767 
768  // 1D fields with the temperature as degree of freedom
769  static Scalar* vaporPressure_;
770 
771  static Scalar* minLiquidDensity__;
772  static Scalar* maxLiquidDensity__;
773 
774  static Scalar* minGasDensity__;
775  static Scalar* maxGasDensity__;
776 
777  // 2D fields with the temperature and pressure as degrees of
778  // freedom
779  static Scalar* gasEnthalpy_;
780  static Scalar* liquidEnthalpy_;
781 
782  static Scalar* gasHeatCapacity_;
783  static Scalar* liquidHeatCapacity_;
784 
785  static Scalar* gasDensity_;
786  static Scalar* liquidDensity_;
787 
788  static Scalar* gasViscosity_;
789  static Scalar* liquidViscosity_;
790 
791  static Scalar* gasThermalConductivity_;
792  static Scalar* liquidThermalConductivity_;
793 
794  // 2D fields with the temperature and density as degrees of
795  // freedom
796  static Scalar* gasPressure_;
797  static Scalar* liquidPressure_;
798 
799  // temperature, pressure and density ranges
800  static Scalar tempMin_;
801  static Scalar tempMax_;
802  static unsigned nTemp_;
803 
804  static Scalar pressMin_;
805  static Scalar pressMax_;
806  static unsigned nPress_;
807 
808  static Scalar densityMin_;
809  static Scalar densityMax_;
810  static unsigned nDensity_;
811 };
812 
813 template <class Scalar, class RawComponent, bool useVaporPressure>
814 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::vaporPressure_;
815 template <class Scalar, class RawComponent, bool useVaporPressure>
816 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::minLiquidDensity__;
817 template <class Scalar, class RawComponent, bool useVaporPressure>
818 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::maxLiquidDensity__;
819 template <class Scalar, class RawComponent, bool useVaporPressure>
820 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::minGasDensity__;
821 template <class Scalar, class RawComponent, bool useVaporPressure>
822 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::maxGasDensity__;
823 template <class Scalar, class RawComponent, bool useVaporPressure>
824 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasEnthalpy_;
825 template <class Scalar, class RawComponent, bool useVaporPressure>
826 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidEnthalpy_;
827 template <class Scalar, class RawComponent, bool useVaporPressure>
828 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasHeatCapacity_;
829 template <class Scalar, class RawComponent, bool useVaporPressure>
830 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidHeatCapacity_;
831 template <class Scalar, class RawComponent, bool useVaporPressure>
832 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasDensity_;
833 template <class Scalar, class RawComponent, bool useVaporPressure>
834 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidDensity_;
835 template <class Scalar, class RawComponent, bool useVaporPressure>
836 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasViscosity_;
837 template <class Scalar, class RawComponent, bool useVaporPressure>
838 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidViscosity_;
839 template <class Scalar, class RawComponent, bool useVaporPressure>
840 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasThermalConductivity_;
841 template <class Scalar, class RawComponent, bool useVaporPressure>
842 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidThermalConductivity_;
843 template <class Scalar, class RawComponent, bool useVaporPressure>
844 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::gasPressure_;
845 template <class Scalar, class RawComponent, bool useVaporPressure>
846 Scalar* TabulatedComponent<Scalar, RawComponent, useVaporPressure>::liquidPressure_;
847 template <class Scalar, class RawComponent, bool useVaporPressure>
848 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::tempMin_;
849 template <class Scalar, class RawComponent, bool useVaporPressure>
850 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::tempMax_;
851 template <class Scalar, class RawComponent, bool useVaporPressure>
852 unsigned TabulatedComponent<Scalar, RawComponent, useVaporPressure>::nTemp_;
853 template <class Scalar, class RawComponent, bool useVaporPressure>
854 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::pressMin_;
855 template <class Scalar, class RawComponent, bool useVaporPressure>
856 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::pressMax_;
857 template <class Scalar, class RawComponent, bool useVaporPressure>
858 unsigned TabulatedComponent<Scalar, RawComponent, useVaporPressure>::nPress_;
859 template <class Scalar, class RawComponent, bool useVaporPressure>
860 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::densityMin_;
861 template <class Scalar, class RawComponent, bool useVaporPressure>
862 Scalar TabulatedComponent<Scalar, RawComponent, useVaporPressure>::densityMax_;
863 template <class Scalar, class RawComponent, bool useVaporPressure>
864 unsigned TabulatedComponent<Scalar, RawComponent, useVaporPressure>::nDensity_;
865 
866 
867 } // namespace Opm
868 
869 #endif
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
A generic class which tabulates all thermodynamic properties of a given component.
Definition: TabulatedComponent.hpp:56
static Evaluation gasThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
The thermal conductivity of gaseous water .
Definition: TabulatedComponent.hpp:495
static bool gasIsCompressible()
Returns true iff the gas phase is assumed to be compressible.
Definition: TabulatedComponent.hpp:402
static Evaluation liquidPressure(const Evaluation &temperature, Scalar density)
The pressure of liquid in at a given density and temperature.
Definition: TabulatedComponent.hpp:388
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the gas .
Definition: TabulatedComponent.hpp:282
static Evaluation liquidHeatCapacity(const Evaluation &temperature, const Evaluation &pressure)
Specific isobaric heat capacity of the liquid .
Definition: TabulatedComponent.hpp:333
static Evaluation liquidInternalEnergy(const Evaluation &temperature, const Evaluation &pressure)
Specific internal energy of the liquid .
Definition: TabulatedComponent.hpp:360
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the tables.
Definition: TabulatedComponent.hpp:72
static Scalar criticalTemperature()
Returns the critical temperature in of the component.
Definition: TabulatedComponent.hpp:227
static Scalar criticalPressure()
Returns the critical pressure in of the component.
Definition: TabulatedComponent.hpp:233
static Scalar molarMass()
The molar mass in of the component.
Definition: TabulatedComponent.hpp:221
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition: TabulatedComponent.hpp:408
static Evaluation gasPressure(const Evaluation &temperature, Scalar density)
The pressure of gas in at a given density and temperature.
Definition: TabulatedComponent.hpp:370
static Evaluation liquidViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of liquid.
Definition: TabulatedComponent.hpp:478
static Scalar triplePressure()
Returns the pressure in at the component's triple point.
Definition: TabulatedComponent.hpp:257
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of gas at a given pressure and temperature .
Definition: TabulatedComponent.hpp:426
static Evaluation gasInternalEnergy(const Evaluation &temperature, const Evaluation &pressure)
Specific internal energy of the gas .
Definition: TabulatedComponent.hpp:350
static Scalar criticalVolume()
Returns the critical volume in of the component.
Definition: TabulatedComponent.hpp:245
static Evaluation gasHeatCapacity(const Evaluation &temperature, const Evaluation &pressure)
Specific isobaric heat capacity of the gas .
Definition: TabulatedComponent.hpp:316
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition: TabulatedComponent.hpp:414
static const char * name()
A human readable name for the component.
Definition: TabulatedComponent.hpp:215
static Scalar acentricFactor()
Returns the acentric factor of the component.
Definition: TabulatedComponent.hpp:239
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of liquid at a given pressure and temperature .
Definition: TabulatedComponent.hpp:444
static Evaluation vaporPressure(const Evaluation &temperature)
The vapor pressure in of the component at a given temperature.
Definition: TabulatedComponent.hpp:267
static Evaluation liquidThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
The thermal conductivity of liquid water .
Definition: TabulatedComponent.hpp:512
static Evaluation gasViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of gas.
Definition: TabulatedComponent.hpp:461
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of the liquid .
Definition: TabulatedComponent.hpp:299
static Scalar tripleTemperature()
Returns the temperature in at the component's triple point.
Definition: TabulatedComponent.hpp:251