12#ifndef DUMUX_BRINE_FLUID_SYSTEM_HH 
   13#define DUMUX_BRINE_FLUID_SYSTEM_HH 
   15#include <dune/common/math.hh> 
   32template< 
class Scalar, 
class H2OType = Components::TabulatedComponent<Dumux::Components::H2O<Scalar>> >
 
   33class Brine : 
public Base< Scalar, Brine<Scalar, H2OType>>
 
  114        return H2O::liquidIsCompressible();
 
 
  142            case H2OIdx: 
return H2O::name();
 
  145        DUNE_THROW(Dune::InvalidStateException, 
"Invalid component index " << compIdx);
 
 
  156            case H2OIdx: 
return H2O::molarMass();
 
  159        DUNE_THROW(Dune::InvalidStateException, 
"Invalid component index " << compIdx);
 
 
  195        if (H2O::isTabulated)
 
  197            std::cout << 
"Initializing tables for the H2O fluid properties (" 
  201            H2O::init(tempMin, tempMax, nTemp, pressMin, pressMax, nPress);
 
 
  216    template <
class Flu
idState>
 
  220        const Scalar temperature = fluidState.temperature(phaseIdx);
 
  221        const Scalar pressure = fluidState.pressure(phaseIdx);
 
  222        const Scalar xNaCl = fluidState.massFraction(phaseIdx, 
NaClIdx);
 
  225        const Scalar TempC = temperature - 273.15;
 
  226        const Scalar pMPa = pressure/1.0E6;
 
  227        const Scalar salinity = max(0.0, xNaCl);
 
  229        const Scalar rhow = H2O::liquidDensity(temperature, pressure);
 
 
  245    template <
class Flu
idState>
 
  250        assert(0 <= phaseIdx  && phaseIdx < 
numPhases);
 
  253        if (phaseIdx == compIdx)
 
  259        return std::numeric_limits<Scalar>::infinity();
 
 
  274    template <
class Flu
idState>
 
  278        const Scalar temperature = fluidState.temperature(phaseIdx);
 
  279        const Scalar xNaCl = fluidState.massFraction(phaseIdx, 
NaClIdx);
 
  285        const Scalar T = max(temperature, 275.0);
 
  286        const Scalar salinity = max(0.0, xNaCl);
 
  288        const Scalar T_C = T - 273.15;
 
  289        const Scalar A = ((0.42*power((pow(salinity, 0.8)-0.17), 2)) + 0.045)*pow(T_C, 0.8);
 
  290        const Scalar mu_brine = 0.1 + (0.333*salinity) + (1.65+(91.9*salinity*salinity*salinity))*exp(-A); 
 
  291        assert(mu_brine > 0.0);
 
  292        return mu_brine/1000.0; 
 
 
  305    template <
class Flu
idState>
 
  310            const Scalar temperature = fluidState.temperature(
H2OIdx);
 
  312            return H2O::vaporPressure(temperature)*fluidState.moleFraction(
phase0Idx, 
H2OIdx);
 
  315            DUNE_THROW(Dune::NotImplemented, 
"NaCl::vaporPressure(t)");
 
  317            DUNE_THROW(Dune::NotImplemented, 
"Invalid component index " << compIdx);
 
 
  334    template <
class Flu
idState>
 
  338        return enthalpy_(fluidState.pressure(phaseIdx),
 
  339                        fluidState.temperature(phaseIdx),
 
 
  350    template <
class Flu
idState>
 
  360            if (componentIdx == 
H2OIdx)
 
  361                return H2O::liquidEnthalpy(T, p);
 
  362            else if (componentIdx == 
NaClIdx)
 
  363                DUNE_THROW(Dune::NotImplemented, 
"The component enthalpy for NaCl is not implemented.");
 
  364            DUNE_THROW(Dune::InvalidStateException, 
"Invalid component index " << componentIdx);
 
  366        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  371    template <
class Flu
idState>
 
  374        return density(fluidState, phaseIdx)/fluidState.averageMolarMass(phaseIdx);
 
 
  379    template <
class Flu
idState>
 
  382        DUNE_THROW(Dune::NotImplemented, 
"FluidSystems::Brine::diffusionCoefficient()");
 
 
  400    template <
class Flu
idState>
 
  408            if (compIIdx > compJIdx)
 
  411                swap(compIIdx, compJIdx);
 
  417                DUNE_THROW(Dune::NotImplemented, 
"Binary diffusion coefficient of components " 
  418                                                 << compIIdx << 
" and " << compJIdx
 
  419                                                 << 
" in phase " << phaseIdx);
 
  422         DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index: " << phaseIdx);
 
 
  433    template <
class Flu
idState>
 
  438            Scalar tempC = fluidState.temperature(phaseIdx)-273.15;
 
  440            Scalar S = 5844.3 * m / (1000 + 58.443 *m);
 
  441            Scalar contribNaClFactor = 1.0 - (2.3434e-3 - 7.924e-6*tempC + 3.924e-8*tempC*tempC)*S + (1.06e-5 - 2.0e-8*tempC + 1.2e-10*tempC*tempC)*S*S;
 
  442            return contribNaClFactor * H2O::liquidThermalConductivity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
 
  444        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index: " << phaseIdx);
 
 
  449    template <
class Flu
idState>
 
  453            const Scalar eps = fluidState.temperature(phaseIdx)*1e-8;
 
  455            return (enthalpy_(fluidState.pressure(phaseIdx),
 
  456                        fluidState.temperature(phaseIdx) +eps ,
 
  458                        - enthalpy_(fluidState.pressure(phaseIdx),
 
  459                        fluidState.temperature(phaseIdx),
 
  462        DUNE_THROW(Dune::InvalidStateException, 
"Invalid phase index " << phaseIdx);
 
 
  476    static Scalar enthalpy_(
const Scalar p, 
const Scalar T, 
const Scalar xNaCl)
 
  479        static const Scalar f[] = { 2.63500E-1, 7.48368E-6, 1.44611E-6, -3.80860E-10 };
 
  482        static const Scalar a[4][3] = { { +9633.6, -4080.0, +286.49 },
 
  483                                        { +166.58, +68.577, -4.6856 },
 
  484                                        { -0.90963, -0.36524, +0.249667E-1 },
 
  485                                        { +0.17965E-2, +0.71924E-3, -0.4900E-4 } };
 
  487        const Scalar theta = T - 273.15;
 
  488        const Scalar salSat = f[0] + f[1]*theta + f[2]*theta*theta + f[3]*theta*theta*theta;
 
  493        const Scalar salinity = min(max(xNaCl,0.0), salSat);
 
  495        const Scalar hw = H2O::liquidEnthalpy(T, p)/1E3; 
 
  498        const Scalar h_NaCl = (3.6710E4*T + 0.5*(6.2770E1)*T*T - ((6.6670E-2)/3)*T*T*T
 
  499                              + ((2.8000E-5)/4)*(T*T*T*T))/(58.44E3)- 2.045698e+02; 
 
  501        const Scalar m = (1E3/58.44)*(salinity/(1-salinity));
 
  505        for (
int i = 0; i<=3; i++)
 
  506            for (
int j=0; j<=2; j++)
 
  507                d_h = d_h + a[i][j] * power(theta, i) * power(m, j);
 
  510        const Scalar delta_h = (4.184/(1E3 + (58.44 * m)))*d_h;
 
  513        const Scalar h_ls1 =(1-salinity)*hw + salinity*h_NaCl + salinity*delta_h; 
 
 
A class for the NaCl properties.
Definition nacl.hh:34
static std::string name()
A human readable name for the NaCl.
Definition nacl.hh:39
static constexpr Scalar molarMass()
The molar mass of NaCl in .
Definition nacl.hh:47
Fluid system base class.
Definition fluidsystems/base.hh:32
ScalarType Scalar
export the scalar type
Definition fluidsystems/base.hh:35
A compositional single phase fluid system consisting of two components, which are H2O and NaCl.
Definition fluidsystems/brine.hh:34
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition fluidsystems/brine.hh:434
static Scalar viscosity(const FluidState &fluidState, int phaseIdx=liquidPhaseIdx)
Return the viscosity of the phase.
Definition fluidsystems/brine.hh:275
static Scalar componentEnthalpy(const FluidState &fluidState, int phaseIdx, int componentIdx)
Returns the specific enthalpy  of a component in a specific phase.
Definition fluidsystems/brine.hh:351
static bool isIdealGas(int phaseIdx=liquidPhaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition fluidsystems/brine.hh:124
static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase .
Definition fluidsystems/brine.hh:380
static constexpr int phase0Idx
Definition fluidsystems/brine.hh:45
static constexpr bool isGas(int phaseIdx=liquidPhaseIdx)
Return whether a phase is gaseous.
Definition fluidsystems/brine.hh:76
static constexpr int comp0Idx
Definition fluidsystems/brine.hh:50
static const int numComponents
Definition fluidsystems/brine.hh:43
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 fluidsystems/brine.hh:191
static Scalar fugacityCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Calculate the fugacity coefficient  of an individual component in a fluid phase.
Definition fluidsystems/brine.hh:246
static constexpr int comp1Idx
Definition fluidsystems/brine.hh:51
static constexpr int H2OIdx
Definition fluidsystems/brine.hh:48
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity  of a fluid phase .
Definition fluidsystems/brine.hh:450
static Scalar molarMass(int compIdx)
Return the molar mass of a component in .
Definition fluidsystems/brine.hh:152
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given a phase's composition, temperature and pressure, return its specific enthalpy .
Definition fluidsystems/brine.hh:335
static Scalar vaporPressure(const FluidState &fluidState, int compIdx)
Vapor pressure of a component .
Definition fluidsystems/brine.hh:306
static const std::string phaseName(int phaseIdx=liquidPhaseIdx)
Return the human readable name of the phase.
Definition fluidsystems/brine.hh:57
Dumux::Components::NaCl< Scalar > NaCl
Definition fluidsystems/brine.hh:40
static constexpr int NaClIdx
Definition fluidsystems/brine.hh:49
H2OType H2O
Definition fluidsystems/brine.hh:39
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition fluidsystems/brine.hh:67
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx=liquidPhaseIdx)
Calculate the molar density  of a fluid phase.
Definition fluidsystems/brine.hh:372
static constexpr int liquidPhaseIdx
Definition fluidsystems/brine.hh:46
static bool isCompressible(int phaseIdx=liquidPhaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition fluidsystems/brine.hh:111
static std::string componentName(int compIdx)
Return the human readable name of a component.
Definition fluidsystems/brine.hh:138
static Scalar density(const FluidState &fluidState, int phaseIdx=liquidPhaseIdx)
Return the phase density [kg/m^3].
Definition fluidsystems/brine.hh:217
static bool isIdealMixture(int phaseIdx=liquidPhaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition fluidsystems/brine.hh:96
static const int numPhases
Definition fluidsystems/brine.hh:42
static void init()
Initialize the fluid system's static parameters generically.
Definition fluidsystems/brine.hh:170
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 fluidsystems/brine.hh:401
A central place for various physical constants occurring in some equations.
Some exceptions thrown in DuMux
Material properties of pure water .
Material properties of pure salt .
A collection of input/output field names for common physical quantities.
std::string liquidPhase() noexcept
I/O name of liquid phase.
Definition name.hh:107
Tabulates all thermodynamic properties of a given untabulated chemical species.