12#ifndef DUMUX_H2O_AIR_XYLENE_FLUID_SYSTEM_HH
13#define DUMUX_H2O_AIR_XYLENE_FLUID_SYSTEM_HH
39template <
class Scalar,
40 class H2OType = Components::TabulatedComponent<Components::H2O<Scalar> > >
42 :
public Base<Scalar, H2OAirXylene<Scalar, H2OType> >
101 H2O::init(tempMin, tempMax, nTemp,
102 pressMin, pressMax, nPress);
117 static constexpr bool isGas(
int phaseIdx)
119 assert(0 <= phaseIdx && phaseIdx <
numPhases);
148 assert(0 <= phaseIdx && phaseIdx <
numPhases);
166 assert(0 <= phaseIdx && phaseIdx <
numPhases);
172 return H2O::liquidIsCompressible();
184 assert(0 <= phaseIdx && phaseIdx <
numPhases);
191 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
201 case H2OIdx:
return H2O::name();
205 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << compIdx);
215 case H2OIdx:
return H2O::molarMass();
219 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << compIdx);
235 template <
class Flu
idState>
241 return H2O::liquidMolarDensity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx))
262 return H2O::gasDensity(fluidState.temperature(phaseIdx), pH2O)
269 template <
class Flu
idState>
272 Scalar temperature = fluidState.temperature(phaseIdx);
273 Scalar pressure = fluidState.pressure(phaseIdx);
281 return H2O::liquidMolarDensity(temperature, pressure);
285 return H2O::gasMolarDensity(temperature, fluidState.partialPressure(
gPhaseIdx,
H2OIdx))
299 template <
class Flu
idState>
305 return H2O::liquidViscosity(fluidState.temperature(phaseIdx),
306 fluidState.pressure(phaseIdx));
311 fluidState.pressure(phaseIdx));
330 Air::gasViscosity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)),
360 muResult = (xAW*muAW)/(xAW+fluidState.moleFraction(
gPhaseIdx,
NAPLIdx)*phiAWC)
369 template <
class Flu
idState>
374 Scalar temperature = fluidState.temperature(phaseIdx);
375 Scalar pressure = fluidState.pressure(phaseIdx);
386 return (1.- xgw)/(xga/diffAW + xgc/diffWC);
388 return (1.- xgc)/(xgw/diffWC + xga/diffAC);
390 DUNE_THROW(Dune::InvalidStateException,
391 "Diffusivity of Air in the gas phase "
392 "is constraint by sum of diffusive fluxes = 0 !\n");
406 diffCont = (1.- xww)/(xwa/diffAWl + xwc/diffWCl);
409 diffCont = (1.- xwc)/(xww/diffWCl + xwa/diffACl);
412 DUNE_THROW(Dune::InvalidStateException,
413 "Diffusivity of water in the water phase "
414 "is constraint by sum of diffusive fluxes = 0 !\n");
424 template <
class Flu
idState>
430 DUNE_THROW(Dune::NotImplemented,
"FluidSystems::H2OAirXylene::binaryDiffusionCoefficient()");
447 template <
class Flu
idState>
452 assert(0 <= phaseIdx && phaseIdx <
numPhases);
455 Scalar T = fluidState.temperature(phaseIdx);
456 Scalar p = fluidState.pressure(phaseIdx);
460 return H2O::vaporPressure(T)/p;
461 else if (compIdx ==
AirIdx)
476 else if (compIdx ==
AirIdx)
478 else if (compIdx ==
H2OIdx)
488 template <
class Flu
idState>
493 DUNE_THROW(Dune::NotImplemented,
"FluidSystems::H2OAirXylene::kelvinVaporPressure()");
506 template <
class Flu
idState>
511 return H2O::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
518 fluidState.pressure(phaseIdx));
519 Scalar hgw = H2O::gasEnthalpy(fluidState.temperature(phaseIdx),
520 fluidState.pressure(phaseIdx));
531 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
540 template <
class Flu
idState>
543 const Scalar T = fluidState.temperature(phaseIdx);
544 const Scalar p = fluidState.pressure(phaseIdx);
548 if (componentIdx ==
H2OIdx)
549 return H2O::liquidEnthalpy(T, p);
550 else if (componentIdx ==
NAPLIdx)
551 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for NAPL in water is not implemented.");
552 else if (componentIdx ==
AirIdx)
553 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for Air in water is not implemented.");
554 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << componentIdx);
558 if (componentIdx ==
H2OIdx)
559 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for water in NAPL is not implemented.");
560 else if (componentIdx ==
NAPLIdx)
562 else if (componentIdx ==
AirIdx)
563 DUNE_THROW(Dune::NotImplemented,
"The component enthalpy for air in NAPL is not implemented.");
564 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << componentIdx);
568 if (componentIdx ==
H2OIdx)
569 return H2O::gasEnthalpy(T, p);
570 else if (componentIdx ==
NAPLIdx)
572 else if (componentIdx ==
AirIdx)
574 DUNE_THROW(Dune::InvalidStateException,
"Invalid component index " << componentIdx);
576 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
581 template <
class Flu
idState>
585 DUNE_THROW(Dune::NotImplemented,
"FluidSystems::H2OAirXylene::heatCapacity()");
590 template <
class Flu
idState>
594 const Scalar temperature = fluidState.temperature(phaseIdx) ;
595 const Scalar pressure = fluidState.pressure(phaseIdx);
598 return H2O::liquidThermalConductivity(temperature, pressure);
608 DUNE_THROW(Dune::InvalidStateException,
"Invalid phase index " << phaseIdx);
614 Scalar rholH2O = H2O::liquidDensity(T, pw);
615 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:35
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for molecular water and air.
Definition h2o_air.hh:55
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for molecular nitrogen in liquid water.
Definition h2o_air.hh:88
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:33
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of Air at a given pressure and temperature.
Definition air.hh:71
static constexpr Scalar molarMass()
The molar mass in of Air.
Definition air.hh:48
static Scalar gasViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of Air at a given pressure and temperature.
Definition air.hh:179
static Scalar gasThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of air.
Definition air.hh:337
static constexpr bool gasIsIdeal()
Returns true, the gas phase is assumed to be ideal.
Definition air.hh:95
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of Air with 273.15 as basis.
Definition air.hh:262
static std::string name()
A human readable name for Air.
Definition air.hh:40
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar density of air in , depending on pressure and temperature.
Definition air.hh:83
Properties of xylene.
Definition xylene.hh:36
static Scalar gasViscosity(Scalar temp, Scalar pressure)
The dynamic viscosity of xylene vapor.
Definition xylene.hh:315
static std::string name()
A human readable name for the xylene.
Definition xylene.hh:44
static Scalar liquidDensity(Scalar temperature, Scalar pressure)
The density of pure xylene at a given pressure and temperature .
Definition xylene.hh:286
static Scalar liquidEnthalpy(const Scalar temperature, const Scalar pressure)
Specific enthalpy of liquid xylene .
Definition xylene.hh:154
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar gas density of xylene gas at a given pressure and temperature.
Definition xylene.hh:249
static constexpr bool gasIsIdeal()
Returns true if the gas phase is assumed to be ideal.
Definition xylene.hh:300
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:93
static Scalar liquidThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of xylene.
Definition xylene.hh:371
static constexpr bool liquidIsCompressible()
Returns true if the liquid phase is assumed to be compressible.
Definition xylene.hh:306
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of xylene vapor .
Definition xylene.hh:225
static constexpr Scalar molarMass()
The molar mass in of xylene.
Definition xylene.hh:50
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of xylene gas at a given pressure and temperature.
Definition xylene.hh:236
static Scalar liquidMolarDensity(Scalar temp, Scalar pressure)
The molar liquid density of pure xylene at a given pressure and temperature .
Definition xylene.hh:261
static Scalar liquidViscosity(Scalar temp, Scalar pressure)
The dynamic viscosity of pure xylene.
Definition xylene.hh:343
Fluid system base class.
Definition fluidsystems/base.hh:32
ScalarType Scalar
Definition fluidsystems/base.hh:35
A three-phase fluid system featuring gas, NAPL and water as phases and distilled water and air (Pseu...
Definition h2oairxylene.hh:43
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:370
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:236
static constexpr bool isGas(int phaseIdx)
Return whether a phase is gaseous.
Definition h2oairxylene.hh:117
static const int gCompIdx
Definition h2oairxylene.hh:67
static Scalar viscosity(const FluidState &fluidState, int phaseIdx)
Return the viscosity of a phase .
Definition h2oairxylene.hh:300
static void init()
Initialize the fluid system's static parameters generically.
Definition h2oairxylene.hh:75
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity of a fluid phase .
Definition h2oairxylene.hh:582
static const int nCompIdx
Definition h2oairxylene.hh:66
static Scalar molarMass(int compIdx)
Return the molar mass of a component in .
Definition h2oairxylene.hh:212
static const int AirIdx
Definition h2oairxylene.hh:60
static const int gPhaseIdx
Definition h2oairxylene.hh:56
static std::string componentName(int compIdx)
Return the human readable name of a component (used in indices)
Definition h2oairxylene.hh:198
static std::string phaseName(int phaseIdx)
Return the human readable name of a phase (used in indices)
Definition h2oairxylene.hh:182
static const int numComponents
Definition h2oairxylene.hh:52
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:425
static const int H2OIdx
Definition h2oairxylene.hh:58
static const int wCompIdx
Definition h2oairxylene.hh:65
static constexpr bool isCompressible(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition h2oairxylene.hh:164
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition h2oairxylene.hh:109
Components::Xylene< Scalar > NAPL
Definition h2oairxylene.hh:48
static bool isIdealGas(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition h2oairxylene.hh:129
H2OType H2O
Definition h2oairxylene.hh:47
static const int nPhaseIdx
Definition h2oairxylene.hh:55
static const int numPhases
Definition h2oairxylene.hh:51
static const int NAPLIdx
Definition h2oairxylene.hh:59
static Scalar fugacityCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Returns the fugacity coefficient of a component in a phase.
Definition h2oairxylene.hh:448
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
Calculate the molar density of a fluid phase.
Definition h2oairxylene.hh:270
static Scalar componentEnthalpy(const FluidState &fluidState, int phaseIdx, int componentIdx)
Returns the specific enthalpy of a component in a specific phase.
Definition h2oairxylene.hh:541
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given all mole fractions in a phase, return the specific phase enthalpy .
Definition h2oairxylene.hh:507
static const int wPhaseIdx
Definition h2oairxylene.hh:54
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition h2oairxylene.hh:591
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:96
static bool isIdealMixture(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition h2oairxylene.hh:146
static Scalar kelvinVaporPressure(const FluidState &fluidState, const int phaseIdx, const int compIdx)
Definition h2oairxylene.hh:489
Dumux::Components::Air< Scalar > Air
Definition h2oairxylene.hh:49
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.
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.