version 3.10.0
Loading...
Searching...
No Matches
h2oheavyoil.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// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
12#ifndef DUMUX_H2O_HEAVYOIL_FLUID_SYSTEM_HH
13#define DUMUX_H2O_HEAVYOIL_FLUID_SYSTEM_HH
14
19
21
23
24#include <dumux/io/name.hh>
25
26namespace Dumux::FluidSystems {
27
33template <class Scalar,
34 class H2OType = Dumux::Components::TabulatedComponent<Dumux::Components::H2O<Scalar> > >
36 : public Base<Scalar, H2OHeavyOil<Scalar, H2OType> >
37{
38 using ThisType = H2OHeavyOil<Scalar, H2OType>;
39
40public:
42 using H2O = H2OType;
43
44
45 static const int numPhases = 3;
46 static const int numComponents = 2;
47
48 static const int wPhaseIdx = 0; // index of the water phase
49 static const int nPhaseIdx = 1; // index of the NAPL phase
50 static const int gPhaseIdx = 2; // index of the gas phase
51
52 static const int H2OIdx = 0;
53 static const int NAPLIdx = 1;
54
55 // export component indices to indicate the main component
56 // of the corresponding phase at atmospheric pressure 1 bar
57 // and room temperature 20°C:
58 static const int wCompIdx = H2OIdx;
59 static const int nCompIdx = NAPLIdx;
60
67 static void init()
68 {
69 init(/*tempMin=*/273.15,
70 /*tempMax=*/623.15,
71 /*numTemp=*/100,
72 /*pMin=*/0.0,
73 /*pMax=*/20e6,
74 /*numP=*/200);
75 }
76
88 static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
89 Scalar pressMin, Scalar pressMax, unsigned nPress)
90 {
91 if (H2O::isTabulated)
92 {
93 H2O::init(tempMin, tempMax, nTemp,
94 pressMin, pressMax, nPress);
95 }
96 }
97
102 static constexpr int getMainComponent(int phaseIdx)
103 {
104 // For the gas phase, choosing a main component appears to be
105 // rather arbitrary. Motivated by the fact that the thermal conductivity
106 // of the gas phase is set to the thermal conductivity of pure water,
107 // water is chosen for now.
108 if (phaseIdx == nPhaseIdx)
109 return nCompIdx;
110 else
111 return wCompIdx;
112 }
113
117 static constexpr bool isMiscible()
118 { return true; }
119
125 static constexpr bool isGas(int phaseIdx)
126 {
127 assert(0 <= phaseIdx && phaseIdx < numPhases);
128 return phaseIdx == gPhaseIdx;
129 }
130
137 static bool isIdealGas(int phaseIdx)
138 { return phaseIdx == gPhaseIdx && H2O::gasIsIdeal() && HeavyOil::gasIsIdeal(); }
139
154 static bool isIdealMixture(int phaseIdx)
155 {
156 assert(0 <= phaseIdx && phaseIdx < numPhases);
157 // we assume Henry's and Raoult's laws for the water phase and
158 // and no interaction between gas molecules of different
159 // components, so all phases are ideal mixtures!
160 return true;
161 }
162
172 static constexpr bool isCompressible(int phaseIdx)
173 {
174 assert(0 <= phaseIdx && phaseIdx < numPhases);
175 // gases are always compressible
176 if (phaseIdx == gPhaseIdx)
177 return true;
178 else if (phaseIdx == wPhaseIdx)
179 // the water component decides for the water phase...
180 return H2O::liquidIsCompressible();
181
182 // the NAPL component decides for the napl phase...
184 }
185
189 static std::string phaseName(int phaseIdx)
190 {
191 assert(0 <= phaseIdx && phaseIdx < numPhases);
192 switch (phaseIdx)
193 {
194 case wPhaseIdx: return IOName::aqueousPhase();
195 case nPhaseIdx: return IOName::naplPhase();
196 case gPhaseIdx: return IOName::gaseousPhase();
197 }
198 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
199 }
200
204 static std::string componentName(int compIdx)
205 {
206 switch (compIdx) {
207 case H2OIdx: return H2O::name();
208 case NAPLIdx: return HeavyOil::name();
209 };
210 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
211 }
212
216 static Scalar molarMass(int compIdx)
217 {
218 switch (compIdx) {
219 case H2OIdx: return H2O::molarMass();
220 case NAPLIdx: return HeavyOil::molarMass();
221 };
222 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
223 }
224
225 using Base<Scalar, ThisType>::density;
227 template <class FluidState>
228 static Scalar density(const FluidState &fluidState, int phaseIdx)
229 {
230 if (phaseIdx == wPhaseIdx) {
231 // See: doctoral thesis of Steffen Ochs 2007
232 // Steam injection into saturated porous media : process analysis including experimental and numerical investigations
233 // http://elib.uni-stuttgart.de/bitstream/11682/271/1/Diss_Ochs_OPUS.pdf
234
235 // This assumes each gas molecule displaces exactly one
236 // molecule in the liquid.
237
238 return H2O::liquidMolarDensity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx))
239 * (H2O::molarMass()*fluidState.moleFraction(wPhaseIdx, H2OIdx)
240 + HeavyOil::molarMass()*fluidState.moleFraction(wPhaseIdx, NAPLIdx));
241 }
242 else if (phaseIdx == nPhaseIdx) {
243 // assume pure NAPL for the NAPL phase
244 Scalar pressure = HeavyOil::liquidIsCompressible()?fluidState.pressure(phaseIdx):1e100;
245 return HeavyOil::liquidDensity(fluidState.temperature(phaseIdx), pressure);
246 }
247
248 assert (phaseIdx == gPhaseIdx);
249 Scalar pH2O =
250 fluidState.moleFraction(gPhaseIdx, H2OIdx) *
251 fluidState.pressure(gPhaseIdx);
252 Scalar pNAPL =
253 fluidState.moleFraction(gPhaseIdx, NAPLIdx) *
254 fluidState.pressure(gPhaseIdx);
255 return H2O::gasDensity(fluidState.temperature(phaseIdx), pH2O)
256 + HeavyOil::gasDensity(fluidState.temperature(phaseIdx), pNAPL);
257 }
258
259 using Base<Scalar, ThisType>::molarDensity;
261 template <class FluidState>
262 static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
263 {
264 Scalar temperature = fluidState.temperature(phaseIdx);
265 Scalar pressure = fluidState.pressure(phaseIdx);
266 if (phaseIdx == nPhaseIdx)
267 {
268 return HeavyOil::liquidMolarDensity(temperature, pressure);
269 }
270 else if (phaseIdx == wPhaseIdx)
271 { // This assumes each gas molecule displaces exactly one
272 // molecule in the liquid.
273 return H2O::liquidMolarDensity(temperature, pressure);
274 }
275 else
276 {
277 return H2O::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, H2OIdx))
278 + HeavyOil::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, NAPLIdx));
279 }
280 }
281
282 using Base<Scalar, ThisType>::viscosity;
284 template <class FluidState>
285 static Scalar viscosity(const FluidState &fluidState,
286 int phaseIdx)
287 {
288 if (phaseIdx == wPhaseIdx) {
289 // assume pure water viscosity
290 return H2O::liquidViscosity(fluidState.temperature(phaseIdx),
291 fluidState.pressure(phaseIdx));
292 }
293 else if (phaseIdx == nPhaseIdx) {
294 // assume pure NAPL viscosity
295 return HeavyOil::liquidViscosity(fluidState.temperature(phaseIdx),
296 fluidState.pressure(phaseIdx));
297 }
298
299 assert (phaseIdx == gPhaseIdx);
300
301 /* Wilke method. See:
302 *
303 * See: R. Reid, et al.: The Properties of Gases and Liquids,
304 * 4th edition, McGraw-Hill, 1987, 407-410
305 * 5th edition, McGraw-Hill, 20001, p. 9.21/22
306 *
307 * in this case, we use a simplified version in order to avoid
308 * computationally costly evaluation of sqrt and pow functions and
309 * divisions
310 * -- compare e.g. with Promo Class p. 32/33
311 */
312 const Scalar mu[numComponents] = {
313 h2oGasViscosityInMixture(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)),
314 HeavyOil::gasViscosity(fluidState.temperature(phaseIdx), HeavyOil::vaporPressure(fluidState.temperature(phaseIdx)))
315 };
316
317 return mu[H2OIdx]*fluidState.moleFraction(gPhaseIdx, H2OIdx)
318 + mu[NAPLIdx]*fluidState.moleFraction(gPhaseIdx, NAPLIdx);
319 }
320
336 template <class FluidState>
337 static Scalar binaryDiffusionCoefficient(const FluidState &fluidState,
338 int phaseIdx,
339 int compIIdx,
340 int compJIdx)
341
342 {
343 const Scalar T = fluidState.temperature(phaseIdx);
344 const Scalar p = fluidState.pressure(phaseIdx);
345
346 // liquid phase
347 if (phaseIdx == wPhaseIdx)
349
350 // gas phase
351 else if (phaseIdx == gPhaseIdx)
353 else
354 DUNE_THROW(Dune::InvalidStateException,
355 "Non-existent binary diffusion coefficient for phase index "
356 << phaseIdx);
357 }
358
359 using Base<Scalar, ThisType>::diffusionCoefficient;
368 template <class FluidState>
369 static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx)
370 {
371 // liquid phase
372 if (phaseIdx == wPhaseIdx)
373 return binaryDiffusionCoefficient(fluidState, phaseIdx, H2OIdx, NAPLIdx);
374 // gas phase
375 else if (phaseIdx == gPhaseIdx)
376 return binaryDiffusionCoefficient(fluidState, phaseIdx, NAPLIdx, H2OIdx);
377 else
378 DUNE_THROW(Dune::InvalidStateException,
379 "Non-existent diffusion coefficient for phase index "<< phaseIdx);
380 }
381
388 template <class FluidState>
389 static Scalar henryCoefficient(const FluidState &fluidState,
390 int phaseIdx,
391 int compIdx)
392 {
393 assert(0 <= phaseIdx && phaseIdx < numPhases);
394 assert(0 <= compIdx && compIdx < numComponents);
395
396 const Scalar T = fluidState.temperature(phaseIdx);
397
398 if (compIdx == NAPLIdx && phaseIdx == wPhaseIdx)
400
401 else if (phaseIdx == nPhaseIdx && compIdx == H2OIdx)
403
404 else
405 DUNE_THROW(Dune::InvalidStateException, "non-existent henry coefficient for phase index " << phaseIdx
406 << " and component index " << compIdx);
407 }
408
415 template <class FluidState>
416 static Scalar partialPressureGas(const FluidState &fluidState, int phaseIdx,
417 int compIdx)
418 {
419 assert(0 <= compIdx && compIdx < numComponents);
420
421 const Scalar T = fluidState.temperature(phaseIdx);
422 if (compIdx == NAPLIdx)
423 return HeavyOil::vaporPressure(T);
424 else if (compIdx == H2OIdx)
425 return H2O::vaporPressure(T);
426 else
427 DUNE_THROW(Dune::InvalidStateException, "non-existent component index " << compIdx);
428 }
429
436 template <class FluidState>
437 static Scalar inverseVaporPressureCurve(const FluidState &fluidState,
438 int phaseIdx,
439 int compIdx)
440 {
441 assert(0 <= compIdx && compIdx < numComponents);
442
443 const Scalar pressure = fluidState.pressure(phaseIdx);
444 if (compIdx == NAPLIdx)
445 return HeavyOil::vaporTemperature(pressure);
446 else if (compIdx == H2OIdx)
447 return H2O::vaporTemperature(pressure);
448 else
449 DUNE_THROW(Dune::InvalidStateException, "non-existent component index " << compIdx);
450 }
451
452 using Base<Scalar, ThisType>::enthalpy;
462 template <class FluidState>
463 static Scalar enthalpy(const FluidState &fluidState,
464 int phaseIdx)
465 {
466 if (phaseIdx == wPhaseIdx) {
467 return H2O::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
468 }
469 else if (phaseIdx == nPhaseIdx) {
470 return HeavyOil::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
471 }
472 else if (phaseIdx == gPhaseIdx) { // gas phase enthalpy depends strongly on composition
473 Scalar hgc = HeavyOil::gasEnthalpy(fluidState.temperature(phaseIdx),
474 fluidState.pressure(phaseIdx));
475 Scalar hgw = H2O::gasEnthalpy(fluidState.temperature(phaseIdx),
476 fluidState.pressure(phaseIdx));
477
478 Scalar result = 0;
479 result += hgw * fluidState.massFraction(gPhaseIdx, H2OIdx);
480 result += hgc * fluidState.massFraction(gPhaseIdx, NAPLIdx);
481
482 return result;
483 }
484 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
485 }
486
493 template <class FluidState>
494 static Scalar componentEnthalpy(const FluidState& fluidState, int phaseIdx, int componentIdx)
495 {
496 const Scalar T = fluidState.temperature(phaseIdx);
497 const Scalar p = fluidState.pressure(phaseIdx);
498
499 if (phaseIdx == wPhaseIdx)
500 {
501 if (componentIdx == H2OIdx)
502 return H2O::liquidEnthalpy(T, p);
503 else if (componentIdx == NAPLIdx)
504 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for NAPL in water is not implemented.");
505 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
506 }
507 else if (phaseIdx == nPhaseIdx)
508 {
509 if (componentIdx == H2OIdx)
510 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for water in NAPL is not implemented.");
511 else if (componentIdx == NAPLIdx)
512 return HeavyOil::liquidEnthalpy(T, p);
513 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
514 }
515 else if (phaseIdx == gPhaseIdx)
516 {
517 if (componentIdx == H2OIdx)
518 return H2O::gasEnthalpy(T, p);
519 else if (componentIdx == NAPLIdx)
520 return HeavyOil::gasEnthalpy(T, p);
521 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
522 }
523 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
524 }
525
526 using Base<Scalar, ThisType>::heatCapacity;
528 template <class FluidState>
529 static Scalar heatCapacity(const FluidState &fluidState,
530 int phaseIdx)
531 {
532 DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2ONAPL::heatCapacity()");
533 }
534
535 using Base<Scalar, ThisType>::thermalConductivity;
544 template <class FluidState>
545 static Scalar thermalConductivity(const FluidState &fluidState,
546 int phaseIdx)
547 {
548 const Scalar temperature = fluidState.temperature(phaseIdx) ;
549 const Scalar pressure = fluidState.pressure(phaseIdx);
550 if (phaseIdx == wPhaseIdx)
551 {
552 return H2O::liquidThermalConductivity(temperature, pressure);
553 }
554 else if (phaseIdx == nPhaseIdx)
555 {
556 return HeavyOil::liquidThermalConductivity(temperature, pressure);
557 }
558 else if (phaseIdx == gPhaseIdx)
559 {
560 return H2O::gasThermalConductivity(temperature, pressure);
561 }
562 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
563 }
564
565};
566
567} // end namespace Dumux::FluidSystems
568
569#endif
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient [m^2/s] for heavy oil in liquid water.
Definition h2o_heavyoil.hh:77
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient [m^2/s] for molecular water and heavy oil.
Definition h2o_heavyoil.hh:64
static Scalar henryOilInWater(Scalar temperature)
Henry coefficient for heavy oil in liquid water.
Definition h2o_heavyoil.hh:34
static Scalar henryWaterInOil(Scalar temperature)
Henry coefficient for water in liquid heavy oil.
Definition h2o_heavyoil.hh:49
Properties of the component heavyoil.
Definition heavyoil.hh:38
static Scalar liquidViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of pure heavyoil.
Definition heavyoil.hh:422
static constexpr Scalar molarMass()
The molar mass in of heavyoil.
Definition heavyoil.hh:51
static Scalar vaporPressure(Scalar temperature)
The saturation vapor pressure in of.
Definition heavyoil.hh:219
static constexpr bool gasIsIdeal()
Returns true if the gas phase is assumed to be ideal.
Definition heavyoil.hh:378
static Scalar liquidMolarDensity(Scalar temperature, Scalar pressure)
The molar density of pure heavyoil in at a given pressure and temperature.
Definition heavyoil.hh:366
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar density of pure heavyoil in , depending on pressure and temperature.
Definition heavyoil.hh:339
static Scalar gasViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of heavyoil vapor.
Definition heavyoil.hh:393
static Scalar vaporTemperature(Scalar pressure)
Definition heavyoil.hh:231
static Scalar liquidThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of heavy oil.
Definition heavyoil.hh:463
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The (ideal) gas density of heavyoil vapor at a given temperature and pressure .
Definition heavyoil.hh:328
static constexpr bool liquidIsCompressible()
Returns true if the liquid phase is assumed to be compressible.
Definition heavyoil.hh:384
static Scalar liquidEnthalpy(const Scalar temperature, const Scalar pressure)
Specific enthalpy of liquid heavyoil .
Definition heavyoil.hh:249
static std::string name()
A human readable name for heavyoil.
Definition heavyoil.hh:45
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of heavyoil vapor .
Definition heavyoil.hh:317
static Scalar liquidDensity(Scalar temperature, Scalar pressure)
The density of pure heavyoil at a given pressure and temperature .
Definition heavyoil.hh:348
Fluid system base class.
Definition fluidsystems/base.hh:32
ScalarType Scalar
Definition fluidsystems/base.hh:35
A compositional fluid system with water and heavy oil components in both the liquid and the gas phase...
Definition h2oheavyoil.hh:37
static Scalar henryCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Henry coefficients of a component in a phase.
Definition h2oheavyoil.hh:389
static constexpr bool isCompressible(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition h2oheavyoil.hh:172
static const int gPhaseIdx
Definition h2oheavyoil.hh:50
static std::string componentName(int compIdx)
Return the human readable name of a component (used in indices)
Definition h2oheavyoil.hh:204
static const int wCompIdx
Definition h2oheavyoil.hh:58
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
Calculate the molar density of a fluid phase.
Definition h2oheavyoil.hh:262
static bool isIdealMixture(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition h2oheavyoil.hh:154
static const int numPhases
Definition h2oheavyoil.hh:45
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition h2oheavyoil.hh:117
static Scalar viscosity(const FluidState &fluidState, int phaseIdx)
Calculate the dynamic viscosity of a fluid phase .
Definition h2oheavyoil.hh:285
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given all mole fractions in a phase, return the specific phase enthalpy .
Definition h2oheavyoil.hh:463
static const int H2OIdx
Definition h2oheavyoil.hh:52
static Scalar partialPressureGas(const FluidState &fluidState, int phaseIdx, int compIdx)
Partial pressures in the gas phase, taken from saturation vapor pressures.
Definition h2oheavyoil.hh:416
static Scalar componentEnthalpy(const FluidState &fluidState, int phaseIdx, int componentIdx)
Returns the specific enthalpy of a component in a specific phase.
Definition h2oheavyoil.hh:494
static const int nPhaseIdx
Definition h2oheavyoil.hh:49
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition h2oheavyoil.hh:545
H2OType H2O
Definition h2oheavyoil.hh:42
static std::string phaseName(int phaseIdx)
Return the human readable name of a phase (used in indices)
Definition h2oheavyoil.hh:189
static const int numComponents
Definition h2oheavyoil.hh:46
static constexpr int getMainComponent(int phaseIdx)
Get the main component of a given phase.
Definition h2oheavyoil.hh:102
static constexpr bool isGas(int phaseIdx)
Return whether a phase is gaseous.
Definition h2oheavyoil.hh:125
Dumux::Components::HeavyOil< Scalar > HeavyOil
Definition h2oheavyoil.hh:41
static const int NAPLIdx
Definition h2oheavyoil.hh:53
static Scalar molarMass(int compIdx)
Return the molar mass of a component in [kg/mol].
Definition h2oheavyoil.hh:216
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the fluid system's static parameters using problem specific temperature and pressure range...
Definition h2oheavyoil.hh:88
static Scalar density(const FluidState &fluidState, int phaseIdx)
Calculate the density of a fluid phase.
Definition h2oheavyoil.hh:228
static Scalar inverseVaporPressureCurve(const FluidState &fluidState, int phaseIdx, int compIdx)
Inverse vapor pressures, taken from inverse saturation vapor pressures.
Definition h2oheavyoil.hh:437
static const int nCompIdx
Definition h2oheavyoil.hh:59
static bool isIdealGas(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition h2oheavyoil.hh:137
static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx)
Calculate the molecular diffusion coefficient for a component in a fluid phase .
Definition h2oheavyoil.hh:369
static const int wPhaseIdx
Definition h2oheavyoil.hh:48
static Scalar binaryDiffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIIdx, int compJIdx)
Given a phase's composition, temperature and pressure, return the binary diffusion coefficient for c...
Definition h2oheavyoil.hh:337
static void init()
Initialize the fluid system's static parameters generically.
Definition h2oheavyoil.hh:67
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity of a fluid phase .
Definition h2oheavyoil.hh:529
Fluid system base class.
Material properties of pure water .
Binary coefficients for water and heavy oil.
Properties of the component heavyoil.
Relations valid for an ideal gas.
A collection of input/output field names for common physical quantities.
Definition h2o.hh:901
Scalar h2oGasViscosityInMixture(Scalar temperature, Scalar pressure)
The dynamic viscosity of steam in a gas mixture.
Definition h2o.hh:964
std::string gaseousPhase() noexcept
I/O name of gaseous phase.
Definition name.hh:111
std::string naplPhase() noexcept
I/O name of napl phase.
Definition name.hh:119
std::string aqueousPhase() noexcept
I/O name of aqueous phase.
Definition name.hh:115
Tabulates all thermodynamic properties of a given untabulated chemical species.