12#ifndef DUMUX_H2O_AIR_XYLENE_FLUID_SYSTEM_HH
13#define DUMUX_H2O_AIR_XYLENE_FLUID_SYSTEM_HH
30namespace FluidSystems {
40template <
class Scalar,
41 class H2OType = Components::TabulatedComponent<Components::H2O<Scalar> > >
43 :
public Base<Scalar, H2OAirXylene<Scalar, H2OType> >
100 if (H2O::isTabulated)
102 H2O::init(tempMin, tempMax, nTemp,
103 pressMin, pressMax, nPress);
118 static constexpr bool isGas(
int phaseIdx)
120 assert(0 <= phaseIdx && phaseIdx <
numPhases);
149 assert(0 <= phaseIdx && phaseIdx <
numPhases);
167 assert(0 <= phaseIdx && phaseIdx <
numPhases);
173 return H2O::liquidIsCompressible();
185 assert(0 <= phaseIdx && phaseIdx <
numPhases);
192 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
202 case H2OIdx:
return H2O::name();
206 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << compIdx);
216 case H2OIdx:
return H2O::molarMass();
220 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << compIdx);
236 template <
class Flu
idState>
242 return H2O::liquidMolarDensity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx))
263 return H2O::gasDensity(fluidState.temperature(phaseIdx), pH2O)
270 template <
class Flu
idState>
273 Scalar temperature = fluidState.temperature(phaseIdx);
274 Scalar pressure = fluidState.pressure(phaseIdx);
282 return H2O::liquidMolarDensity(temperature, pressure);
286 return H2O::gasMolarDensity(temperature, fluidState.partialPressure(
gPhaseIdx,
H2OIdx))
300 template <
class Flu
idState>
306 return H2O::liquidViscosity(fluidState.temperature(phaseIdx),
307 fluidState.pressure(phaseIdx));
312 fluidState.pressure(phaseIdx));
330 h2oGasViscosityInMixture(fluidState.temperature(phaseIdx),fluidState.pressure(phaseIdx)),
331 Air::gasViscosity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)),
361 muResult = (xAW*muAW)/(xAW+fluidState.moleFraction(
gPhaseIdx,
NAPLIdx)*phiAWC)
370 template <
class Flu
idState>
377 Scalar temperature = fluidState.temperature(phaseIdx);
378 Scalar pressure = fluidState.pressure(phaseIdx);
389 return (1.- xgw)/(xga/diffAW + xgc/diffWC);
391 return (1.- xgc)/(xgw/diffWC + xga/diffAC);
393 DUNE_THROW(Dune::InvalidStateException,
394 "Diffusivity of Air in the gas phase "
395 "is constraint by sum of diffusive fluxes = 0 !\n");
407 diffCont = (1.- xww)/(xwa/diffAWl + xwc/diffWCl);
410 diffCont = (1.- xwc)/(xww/diffWCl + xwa/diffACl);
413 DUNE_THROW(Dune::InvalidStateException,
414 "Diffusivity of water in the water phase "
415 "is constraint by sum of diffusive fluxes = 0 !\n");
425 template <
class Flu
idState>
431 DUNE_THROW(Dune::NotImplemented,
"FluidSystems::H2OAirXylene::binaryDiffusionCoefficient()");
448 template <
class Flu
idState>
453 assert(0 <= phaseIdx && phaseIdx <
numPhases);
456 Scalar T = fluidState.temperature(phaseIdx);
457 Scalar p = fluidState.pressure(phaseIdx);
461 return H2O::vaporPressure(T)/p;
462 else if (compIdx ==
AirIdx)
477 else if (compIdx ==
AirIdx)
479 else if (compIdx ==
H2OIdx)
489 template <
class Flu
idState>
494 DUNE_THROW(Dune::NotImplemented,
"FluidSystems::H2OAirXylene::kelvinVaporPressure()");
507 template <
class Flu
idState>
512 return H2O::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
519 fluidState.pressure(phaseIdx));
520 Scalar hgw = H2O::gasEnthalpy(fluidState.temperature(phaseIdx),
521 fluidState.pressure(phaseIdx));
532 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
541 template <
class Flu
idState>
544 const Scalar T = fluidState.temperature(phaseIdx);
545 const Scalar p = fluidState.pressure(phaseIdx);
549 if (componentIdx ==
H2OIdx)
550 return H2O::liquidEnthalpy(T, p);
551 else if (componentIdx ==
NAPLIdx)
552 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for NAPL in water is not implemented.");
553 else if (componentIdx ==
AirIdx)
554 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for Air in water is not implemented.");
555 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << componentIdx);
559 if (componentIdx ==
H2OIdx)
560 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for water in NAPL is not implemented.");
561 else if (componentIdx ==
NAPLIdx)
563 else if (componentIdx ==
AirIdx)
564 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for air in NAPL is not implemented.");
565 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << componentIdx);
569 if (componentIdx ==
H2OIdx)
570 return H2O::gasEnthalpy(T, p);
571 else if (componentIdx ==
NAPLIdx)
573 else if (componentIdx ==
AirIdx)
575 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << componentIdx);
577 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
582 template <
class Flu
idState>
586 DUNE_THROW(Dune::NotImplemented,
"FluidSystems::H2OAirXylene::heatCapacity()");
591 template <
class Flu
idState>
595 const Scalar temperature = fluidState.temperature(phaseIdx) ;
596 const Scalar pressure = fluidState.pressure(phaseIdx);
599 return H2O::liquidThermalConductivity(temperature, pressure);
609 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
615 Scalar rholH2O = H2O::liquidDensity(T, pw);
616 Scalar clH2O = rholH2O/H2O::molarMass();
A simple class for the air fluid properties.
Binary coefficients for air and xylene.
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for air and xylene. method according to Wilke and Lee see W....
Definition air_xylene.hh:50
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for air and xylene in liquid water.
Definition air_xylene.hh:96
static Scalar henry(Scalar temperature)
Henry coefficient for air in liquid water.
Definition h2o_air.hh:36
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for molecular water and air.
Definition h2o_air.hh:56
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for molecular nitrogen in liquid water.
Definition h2o_air.hh:89
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for xylene in liquid water.
Definition h2o_xylene.hh:105
static Scalar henry(Scalar temperature)
Henry coefficient for xylene in liquid water.
Definition h2o_xylene.hh:40
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for molecular water and xylene.
Definition h2o_xylene.hh:57
A class for the air fluid properties.
Definition air.hh:35
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of Air at a given pressure and temperature.
Definition air.hh:73
static constexpr Scalar molarMass()
The molar mass in of Air.
Definition air.hh:50
static Scalar gasViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of Air at a given pressure and temperature.
Definition air.hh:181
static Scalar gasThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of air.
Definition air.hh:339
static constexpr bool gasIsIdeal()
Returns true, the gas phase is assumed to be ideal.
Definition air.hh:97
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of Air with 273.15 as basis.
Definition air.hh:264
static std::string name()
A human readable name for Air.
Definition air.hh:42
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar density of air in , depending on pressure and temperature.
Definition air.hh:85
Properties of xylene.
Definition xylene.hh:37
static Scalar gasViscosity(Scalar temp, Scalar pressure)
The dynamic viscosity of xylene vapor.
Definition xylene.hh:316
static std::string name()
A human readable name for the xylene.
Definition xylene.hh:45
static Scalar liquidDensity(Scalar temperature, Scalar pressure)
The density of pure xylene at a given pressure and temperature .
Definition xylene.hh:287
static Scalar liquidEnthalpy(const Scalar temperature, const Scalar pressure)
Specific enthalpy of liquid xylene .
Definition xylene.hh:155
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar gas density of xylene gas at a given pressure and temperature.
Definition xylene.hh:250
static constexpr bool gasIsIdeal()
Returns true if the gas phase is assumed to be ideal.
Definition xylene.hh:301
static Scalar vaporPressure(Scalar temperature)
The saturation vapor pressure in of pure xylene at a given temperature according to Antoine after Be...
Definition xylene.hh:94
static Scalar liquidThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of xylene.
Definition xylene.hh:372
static constexpr bool liquidIsCompressible()
Returns true if the liquid phase is assumed to be compressible.
Definition xylene.hh:307
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of xylene vapor .
Definition xylene.hh:226
static constexpr Scalar molarMass()
The molar mass in of xylene.
Definition xylene.hh:51
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of xylene gas at a given pressure and temperature.
Definition xylene.hh:237
static Scalar liquidMolarDensity(Scalar temp, Scalar pressure)
The molar liquid density of pure xylene at a given pressure and temperature .
Definition xylene.hh:262
static Scalar liquidViscosity(Scalar temp, Scalar pressure)
The dynamic viscosity of pure xylene.
Definition xylene.hh:344
Fluid system base class.
Definition fluidsystems/base.hh:33
ScalarType Scalar
export the scalar type
Definition fluidsystems/base.hh:36
A three-phase fluid system featuring gas, NAPL and water as phases and distilled water and air (Pseu...
Definition h2oairxylene.hh:44
static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase .
Definition h2oairxylene.hh:371
static Scalar density(const FluidState &fluidState, int phaseIdx)
Given a phase's composition, temperature, pressure, and the partial pressures of all components,...
Definition h2oairxylene.hh:237
static constexpr bool isGas(int phaseIdx)
Return whether a phase is gaseous.
Definition h2oairxylene.hh:118
static const int gCompIdx
Definition h2oairxylene.hh:68
static Scalar viscosity(const FluidState &fluidState, int phaseIdx)
Return the viscosity of a phase .
Definition h2oairxylene.hh:301
static void init()
Initialize the fluid system's static parameters generically.
Definition h2oairxylene.hh:76
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity of a fluid phase .
Definition h2oairxylene.hh:583
static const int nCompIdx
Definition h2oairxylene.hh:67
static Scalar molarMass(int compIdx)
Return the molar mass of a component in .
Definition h2oairxylene.hh:213
static const int AirIdx
Definition h2oairxylene.hh:61
static const int gPhaseIdx
Definition h2oairxylene.hh:57
static std::string componentName(int compIdx)
Return the human readable name of a component (used in indices)
Definition h2oairxylene.hh:199
static std::string phaseName(int phaseIdx)
Return the human readable name of a phase (used in indices)
Definition h2oairxylene.hh:183
static const int numComponents
Definition h2oairxylene.hh:53
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 h2oairxylene.hh:426
static const int H2OIdx
Definition h2oairxylene.hh:59
static const int wCompIdx
Definition h2oairxylene.hh:66
static constexpr bool isCompressible(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition h2oairxylene.hh:165
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition h2oairxylene.hh:110
static bool isIdealGas(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition h2oairxylene.hh:130
H2OType H2O
Definition h2oairxylene.hh:48
static const int nPhaseIdx
Definition h2oairxylene.hh:56
static const int numPhases
Definition h2oairxylene.hh:52
static const int NAPLIdx
Definition h2oairxylene.hh:60
static Scalar fugacityCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Returns the fugacity coefficient of a component in a phase.
Definition h2oairxylene.hh:449
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
Calculate the molar density of a fluid phase.
Definition h2oairxylene.hh:271
static Scalar componentEnthalpy(const FluidState &fluidState, int phaseIdx, int componentIdx)
Returns the specific enthalpy of a component in a specific phase.
Definition h2oairxylene.hh:542
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given all mole fractions in a phase, return the specific phase enthalpy .
Definition h2oairxylene.hh:508
static const int wPhaseIdx
Definition h2oairxylene.hh:55
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition h2oairxylene.hh:592
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 h2oairxylene.hh:97
static bool isIdealMixture(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition h2oairxylene.hh:147
static Scalar kelvinVaporPressure(const FluidState &fluidState, const int phaseIdx, const int compIdx)
Definition h2oairxylene.hh:490
Material properties of pure water .
Binary coefficients for water and air.
Binary coefficients for water and xylene.
Relations valid for an ideal gas.
A collection of input/output field names for common physical quantities.
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.