Commit 0431ae6d authored by Shirin's avatar Shirin
Browse files

el2p_Decoupled for relative Pn=0, and corrected for avg sw for initial stress values.

parent dfac8ae9
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
*
* \brief Binary coefficients for water and nitrogen.
*/
#ifndef DUMUX_BINARY_COEFF_H2O_AIR_RELPRESS_HH
#define DUMUX_BINARY_COEFF_H2O_AIR_RELPRESS_HH
#include <cmath>
namespace Dumux
{
namespace BinaryCoeff
{
/*!
* \ingroup Binarycoefficients
* \brief Binary coefficients for water and nitrogen.
*/
class H2O_Air
{
public:
/*!
* \brief Henry coefficient \f$\mathrm{[Pa]}\f$ for air in liquid water.
* \param temperature the temperature \f$\mathrm{[K]}\f$
*
* Henry coefficient See:
* Stefan Finsterle (1993, page 33 Formula (2.9)) \cite finsterle1993 <BR>
* (fitted to data from Tchobanoglous & Schroeder, 1985 \cite tchobanoglous1985 )
*/
template <class Scalar>
static Scalar henry(Scalar temperature)
{
using std::exp;
Scalar r = (0.8942+1.47*exp(-0.04394*(temperature-273.15)))*1.E-10;
return 1./r;
}
/*!
* \brief Binary diffusion coefficient \f$\mathrm{[m^2/s]}\f$ for molecular water and air
*
* \param temperature the temperature \f$\mathrm{[K]}\f$
* \param pressure the phase pressure \f$\mathrm{[Pa]}\f$
* Vargaftik: Tables on the thermophysical properties of liquids and gases.
* John Wiley & Sons, New York, 1975. \cite vargaftik1975 <BR>
* Walker, Sabey, Hampton: Studies of heat transfer and water migration in soils.
* Dep. of Agricultural and Chemical Engineering, Colorado State University,
* Fort Collins, 1981. \cite walker1981
*/
template <class Scalar>
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
{
const Scalar Theta=1.8;
const Scalar Daw=2.13e-5; /* reference value */
const Scalar pg0=1.e5; /* reference pressure */
const Scalar T0=273.15; /* reference temperature */
Scalar Dgaw;
using std::pow;
Dgaw=Daw*(pg0/(pressure+1e5))*pow((temperature/T0),Theta);
return Dgaw;
}
/*!
* Lacking better data on water-air diffusion in liquids, we use at the
* moment the diffusion coefficient of the air's main component nitrogen!!
* \brief Diffusion coefficient \f$\mathrm{[m^2/s]}\f$ for molecular nitrogen in liquid water.
*
* \param temperature the temperature \f$\mathrm{[K]}\f$
* \param pressure the phase pressure \f$\mathrm{[Pa]}\f$
*
* The empirical equations for estimating the diffusion coefficient in
* infinite solution which are presented in Reid, 1987 all show a
* linear dependency on temperature. We thus simply scale the
* experimentally obtained diffusion coefficient of Ferrell and
* Himmelblau by the temperature.
*
* See:
* R. Reid et al. (1987, pp. 599) \cite reid1987 <BR>
* R. Ferrell, D. Himmelblau (1967, pp. 111-115) \cite ferrell1967
*/
template <class Scalar>
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
{
const Scalar Texp = 273.15 + 25; // [K]
const Scalar Dexp = 2.01e-9; // [m^2/s]
return Dexp * temperature/Texp;
}
};
}
} // end namespace
#endif
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
*
* \ingroup Components
*
* \brief A simple class for the air fluid properties
*/
#ifndef DUMUX_AIR_RELPRESS_HH
#define DUMUX_AIR_RELPRESS_HH
#include <dumux/common/exceptions.hh>
#include <dumux/material/components/component.hh>
#include <dumux/material/idealgas.hh>
namespace Dumux
{
/*!
* \ingroup Components
*
* \brief A class for the air fluid properties
*
* \tparam Scalar The type used for scalar values
*/
template <class Scalar>
class AirRelPress : public Component<Scalar, AirRelPress<Scalar> >
{
typedef Dumux::IdealGas<Scalar> IdealGas;
public:
/*!
* \brief A human readable name for Air.
*/
static std::string name()
{ return "Air"; }
/*!
* \brief The molar mass in \f$\mathrm{[kg/mol]}\f$ of Air.
*
* Taken from constrelair.hh.
*/
static Scalar molarMass()
{ return 0.02896; /* [kg/mol] */ }
/*!
* \brief Returns the critical temperature \f$\mathrm{[K]}\f$ of Air.
*/
static Scalar criticalTemperature()
{ return 132.6312; /* [K] */ }
/*!
* \brief Returns the critical pressure \f$\mathrm{[Pa]}\f$ of Air.
*/
static Scalar criticalPressure()
{ return 37.86e5-1e5; /* [Pa] */ }
/*!
* \brief The density \f$\mathrm{[kg/m^3]}\f$ of Air at a given pressure and temperature.
*
* Ideal gas is assumed.
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of phase in \f$\mathrm{[Pa]}\f$
*/
static Scalar gasDensity(Scalar temperature, Scalar pressure)
{
// Assume an ideal gas
return IdealGas::density(molarMass(), temperature, pressure+1e5);
}
/*!
* \brief Returns true, the gas phase is assumed to be compressible
*/
static bool gasIsCompressible()
{ return true; }
/*!
* \brief Returns true, the gas phase is assumed to be ideal
*/
static bool gasIsIdeal()
{ return true; }
/*!
* \brief The pressure \f$\mathrm{[Pa]}\f$ of gaseous Air at a given density and temperature.
*
* Ideal gas is assumed.
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param density density of component in \f$\mathrm{[kg/m^3]}\f$
*/
static Scalar gasPressure(Scalar temperature, Scalar density)
{
// Assume an ideal gas
return IdealGas::pressure(temperature, density/molarMass()) - 1e5;
}
/*!
* \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of Air at a given pressure and temperature.
*
* Criticial specific volume calculated by \f$V_c = (R*T_c)/p_c\f$.
*
* Reid et al. (1987, pp 396-397, 667) \cite reid1987 <BR>
* Poling et al. (2001, pp 9.7-9.8) \cite poling2001 <BR>
*
* Accentric factor taken from: <BR>
* Adebiyi (2003) \cite adebiyi2003
*
* air is a non-polar substance,
* thus dipole moment mu is zero, as well the dimensionless dipole moment mu_r
* therefore not considered below
* the same holds for the correction value kappa for highly polar substances
*
* This calculation was introduced into Dumux in 2012 although the method here
* is designed for general polar substances. Air, however, is (a) non-polar,
* and (b) there are more precise methods available
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
static Scalar oldGasViscosity(Scalar temperature, Scalar pressure)
{
const Scalar Tc = criticalTemperature();
const Scalar Vc = 84.525138; // critical specific volume [cm^3/mol]
const Scalar omega = 0.078; // accentric factor
const Scalar M = molarMass() * 1e3; // molar mas [g/mol]
const Scalar Fc = 1.0 - 0.2756*omega;
const Scalar Tstar = 1.2593*temperature/Tc;
using std::exp;
using std::pow;
const Scalar Omega_v = 1.16145*pow(Tstar, -0.14874)
+ 0.52487*exp(-0.77320*Tstar)
+ 2.16178*exp(-2.43787*Tstar);
using std::cbrt;
using std::sqrt;
const Scalar mu = 40.785 * Fc * sqrt(M * temperature)/(cbrt(Vc * Vc) * Omega_v);
// convertion from micro poise to Pa s
return mu/1.0e6/10.0;
}
/*!
* \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of Air at a given pressure and temperature.
*
* Simple method, already implemented in MUFTE-UG, but pretty accurate.
*
* The pressure correction is even simpler and developed and tested by
* Holger Class in 2016 against the results of the Lemmon and Jacobsen (2004)
* approach \cite LemmonJacobsen2004
* It shows very reasonable results throughout realistic pressure and
* temperature ranges up to several hundred Kelvin and up to 500 bar
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
static Scalar gasViscosity(Scalar temperature, Scalar pressure)
{
// above 1200 K, the function becomes inaccurate
// since this should realistically never happen, we can live with it
const Scalar tempCelsius = temperature - 273.15;
const Scalar pressureCorrectionFactor = 9.7115e-9*tempCelsius*tempCelsius - 5.5e-6*tempCelsius + 0.0010809;
using std::sqrt;
const Scalar mu = 1.496e-6 * sqrt(temperature * temperature * temperature) / (temperature + 120.0)
* (1.0 + ((pressure+1e5)/1.0e5 - 1.0)*pressureCorrectionFactor);
return mu;
}
/*!
* \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of Air at a given pressure and temperature.
*
* Simple method, already implemented in MUFTE-UG, but pretty accurate
* at atmospheric pressures.
* Gas viscosity is not very dependent on pressure. Thus, for
* low pressures one might switch the pressure correction off
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
static Scalar simpleGasViscosity(Scalar temperature, Scalar pressure)
{
// above 1200 K, the function becomes inaccurate
// since this should realistically never happen, we can live with it
using std::sqrt;
return 1.496e-6 * sqrt(temperature * temperature * temperature) / (temperature + 120.0);
}
/*!
* \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of Air at a given pressure and temperature.
*
* This is a very exact approach by Lemmon and Jacobsen (2004) \cite LemmonJacobsen2004
* All the values and parameters used below are explained in their paper
* Since they use ''eta'' for dyn. viscosity, we do it as well for easier
* comparison with the paper
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
static Scalar exactGasViscosity(Scalar temperature, Scalar pressure)
{
const Scalar epsk = 103.3; // [K]
using namespace std;
const Scalar logTstar = log(temperature/epsk);
const Scalar Omega = exp(0.431
- 0.4623*logTstar
+ 0.08406*logTstar*logTstar
+ 0.005341*logTstar*logTstar*logTstar
- 0.00331*logTstar*logTstar*logTstar*logTstar);
const Scalar sigma = 0.36; // [nm]
const Scalar eta0 = 0.0266958*sqrt(1000.0*molarMass()*temperature)/(sigma*sigma*Omega);
const Scalar tau = criticalTemperature()/temperature;
const Scalar rhoc = 10.4477; // [mol/m^3]
const Scalar delta = 0.001*pressure/(temperature*8.3144598)/rhoc;
const Scalar etaR = 10.72 * pow(tau, 0.2) * delta
+ 1.122 * pow(tau, 0.05) * pow(delta, 4)
+ 0.002019 * pow(tau, 2.4) * pow(delta, 9)
- 8.876 * pow(tau, 0.6) * delta * exp(-delta)
- 0.02916 * pow(tau, 3.6) * pow(delta, 8) * exp(-delta);
return (eta0 + etaR)*1e-6;
}
/*!
* \brief Specific enthalpy of Air \f$\mathrm{[J/kg]}\f$
* with 273.15 \f$ K \f$ as basis. <BR>
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*
* Kays et al. (2005, 431ff) \cite kays2005 <BR>
*/
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
{
return 1005*(temperature-273.15);
}
/*!
* \brief Specific internal energy of Air \f$\mathrm{[J/kg]}\f$.
*
* Definition of enthalpy: \f$h= u + pv = u + p / \rho\f$.
* Rearranging for internal energy yields: \f$u = h - pv\f$.
* Exploiting the Ideal Gas assumption
* (\f$pv = R_{\textnormal{specific}} T\f$) gives: \f$u = h - R / M T \f$.
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
static const Scalar gasInternalEnergy(Scalar temperature,
Scalar pressure)
{
return gasEnthalpy(temperature, pressure+1e5)
- IdealGas::R * temperature // = pressure * molar volume for an ideal gas
/ molarMass(); // conversion from [J/(mol K)] to [J/(kg K)]
}
/*!
* \brief Specific isobaric heat capacity \f$\mathrm{[J/(kg*K)]}\f$ of pure
* air.
*
* This methods uses the formula for "zero-pressure" heat capacity that
* is only dependent on temperature, because the pressure dependence is rather small.
* This one should be accurate for a pressure of 1 atm.
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*
* Values taken from Hollis (1996) \cite hollis1996 <BR>
* "Tables of Thermal Properties of Gases"
*/
static const Scalar gasHeatCapacity(Scalar temperature,
Scalar pressure)
{
// scale temperature with reference temp of 100K
Scalar phi = temperature/100;
using std::pow;
Scalar c_p = 0.661738E+01
-0.105885E+01 * phi
+0.201650E+00 * pow(phi,2)
-0.196930E-01 * pow(phi,3)
+0.106460E-02 * pow(phi,4)
-0.303284E-04 * pow(phi,5)
+0.355861E-06 * pow(phi,6);
c_p += -0.549169E+01 * pow(phi,-1)
+0.585171E+01 * pow(phi,-2)
-0.372865E+01 * pow(phi,-3)
+0.133981E+01 * pow(phi,-4)
-0.233758E+00 * pow(phi,-5)
+0.125718E-01 * pow(phi,-6);
c_p *= IdealGas::R / (molarMass() * 1000); // in J/mol/K * mol / kg / 1000 = kJ/kg/K
return c_p;
}
/*!
* \brief Thermal conductivity \f$\mathrm{[[W/(m*K)]}\f$ of air.
*
* Isobaric Properties for Nitrogen in: NIST Standard \cite NIST <BR>
* evaluated at p=.1 MPa, T=20°C <BR>
* Nitrogen: 0.025398 <BR>
* Oxygen: 0.026105 <BR>
* lambda_air is approximately 0.78*lambda_N2+0.22*lambda_O2
*
* \param temperature absolute temperature in \f$\mathrm{[K]}\f$
* \param pressure of the phase in \f$\mathrm{[Pa]}\f$
*/
static Scalar gasThermalConductivity(Scalar temperature, Scalar pressure)
{
return 0.0255535;
}
};
} // end namespace
#endif
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup Components
* \brief A much simpler (and thus potentially less buggy) version of
* pure water.
*/
#ifndef DUMUX_SIMPLE_H2O_RELPRESS_HH
#define DUMUX_SIMPLE_H2O_RELPRESS_HH
#include <dumux/material/idealgas.hh>
#include "component.hh"
#include <cmath>
namespace Dumux
{
/*!
* \ingroup Components
*
* \brief A much simpler (and thus potentially less buggy) version of
* pure water.
*
* \tparam Scalar The type used for scalar values
*/
template <class Scalar>
class SimpleH2ORelPress : public Component<Scalar, SimpleH2ORelPress<Scalar> >
{
typedef Dumux::IdealGas<Scalar> IdealGas;
static const Scalar R; // specific gas constant of water
public:
/*!
* \brief A human readable name for the water.
*/
static std::string name()
{ return "H2O"; }
/*!
* \brief The molar mass in \f$\mathrm{[kg/mol]}\f$ of water.
*/
static Scalar molarMass()
{ return 18e-3; }
/*!
* \brief Returns the critical temperature \f$\mathrm{[K]}\f$ of water.
*/
static Scalar criticalTemperature()
{ return 647.096; /* [K] */ }
/*!
* \brief Returns the critical pressure \f$\mathrm{[Pa]}\f$ of water.
*/
static Scalar criticalPressure()
{ return 22.064e6-1e5; /* [N/m^2] */ }
/*!
* \brief Returns the temperature \f$\mathrm{[K]}\f$ at water's triple point.
*/
static Scalar tripleTemperature()
{ return 273.16; /* [K] */ }
/*!
* \brief Returns the pressure \f$\mathrm{[Pa]}\f$ at water's triple point.
*/
static Scalar triplePressure()
{ return 611.657-1e5; /* [N/m^2] */ }
/*!
* \brief The vapor pressure in \f$\mathrm{[Pa]}\f$ of pure water
* at a given temperature.
*
*\param T temperature of component in \f$\mathrm{[K]}\f$
*
* See:
*
* IAPWS: "Revised Release on the IAPWS Industrial Formulation
* 1997 for the Thermodynamic Properties of Water and Steam",
* http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997
*/
static Scalar vaporPressure(Scalar T)
{
if (T > criticalTemperature())
return criticalPressure();
if (T < tripleTemperature())
return 0; // water is solid: We don't take sublimation into account
static const Scalar n[10] = {
0.11670521452767e4, -0.72421316703206e6, -0.17073846940092e2,
0.12020824702470e5, -0.32325550322333e7, 0.14915108613530e2,
-0.48232657361591e4, 0.40511340542057e6, -0.23855557567849,
0.65017534844798e3
};
Scalar sigma = T + n[8]/(T - n[9]);
Scalar A = (sigma + n[0])*sigma + n[1];
Scalar B = (n[2]*sigma + n[3])*sigma + n[4];
Scalar C = (n[5]*sigma + n[6])*sigma + n[7];
using std::sqrt;
Scalar tmp = Scalar(2.0)*C/(sqrt(B*B - Scalar(4.0)*A*C) - B);
tmp *= tmp;
tmp *= tmp;
return Scalar(1e6)*tmp - 1e5;
}
/*!
* \brief Specific enthalpy of water steam \f$\mathrm{[J/kg]}\f$.
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
static const Scalar gasEnthalpy(Scalar temperature,
Scalar pressure)
{ return 1976*(temperature - 293.15) + 2.45e6; }
/*!
* \brief Specific enthalpy of liquid water \f$\mathrm{[J/kg]}\f$.
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
static const Scalar liquidEnthalpy(Scalar temperature,
Scalar pressure)
{
return 4180*(temperature - 293.15);
}