13#ifndef DUMUX_H2O_AIR_SYSTEM_HH 
   14#define DUMUX_H2O_AIR_SYSTEM_HH 
   37template<
bool fastButSimplifiedRelations = false>
 
   56template <
class Scalar,
 
   58          class Policy = H2OAirDefaultPolicy<>,
 
   59          bool useKelvinVaporPressure = 
false>
 
   61: 
public Base<Scalar, H2OAir<Scalar, H2Otype, Policy> >
 
   92        assert(0 <= phaseIdx && phaseIdx < 
numPhases);
 
   98        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  112    static constexpr bool isGas(
int phaseIdx)
 
  114        assert(0 <= phaseIdx && phaseIdx < 
numPhases);
 
 
  134        assert(0 <= phaseIdx && phaseIdx < 
numPhases);
 
 
  152        assert(0 <= phaseIdx && phaseIdx < 
numPhases);
 
  157        return H2O::liquidIsCompressible();
 
 
  170            return H2O::liquidViscosityIsConstant();
 
  172        else if (phaseIdx == 
gasPhaseIdx && Policy::useAirViscosityAsGasMixtureViscosity())
 
 
  187        assert(0 <= phaseIdx && phaseIdx < 
numPhases);
 
 
  207            case H2OIdx: 
return H2O::name();
 
  210        DUNE_THROW(Dune::InvalidStateException, 
"Invalid component index " << compIdx);
 
 
  222            case H2OIdx: 
return H2O::molarMass();
 
  225        DUNE_THROW(Dune::InvalidStateException, 
"Invalid component index " << compIdx);
 
 
  235        static const Scalar TCrit[] = {
 
  236            H2O::criticalTemperature(),
 
  241        return TCrit[compIdx];
 
 
  251        static const Scalar pCrit[] = {
 
  252            H2O::criticalPressure(),
 
  257        return pCrit[compIdx];
 
 
  266    template <
class Flu
idState>
 
  271            const auto t = fluidState.temperature(
H2OIdx);
 
  273            if constexpr (!useKelvinVaporPressure)
 
  274                return H2O::vaporPressure(t);
 
  277                const auto pc =  (fluidState.wettingPhase() == (
int) 
H2OIdx)
 
  278                                 ? fluidState.pressure(
AirIdx)-fluidState.pressure(
H2OIdx)
 
  279                                 : fluidState.pressure(
H2OIdx)-fluidState.pressure(
AirIdx);
 
  285        else if (compIdx == 
AirIdx)
 
  287            DUNE_THROW(Dune::NotImplemented, 
"Air::vaporPressure(t)");
 
  289            DUNE_THROW(Dune::NotImplemented, 
"Invalid component index " << compIdx);
 
 
  299        DUNE_THROW(Dune::NotImplemented,
 
  300                   "H2OAirFluidSystem::criticalMolarVolume()");
 
 
  310        static const Scalar accFac[] = {
 
  311            H2O::acentricFactor(),
 
  312            Air::acentricFactor()
 
  316        return accFac[compIdx];
 
 
  353        std::cout << 
"The H2O-air fluid system was configured with the following policy:\n";
 
  354        std::cout << 
" - use H2O density as liquid mixture density: " << std::boolalpha << Policy::useH2ODensityAsLiquidMixtureDensity() << 
"\n";
 
  355        std::cout << 
" - use ideal gas density: " << std::boolalpha << Policy::useIdealGasDensity() << 
"\n";
 
  356        std::cout << 
" - use air viscosity as gas mixture viscosity: " << std::boolalpha << Policy::useAirViscosityAsGasMixtureViscosity() << std::endl;
 
  358        if (H2O::isTabulated)
 
  360            H2O::init(tempMin, tempMax, nTemp,
 
  361                      pressMin, pressMax, nPress);
 
 
  379    template <
class Flu
idState>
 
  383        assert(0 <= phaseIdx  && phaseIdx < 
numPhases);
 
  385        const Scalar T = fluidState.temperature(phaseIdx);
 
  386        const Scalar p = fluidState.pressure(phaseIdx);
 
  390            if (Policy::useH2ODensityAsLiquidMixtureDensity())
 
  392                return H2O::liquidDensity(T, p);
 
  398                return H2O::liquidMolarDensity(T, p)
 
  405            if (Policy::useIdealGasDensity())
 
  415        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  420    template <
class Flu
idState>
 
  423        assert(0 <= phaseIdx  && phaseIdx < 
numPhases);
 
  425        const Scalar T = fluidState.temperature(phaseIdx);
 
  426        const Scalar p = fluidState.pressure(phaseIdx);
 
  432            return H2O::liquidMolarDensity(T, p);
 
  436            if (Policy::useIdealGasDensity())
 
  440            return H2O::gasMolarDensity(T, fluidState.partialPressure(
phase1Idx, 
H2OIdx))
 
  443        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  459    template <
class Flu
idState>
 
  463        assert(0 <= phaseIdx  && phaseIdx < 
numPhases);
 
  465        Scalar T = fluidState.temperature(phaseIdx);
 
  466        Scalar p = fluidState.pressure(phaseIdx);
 
  471            return H2O::liquidViscosity(T, p);
 
  475            if(Policy::useAirViscosityAsGasMixtureViscosity()){
 
  500                        Scalar phiIJ = 1 + sqrt(mu[i]/mu[j] * sqrt(M[j]/M[i]));
 
  501                        phiIJ = phiIJ * phiIJ / sqrt(8*(1 + M[i]/M[j]));
 
  502                        divisor += fluidState.moleFraction(phaseIdx, j)*phiIJ;
 
  504                    muResult += fluidState.moleFraction(phaseIdx, i)*mu[i] / divisor;
 
  509        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  531    template <
class Flu
idState>
 
  536        assert(0 <= phaseIdx  && phaseIdx < 
numPhases);
 
  539        Scalar T = fluidState.temperature(phaseIdx);
 
  540        Scalar p = fluidState.pressure(phaseIdx);
 
 
  559    template <
class Flu
idState>
 
  563               / H2O::vaporPressure(fluidState.temperature(
gasPhaseIdx));
 
 
  568    template <
class Flu
idState>
 
  573        DUNE_THROW(Dune::NotImplemented, 
"FluidSystems::H2OAir::diffusionCoefficient()");
 
 
  578    template <
class Flu
idState>
 
  585        if (compIIdx > compJIdx)
 
  586            swap(compIIdx, compJIdx);
 
  588        Scalar T = fluidState.temperature(phaseIdx);
 
  589        Scalar p = fluidState.pressure(phaseIdx);
 
  597                DUNE_THROW(Dune::InvalidStateException,
 
  598                           "Binary diffusion coefficient of components " 
  599                            << compIIdx << 
" and " << compJIdx
 
  600                            << 
" in phase " << phaseIdx << 
" is undefined!\n");
 
  609                DUNE_THROW(Dune::InvalidStateException,
 
  610                           "Binary diffusion coefficient of components " 
  611                           << compIIdx << 
" and " << compJIdx
 
  612                           << 
" in phase " << phaseIdx << 
" is undefined!\n");
 
  615        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  636    template <
class Flu
idState>
 
  640        const Scalar T = fluidState.temperature(phaseIdx);
 
  641        const Scalar p = fluidState.pressure(phaseIdx);
 
  644            return H2O::liquidEnthalpy(T, p);
 
  650        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  660    template <
class Flu
idState>
 
  665        const Scalar T = fluidState.temperature(phaseIdx);
 
  666        const Scalar p = fluidState.pressure(phaseIdx);
 
  671            return H2O::liquidEnthalpy(T, p);
 
  675            if (componentIdx == 
H2OIdx)
 
  677                return H2O::gasEnthalpy(T, p);
 
  679            else if (componentIdx == 
AirIdx)
 
  683            DUNE_THROW(Dune::InvalidStateException, 
"Invalid component index " << componentIdx);
 
  685        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  698    template <
class Flu
idState>
 
  702        assert(0 <= phaseIdx  && phaseIdx < 
numPhases);
 
  704        const Scalar temperature  = fluidState.temperature(phaseIdx) ;
 
  705        const Scalar pressure = fluidState.pressure(phaseIdx);
 
  708            return H2O::liquidThermalConductivity(temperature, pressure);
 
  715            DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  729    template <
class Flu
idState>
 
  733        const Scalar temperature  = fluidState.temperature(phaseIdx);
 
  734        const Scalar pressure = fluidState.pressure(phaseIdx);
 
  738            return H2O::liquidHeatCapacity(temperature, pressure);
 
  743                   + H2O::gasHeatCapacity(temperature, pressure) * fluidState.moleFraction(
gasPhaseIdx, 
H2OIdx);
 
  746            DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
 
A simple class for the air fluid properties.
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
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 const Scalar gasHeatCapacity(Scalar temperature, Scalar pressure)
Specific isobaric heat capacity  of pure air.
Definition air.hh:299
static Scalar criticalPressure()
Returns the critical pressure  of Air.
Definition air.hh:60
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 gasViscosityIsConstant()
Returns true if the gas phase viscosity is constant.
Definition air.hh:101
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 Scalar criticalTemperature()
Returns the critical temperature  of Air.
Definition air.hh:54
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
Tabulates all thermodynamic properties of a given component.
Definition tabulatedcomponent.hh:672
static constexpr Scalar R
The ideal gas constant .
Definition constants.hh:32
Fluid system base class.
Definition fluidsystems/base.hh:32
ScalarType Scalar
Definition fluidsystems/base.hh:35
A compositional two-phase fluid system with water and air as components in both, the liquid and the g...
Definition h2oair.hh:62
static Scalar viscosity(const FluidState &fluidState, int phaseIdx)
Calculate the dynamic viscosity of a fluid phase .
Definition h2oair.hh:460
static constexpr bool isGas(int phaseIdx)
Return whether a phase is gaseous.
Definition h2oair.hh:112
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given a phase's composition, temperature and pressure, return its specific enthalpy .
Definition h2oair.hh:637
static Scalar criticalMolarVolume(int compIdx)
Molar volume of a component at the critical point .
Definition h2oair.hh:297
static Scalar vaporPressure(const FluidState &fluidState, int compIdx)
Vapor pressure of a component .
Definition h2oair.hh:267
static constexpr int comp1Idx
Definition h2oair.hh:81
static constexpr bool isIdealGas(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition h2oair.hh:185
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 h2oair.hh:350
static std::string componentName(int compIdx)
Return the human readable name of a component.
Definition h2oair.hh:203
static constexpr int numComponents
Definition h2oair.hh:71
static constexpr int liquidPhaseIdx
Definition h2oair.hh:73
static constexpr int gasPhaseIdx
Definition h2oair.hh:74
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity of a fluid phase. .
Definition h2oair.hh:730
static std::string phaseName(int phaseIdx)
Return the human readable name of a phase.
Definition h2oair.hh:90
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition h2oair.hh:699
Dumux::Components::Air< Scalar > Air
Definition h2oair.hh:68
static Scalar componentEnthalpy(const FluidState &fluidState, int phaseIdx, int componentIdx)
Returns the specific enthalpy  of a component in a specific phase.
Definition h2oair.hh:661
static constexpr bool isIdealMixture(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition h2oair.hh:132
static Scalar criticalPressure(int compIdx)
Critical pressure of a component .
Definition h2oair.hh:249
static constexpr int phase1Idx
Definition h2oair.hh:76
static constexpr int comp0Idx
Definition h2oair.hh:80
static constexpr int numPhases
Definition h2oair.hh:70
static constexpr bool isCompressible(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition h2oair.hh:150
static constexpr int H2OIdx
Definition h2oair.hh:78
static constexpr bool viscosityIsConstant(int phaseIdx)
Returns true if and only if a fluid phase is assumed to have a constant viscosity.
Definition h2oair.hh:166
H2Otype H2O
Definition h2oair.hh:67
static Scalar fugacityCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Returns the fugacity coefficient  of a component in a phase.
Definition h2oair.hh:532
static Scalar relativeHumidity(const FluidState &fluidState)
Returns the relative humidity of the gas phase.
Definition h2oair.hh:560
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
Calculate the molar density  of a fluid phase.
Definition h2oair.hh:421
static constexpr int liquidCompIdx
Definition h2oair.hh:82
static void init()
Initialize the fluid system's static parameters generically.
Definition h2oair.hh:329
static Scalar density(const FluidState &fluidState, const int phaseIdx)
Given a phase's composition, temperature, pressure, and the partial pressures of all components,...
Definition h2oair.hh:380
static constexpr int AirIdx
Definition h2oair.hh:79
static Scalar acentricFactor(int compIdx)
The acentric factor of a component .
Definition h2oair.hh:308
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 h2oair.hh:579
static constexpr int phase0Idx
Definition h2oair.hh:75
static constexpr int gasCompIdx
Definition h2oair.hh:83
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition h2oair.hh:104
static Scalar criticalTemperature(int compIdx)
Critical temperature of a component .
Definition h2oair.hh:233
static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase .
Definition h2oair.hh:569
static Scalar molarMass(int compIdx)
Return the molar mass of a component .
Definition h2oair.hh:218
Relations valid for an ideal gas.
Definition idealgas.hh:25
static constexpr Scalar density(Scalar avgMolarMass, Scalar temperature, Scalar pressure)
The density of the gas in , depending on pressure, temperature and average molar mass of the gas.
Definition idealgas.hh:37
static constexpr Scalar molarDensity(Scalar temperature, Scalar pressure)
The molar density of the gas , depending on pressure and temperature.
Definition idealgas.hh:58
Some exceptions thrown in DuMux
Material properties of pure water .
Binary coefficients for water and air.
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 liquidPhase() noexcept
I/O name of liquid phase.
Definition name.hh:107
Policy for the H2O-air fluid system.
Definition h2oair.hh:39
static constexpr bool useAirViscosityAsGasMixtureViscosity()
Definition h2oair.hh:42
static constexpr bool useIdealGasDensity()
Definition h2oair.hh:41
static constexpr bool useH2ODensityAsLiquidMixtureDensity()
Definition h2oair.hh:40
Tabulates all thermodynamic properties of a given untabulated chemical species.