2p1c.hh

Go to the documentation of this file.

// -*- 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 3 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/>.   *
 *****************************************************************************/
#ifndef DUMUX_2P_1C_FLUID_SYSTEM_HH
#define DUMUX_2P_1C_FLUID_SYSTEM_HH
 
#include <limits>
#include <cassert>
#include <iostream>
 
#include <dune/common/exceptions.hh>
 
#include <dumux/io/name.hh>
 
#include "base.hh"
 
namespace Dumux {
namespace FluidSystems {
template <class Scalar, class ComponentType>
class TwoPOneC
    : public Base<Scalar, TwoPOneC<Scalar, ComponentType> >
{
    using ThisType = TwoPOneC<Scalar, ComponentType>;
    using Base = Dumux::FluidSystems::Base<Scalar, ThisType>;
    using Component = ComponentType;
 
public:
    static constexpr int numPhases = 2; 
    static constexpr int numComponents = 1; 
 
    static constexpr int liquidPhaseIdx = 0; 
    static constexpr int gasPhaseIdx = 1; 
    static constexpr int phase0Idx = liquidPhaseIdx; 
    static constexpr int phase1Idx = gasPhaseIdx; 
 
    static constexpr int comp0Idx = 0; 
 
    /****************************************
    * Fluid phase related static parameters
    ****************************************/
    static std::string phaseName(int phaseIdx)
    {
        static std::string name[] = {
            std::string(IOName::liquidPhase()),
            std::string(IOName::gaseousPhase()),
        };
 
        assert(0 <= phaseIdx && phaseIdx < numPhases);
        return name[phaseIdx];
    }
 
    static constexpr bool isMiscible()
    { return false; }
 
    static constexpr bool isGas(int phaseIdx)
    {
        assert(0 <= phaseIdx && phaseIdx < numPhases);
        return phaseIdx == gasPhaseIdx;
    }
 
    static constexpr bool isIdealMixture(int phaseIdx)
    {
        assert(0 <= phaseIdx && phaseIdx < numPhases);
        // we assume Henry's and Raoult's laws for the water phase and
        // and no interaction between gas molecules of different
        // components, so all phases are ideal mixtures!
        return true;
    }
 
    static constexpr bool isCompressible(int phaseIdx)
    {
        assert(0 <= phaseIdx && phaseIdx < numPhases);
        // gases are always compressible
        if (phaseIdx == gasPhaseIdx)
            return true;
        // the component decides for the liquid phase...
        return Component::liquidIsCompressible();
    }
 
    static constexpr bool isIdealGas(int phaseIdx)
    {
        assert(0 <= phaseIdx && phaseIdx < numPhases);
 
        if (phaseIdx == gasPhaseIdx)
            // let the components decide
            return Component::gasIsIdeal();
        return false; // not a gas
    }
 
    /****************************************
     * Component related static parameters
     ****************************************/
 
    static std::string componentName(int compIdx)
    {
        assert(0 <= compIdx && compIdx < numComponents);
        return Component::name();
    }
 
    static Scalar molarMass(int compIdx)
    {
        assert(0 <= compIdx && compIdx < numComponents);
        return Component::molarMass();
    }
 
    static Scalar criticalTemperature(int compIdx)
    {
        assert(0 <= compIdx && compIdx < numComponents);
        return Component::criticalTemperature();
    }
 
    static Scalar criticalPressure(int compIdx)
    {
        assert(0 <= compIdx && compIdx < numComponents);
        return Component::criticalPressure();
    }
 
    static Scalar criticalMolarVolume(int compIdx)
    {
        DUNE_THROW(Dune::NotImplemented,
                   "TwoPOneC::criticalMolarVolume()");
        return 0;
    }
 
    static Scalar acentricFactor(int compIdx)
    {
        assert(0 <= compIdx && compIdx < numComponents);
        return Component::acentricFactor();
    }
 
    /****************************************
     * thermodynamic relations
     ****************************************/
 
    static void init()
    {
        init(/*tempMin=*/273.15,
             /*tempMax=*/623.15,
             /*numTemp=*/100,
             /*pMin=*/0.0,
             /*pMax=*/20e6,
             /*numP=*/200);
    }
 
    static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
                     Scalar pressMin, Scalar pressMax, unsigned nPress)
    {
        if (Component::isTabulated)
        {
            Component::init(tempMin, tempMax, nTemp,
                            pressMin, pressMax, nPress);
        }
    }
 
    using Base::density;
    template <class FluidState>
    static Scalar density(const FluidState &fluidState,
                          int phaseIdx)
    {
        assert(0 <= phaseIdx  && phaseIdx < numPhases);
 
        Scalar temperature = fluidState.temperature(phaseIdx);
        Scalar pressure = fluidState.pressure(phaseIdx);
 
        // liquid phase
        if (phaseIdx == liquidPhaseIdx) {
            return Component::liquidDensity(temperature, pressure);
        }
        else if (phaseIdx == gasPhaseIdx)// gas phase
        {
            return Component::gasDensity(temperature, pressure);
        }
        else
            DUNE_THROW(Dune::InvalidStateException, "Invalid phase index.");
    }
 
    using Base::molarDensity;
    template <class FluidState>
    static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
    {
        Scalar temperature = fluidState.temperature(phaseIdx);
        Scalar pressure = fluidState.temperature(phaseIdx);
        if (phaseIdx == liquidPhaseIdx)
            return Component::liquidMolarDensity(temperature, pressure);
        else if (phaseIdx == gasPhaseIdx)
            return Component::gasMolarDensity(temperature, pressure);
        else
            DUNE_THROW(Dune::InvalidStateException, "Invalid phase index.");
    }
 
    using Base::viscosity;
    template <class FluidState>
    static Scalar viscosity(const FluidState &fluidState,
                            int phaseIdx)
    {
        assert(0 <= phaseIdx  && phaseIdx < numPhases);
 
        Scalar temperature = fluidState.temperature(phaseIdx);
        Scalar pressure = fluidState.pressure(phaseIdx);
 
        // liquid phase
        if (phaseIdx == liquidPhaseIdx)
            return Component::liquidViscosity(temperature, pressure);
        else if (phaseIdx == gasPhaseIdx) // gas phase
            return Component::gasViscosity(temperature, pressure);
        else
            DUNE_THROW(Dune::InvalidStateException, "Invalid phase index.");
    }
 
    template <class FluidState>
    static Scalar vaporTemperature(const FluidState &fluidState,
                                   const unsigned int phaseIdx)
    {
        assert(0 <= phaseIdx  && phaseIdx < numPhases);
        Scalar pressure = fluidState.pressure(gasPhaseIdx);
 
        return Component::vaporTemperature(pressure);
    }
 
    using Base::fugacityCoefficient;
    template <class FluidState>
    static Scalar fugacityCoefficient(const FluidState &fluidState,
                                      int phaseIdx,
                                      int compIdx)
    {
        assert(0 <= phaseIdx  && phaseIdx < numPhases);
        assert(0 <= compIdx  && compIdx < numComponents);
 
        Scalar temperature = fluidState.temperature(phaseIdx);
        Scalar pressure = fluidState.pressure(phaseIdx);
 
        // liquid phase
        if (phaseIdx == liquidPhaseIdx)
           return Component::vaporPressure(temperature)/pressure;
 
        // for the gas phase, assume an ideal gas when it comes to
        // fugacity (-> fugacity == partial pressure)
        else if (phaseIdx == gasPhaseIdx)
            return 1.0;
 
        else
            DUNE_THROW(Dune::InvalidStateException, "Invalid phase index.");
    }
 
 
    using Base::diffusionCoefficient;
    template <class FluidState>
    static Scalar diffusionCoefficient(const FluidState &fluidState,
                                       int phaseIdx,
                                       int compIdx)
    {
        // TODO!
        DUNE_THROW(Dune::NotImplemented, "Diffusion coefficients");
    }
 
    using Base::binaryDiffusionCoefficient;
    template <class FluidState>
    static Scalar binaryDiffusionCoefficient(const FluidState &fluidState,
                                             int phaseIdx,
                                             int compIIdx,
                                             int compJIdx)
 
    { DUNE_THROW(Dune::NotImplemented, "Binary Diffusion coefficients"); }
 
    using Base::enthalpy;
    template <class FluidState>
    static Scalar enthalpy(const FluidState &fluidState,
                           int phaseIdx)
    {
        assert(0 <= phaseIdx  && phaseIdx < numPhases);
 
        // liquid phase
        if (phaseIdx == liquidPhaseIdx)
            return Component::liquidEnthalpy(fluidState.temperature(phaseIdx),
                                             fluidState.pressure(phaseIdx));
 
        else if (phaseIdx == gasPhaseIdx) // gas phase
            return Component::gasEnthalpy(fluidState.temperature(phaseIdx),
                                          fluidState.pressure(phaseIdx));
 
        else
            DUNE_THROW(Dune::InvalidStateException, "Invalid phase index.");
    }
 
    using Base::thermalConductivity;
    template <class FluidState>
    static Scalar thermalConductivity(const FluidState &fluidState,
                                      const int phaseIdx)
    {
        assert(0 <= phaseIdx  && phaseIdx < numPhases);
        // liquid phase
        if (phaseIdx == liquidPhaseIdx)
            return Component::liquidThermalConductivity(fluidState.temperature(phaseIdx),
                                                        fluidState.pressure(phaseIdx)); //0.68 ;
 
        else if (phaseIdx == gasPhaseIdx) // gas phase
            return Component::gasThermalConductivity(fluidState.temperature(phaseIdx),
                                                     fluidState.pressure(phaseIdx)); //0.0248;
 
        else
            DUNE_THROW(Dune::InvalidStateException, "Invalid phase index.");
    }
 
    using Base::heatCapacity;
    template <class FluidState>
    static Scalar heatCapacity(const FluidState &fluidState,
                               int phaseIdx)
    {
        assert(0 <= phaseIdx  && phaseIdx < numPhases);
        // liquid phase
        if (phaseIdx == liquidPhaseIdx)
            return Component::liquidHeatCapacity(fluidState.temperature(phaseIdx),
                                                 fluidState.pressure(phaseIdx));//4.217e3 ;
 
        else if (phaseIdx == gasPhaseIdx) // gas phase
            return Component::gasHeatCapacity(fluidState.temperature(phaseIdx),
                                              fluidState.pressure(phaseIdx));//2.029e3;
 
        else
            DUNE_THROW(Dune::InvalidStateException, "Invalid phase index.");
    }
};
 
} // end namespace FluidSystems
} // end namespace Dumux
 
#endif