diff --git a/opm/material/fluidmatrixinteractions/BrooksCorey.hpp b/opm/material/fluidmatrixinteractions/BrooksCorey.hpp index 896cc2419..38dca5477 100644 --- a/opm/material/fluidmatrixinteractions/BrooksCorey.hpp +++ b/opm/material/fluidmatrixinteractions/BrooksCorey.hpp @@ -85,20 +85,7 @@ public: "this material law!"); /*! - * \brief The capillary pressure-saturation curves according to van Genuchten. - * - * Van Genuchten's empirical capillary pressure <-> saturation - * function is given by - * \f[ - * p_{c,wn} = p_n - p_w = ({S_w}^{-1/m} - 1)^{1/n}/\alpha - * \f] - * - * \param values A random access container which stores the - * relative pressure of each fluid phase. - * \param params The parameter object expressing the coefficients - * required by the van Genuchten law. - * \param fs The fluid state for which the capillary pressure - * ought to be calculated + * \brief The capillary pressure-saturation curves. */ template static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs) @@ -119,12 +106,12 @@ public: } /*! - * \brief The relative permeability-saturation curves according to van Genuchten. + * \brief The relative permeability-saturation curves. * * \param values A random access container which stores the * relative permeability of each fluid phase. * \param params The parameter object expressing the coefficients - * required by the van Genuchten law. + * required by the material law. * \param fs The fluid state for which the relative permeabilities * ought to be calculated */ @@ -135,6 +122,132 @@ public: values[Traits::nPhaseIdx] = krn(params, fs); } + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = 0; + if (satPhaseIdx == Traits::wPhaseIdx) + values[Traits::nPhaseIdx] = dpcwn_dSw(params, state); + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + if (satPhaseIdx == Traits::wPhaseIdx) { + values[Traits::wPhaseIdx] = twoPhaseSatDKrw_dSw(params, state.saturation(Traits::wPhaseIdx)); + values[Traits::nPhaseIdx] = 0; + } + else { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = - twoPhaseSatDKrn_dSw(params, 1 - state.saturation(Traits::nPhaseIdx)); + } + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + /*! * \brief The capillary pressure-saturation curve according to * Brooks and Corey. @@ -279,7 +392,7 @@ public: * \param params The parameters of the capillary pressure curve * (for Brooks-Corey: Entry pressure and shape factor) */ - static Scalar dkrw_dSw(const Params ¶ms, Scalar Sw) + static Scalar twoPhaseSatDKrw_dSw(const Params ¶ms, Scalar Sw) { assert(0 <= Sw && Sw <= 1); @@ -297,7 +410,7 @@ public: */ template static Scalar krn(const Params ¶ms, const FluidState &fs) - { return twoPhaseSatKrn(params, fs.saturation(Traits::wPhaseIdx)); } + { return twoPhaseSatKrn(params, 1.0 - fs.saturation(Traits::nPhaseIdx)); } static Scalar twoPhaseSatKrn(const Params ¶ms, Scalar Sw) { @@ -318,7 +431,7 @@ public: * \param params The parameters of the capillary pressure curve * (for Brooks-Corey: Entry pressure and shape factor) */ - static Scalar dkrn_dSw(const Params ¶ms, Scalar Sw) + static Scalar twoPhaseSatDKrn_dSw(const Params ¶ms, Scalar Sw) { assert(0 <= Sw && Sw <= 1); diff --git a/opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp b/opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp index e70fea6eb..1d5566e11 100644 --- a/opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp +++ b/opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp @@ -152,6 +152,128 @@ public: EffLaw::relativePermeabilities(values, params, overlayFs); } + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + typedef Opm::SaturationOverlayFluidState OverlayFluidState; + + OverlayFluidState overlayFs(state); + for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) { + overlayFs.setSaturation(phaseIdx, + effectiveSaturation(params, + state.saturation(phaseIdx), + phaseIdx)); + } + + EffLaw::dCapillaryPressures_dSaturation(values, params, overlayFs, satPhaseIdx); + + // multiply with dS_eff / dS_abs + Scalar dSeff_dSabs = dSeff_dSabs_(params, satPhaseIdx); + for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) + values[phaseIdx] *= dSeff_dSabs; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { EffLaw::dCapillaryPressures_dPressure(values, params, state, pPhaseIdx); } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { EffLaw::dCapillaryPressures_dTemperature(values, params, state); } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { EffLaw::dCapillaryPressures_dMoleFraction(values, params, state, phaseIdx, compIdx); } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + typedef Opm::SaturationOverlayFluidState OverlayFluidState; + + OverlayFluidState overlayFs(state); + for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) { + overlayFs.setSaturation(phaseIdx, + effectiveSaturation(params, + state.saturation(phaseIdx), + phaseIdx)); + } + + EffLaw::dRelativePermeabilities_dSaturation(values, params, overlayFs, satPhaseIdx); + + // multiply with dS_eff / dS_abs + Scalar dSeff_dSabs = dSeff_dSabs_(params, satPhaseIdx); + for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) + values[phaseIdx] *= dSeff_dSabs; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { EffLaw::dRelativePermeabilities_dPressure(values, params, state, pPhaseIdx); } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { EffLaw::dRelativePermeabilities_dTemperature(values, params, state); } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { EffLaw::dRelativePermeabilities_dMoleFraction(values, params, state, phaseIdx, compIdx); } + /*! * \brief The capillary pressure-saturation curve. * @@ -390,20 +512,6 @@ public: { return S*(1 - params.sumResidualSaturations()) + params.residualSaturation(phaseIdx); } private: - /*! - * \brief Convert an effective wetting saturation to an absolute one. - * - * \param Swe Effective saturation of the non-wetting phase \f$\overline{S}_n\f$. - * \param params A container object that is populated with the appropriate coefficients for the respective law. - * Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container - * is constructed accordingly. Afterwards the values are set there, too. - * \return Absolute saturation of the non-wetting phase. - */ - static Scalar SweToSw_(const Params ¶ms, Scalar Swe) - { - return Swe*(1 - params.Swr() - params.Snr()) + params.Swr(); - } - /*! * \brief Derivative of the effective saturation w.r.t. the absolute saturation. * @@ -412,8 +520,8 @@ private: * is constructed accordingly. Afterwards the values are set there, too. * \return Derivative of the effective saturation w.r.t. the absolute saturation. */ - static Scalar dSwe_dSw_(const Params ¶ms) - { return 1.0/(1 - params.Swr() - params.Snr()); } + static Scalar dSeff_dSabs_(const Params ¶ms, int phaseIdx) + { return 1.0/(1 - params.sumResidualSaturations()); } /*! * \brief Derivative of the absolute saturation w.r.t. the effective saturation. @@ -423,8 +531,8 @@ private: * is constructed accordingly. Afterwards the values are set there, too. * \return Derivative of the absolute saturation w.r.t. the effective saturation. */ - static Scalar dSw_dSwe_(const Params ¶ms) - { return 1 - params.Swr() - params.Snr(); } + static Scalar dSabs_dSeff_(const Params ¶ms, int phaseIdx) + { return 1 - params.sumResidualSaturations(); } }; } // namespace Opm diff --git a/opm/material/fluidmatrixinteractions/LinearMaterial.hpp b/opm/material/fluidmatrixinteractions/LinearMaterial.hpp index bab5e9b13..034216304 100644 --- a/opm/material/fluidmatrixinteractions/LinearMaterial.hpp +++ b/opm/material/fluidmatrixinteractions/LinearMaterial.hpp @@ -115,7 +115,7 @@ public: const Params ¶ms, const FluidState &state) { - OPM_THROW(std::runtime_error, "Not implemented: MpLinearMaterial::saturations()"); + OPM_THROW(std::runtime_error, "Not implemented: LinearMaterial::saturations()"); } /*! @@ -134,6 +134,134 @@ public: } } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + + values[satPhaseIdx] = params.pcMaxSat(satPhaseIdx) - params.pcMinSat(satPhaseIdx); + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + + // -> linear relation between 0 and 1, else constant + if (state.saturation(satPhaseIdx) >= 0 && + state.saturation(satPhaseIdx) <= 1) + { + values[satPhaseIdx] = 1.0; + } + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + /*! * \brief The difference between the pressures of the non-wetting and wetting phase. */ diff --git a/opm/material/fluidmatrixinteractions/MaterialTraits.hpp b/opm/material/fluidmatrixinteractions/MaterialTraits.hpp index 93adadcd0..0c7bb3378 100644 --- a/opm/material/fluidmatrixinteractions/MaterialTraits.hpp +++ b/opm/material/fluidmatrixinteractions/MaterialTraits.hpp @@ -28,6 +28,24 @@ #define OPM_MATERIAL_TRAITS_HH namespace Opm { +/*! + * \ingroup material + * + * \brief A generic traits class which does not provide any indices. + * + * This traits class is intended to be used by the NullMaterial + */ +template +class NullMaterialTraits +{ +public: + //! The type used for scalar floating point values + typedef ScalarT Scalar; + + //! The number of fluid phases + static const int numPhases = numPhasesV; +}; + /*! * \ingroup material * diff --git a/opm/material/fluidmatrixinteractions/NullMaterial.hpp b/opm/material/fluidmatrixinteractions/NullMaterial.hpp new file mode 100644 index 000000000..8da660a59 --- /dev/null +++ b/opm/material/fluidmatrixinteractions/NullMaterial.hpp @@ -0,0 +1,341 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/***************************************************************************** + * Copyright (C) 2011-2012 by Andreas Lauser * + * * + * 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 . * + *****************************************************************************/ +/*! + * \file + * \copydoc Opm::NullMaterial + */ +#ifndef OPM_NULL_MATERIAL_HH +#define OPM_NULL_MATERIAL_HH + +#include "NullMaterialParams.hpp" + +#include +#include + +#include + +namespace Opm +{ +/*! + * \ingroup material + * + * \brief Implements a dummy linear saturation-capillary pressure + * relation which just disables capillary pressure. + */ +template +class NullMaterial : public TraitsT +{ +public: + typedef TraitsT Traits; + typedef NullMaterialParams Params; + typedef typename Traits::Scalar Scalar; + + //! The number of fluid phases + static const int numPhases = Traits::numPhases; + + //! Specify whether this material law implements the two-phase + //! convenience API + static const bool implementsTwoPhaseApi = (numPhases == 2); + + //! Specify whether this material law implements the two-phase + //! convenience API which only depends on the phase saturations + static const bool implementsTwoPhaseSatApi = (numPhases == 2); + + //! Specify whether the quantities defined by this material law + //! are saturation dependent + //! + //! In this law, the relative permeabilities are saturation + //! dependent, even if capillary pressure is always zero + static const bool isSaturationDependent = true; + + //! Specify whether the quantities defined by this material law + //! are dependent on the absolute pressure + static const bool isPressureDependent = false; + + //! Specify whether the quantities defined by this material law + //! are temperature dependent + static const bool isTemperatureDependent = false; + + //! Specify whether the quantities defined by this material law + //! are dependent on the phase composition + static const bool isCompositionDependent = false; + + /*! + * \brief Returns constant 0 for all phases. + * + * \param values Container for the return values + * \param params Parameters + * \param state The fluid state + */ + template + static void capillaryPressures(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) + values[phaseIdx] = 0.0; + } + + /*! + * \brief The inverse of the capillary pressure + */ + template + static void saturations(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { OPM_THROW(std::logic_error, "Not defined: NullMaterial::saturations()"); } + + /*! + * \brief The relative permeability of all phases. + */ + template + static void relativePermeabilities(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) { + values[phaseIdx] = std::max(std::min(state.saturation(phaseIdx),1.0),0.0); + } + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + // -> not saturation dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + + // -> linear relation between 0 and 1, else constant + if (state.saturation(satPhaseIdx) >= 0 && + state.saturation(satPhaseIdx) <= 1) + { + values[satPhaseIdx] = 1.0; + } + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The difference between the pressures of the non-wetting and wetting phase. + */ + template + static typename std::enable_if<(numPhases > 1), ScalarT>::type + pcwn(const Params ¶ms, const FluidState &fs) + { return 0; } + + template + static typename std::enable_if::type + twoPhaseSatPcwn(const Params ¶ms, Scalar Sw) + { return 0; } + + /*! + * \brief Calculate wetting phase saturation given that the rest + * of the fluid state has been initialized + */ + template + static Scalar Sw(const Params ¶ms, const FluidState &fs) + { OPM_THROW(std::logic_error, "Not defined: Sw()"); } + + template + static typename std::enable_if::type + twoPhaseSatSw(const Params ¶ms, Scalar pcwn) + { OPM_THROW(std::logic_error, "Not defined: twoPhaseSatSw()"); } + + /*! + * \brief Calculate non-wetting phase saturation given that the + * rest of the fluid state has been initialized + */ + template + static Scalar Sn(const Params ¶ms, const FluidState &fs) + { OPM_THROW(std::logic_error, "Not defined: Sn()"); } + + template + static typename std::enable_if::type + twoPhaseSatSn(const Params ¶ms, Scalar pcwn) + { OPM_THROW(std::logic_error, "Not defined: twoPhaseSatSn()"); } + + /*! + * \brief Calculate gas phase saturation given that the rest of + * the fluid state has been initialized + * + * This method is only available for at least three fluid phases + */ + template + static typename std::enable_if< (numPhases > 2), ScalarT>::type + Sg(const Params ¶ms, const FluidState &fs) + { OPM_THROW(std::logic_error, "Not defined: Sg()"); } + + /*! + * \brief The relative permability of the wetting phase + */ + template + static typename std::enable_if<(numPhases > 1), ScalarT>::type + krw(const Params ¶ms, const FluidState &fs) + { return std::max(0.0, std::min(1.0, fs.saturation(Traits::wPhaseIdx))); } + + template + static typename std::enable_if::type + twoPhaseSatKrw(const Params ¶ms, Scalar Sw) + { return std::max(0.0, std::min(1.0, Sw)); } + + /*! + * \brief The relative permability of the liquid non-wetting phase + */ + template + static typename std::enable_if<(numPhases > 1), ScalarT>::type + krn(const Params ¶ms, const FluidState &fs) + { return std::max(0.0, std::min(1.0, fs.saturation(Traits::nPhaseIdx))); } + + template + static typename std::enable_if::type + twoPhaseSatKrn(const Params ¶ms, Scalar Sw) + { return std::max(0.0, std::min(1.0, 1 - Sw)); } + + /*! + * \brief The relative permability of the gas phase + * + * This method is only available for at least three fluid phases + */ + template + static typename std::enable_if< (numPhases > 2), ScalarT>::type + krg(const Params ¶ms, const FluidState &fs) + { return std::max(0.0, std::min(1.0, fs.saturation(Traits::gPhaseIdx))); } + + /*! + * \brief The difference between the pressures of the gas and the non-wetting phase. + * + * This method is only available for at least three fluid phases + */ + template + static typename std::enable_if< (Traits::numPhases > 2), ScalarT>::type + pcng(const Params ¶ms, const FluidState &fs) + { return 0; } +}; +} // namespace Opm + +#endif diff --git a/opm/material/fluidmatrixinteractions/NullMaterialLaw.hpp b/opm/material/fluidmatrixinteractions/NullMaterialLaw.hpp deleted file mode 100644 index 8c22704c9..000000000 --- a/opm/material/fluidmatrixinteractions/NullMaterialLaw.hpp +++ /dev/null @@ -1,91 +0,0 @@ -// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- -// vi: set et ts=4 sw=4 sts=4: -/***************************************************************************** - * Copyright (C) 2011-2012 by Andreas Lauser * - * * - * 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 . * - *****************************************************************************/ -/*! - * \file - * \copydoc Opm::NullMaterialLaw - */ -#ifndef OPM_NULL_MATERIAL_LAW_HH -#define OPM_NULL_MATERIAL_LAW_HH - -#include "NullMaterialLawParams.hpp" - -#include -#include - -#include - -namespace Opm -{ -/*! - * \ingroup material - * - * \brief Implements a dummy linear saturation-capillary pressure - * relation which just disables capillary pressure. - */ -template > -class NullMaterialLaw -{ -public: - typedef ParamsT Params; - typedef typename Params::Scalar Scalar; - enum { numPhases = numPhasesV }; - - /*! - * \brief Returns constant 0 for all phases. - * - * \param values Container for the return values - * \param params Parameters - * \param state The fluid state - */ - template - static void capillaryPressures(ContainerT &values, - const Params ¶ms, - const FluidState &state) - { - for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) - values[phaseIdx] = 0.0; - } - - /*! - * \brief The inverse of the capillary pressure - */ - template - static void saturations(ContainerT &values, - const Params ¶ms, - const FluidState &state) - { - OPM_THROW(std::runtime_error, "Not implemented: MpLinearMaterial::saturations()"); - } - - /*! - * \brief The relative permeability of all phases. - */ - template - static void relativePermeabilities(ContainerT &values, - const Params ¶ms, - const FluidState &state) - { - for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) { - values[phaseIdx] = std::max(std::min(state.saturation(phaseIdx),1.0),0.0); - } - } -}; -} // namespace Opm - -#endif diff --git a/opm/material/fluidmatrixinteractions/NullMaterialLawParams.hpp b/opm/material/fluidmatrixinteractions/NullMaterialParams.hpp similarity index 82% rename from opm/material/fluidmatrixinteractions/NullMaterialLawParams.hpp rename to opm/material/fluidmatrixinteractions/NullMaterialParams.hpp index 26f0c7c1b..68d1bdf25 100644 --- a/opm/material/fluidmatrixinteractions/NullMaterialLawParams.hpp +++ b/opm/material/fluidmatrixinteractions/NullMaterialParams.hpp @@ -18,24 +18,29 @@ *****************************************************************************/ /*! * \file - * \copydoc Opm::NullMaterialLawParams + * \copydoc Opm::NullMaterialParams */ -#ifndef OPM_NULL_MATERIAL_LAW_PARAMS_HH -#define OPM_NULL_MATERIAL_LAW_PARAMS_HH +#ifndef OPM_NULL_MATERIAL_PARAMS_HH +#define OPM_NULL_MATERIAL_PARAMS_HH namespace Opm { /*! * \brief Reference implementation of params for the linear M-phase * material material. */ -template -class NullMaterialLawParams +template +class NullMaterialParams { public: - typedef ScalarT Scalar; - enum { numPhases = numPhasesV }; + typedef typename TraitsT::Scalar Scalar; - NullMaterialLawParams() + NullMaterialParams() + { } + + /*! + * \brief Finish the construction of the parameter object. + */ + void finalize() { } }; } // namespace Opm diff --git a/opm/material/fluidmatrixinteractions/ParkerLenhard.hpp b/opm/material/fluidmatrixinteractions/ParkerLenhard.hpp index d7b89e5af..bf7417d4d 100644 --- a/opm/material/fluidmatrixinteractions/ParkerLenhard.hpp +++ b/opm/material/fluidmatrixinteractions/ParkerLenhard.hpp @@ -361,6 +361,132 @@ public: values[Traits::nPhaseIdx] = krn(params, fs); } + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = 0; + if (satPhaseIdx == Traits::wPhaseIdx) + values[Traits::nPhaseIdx] = twoPhaseSatDpcwn_dSw(params, state.saturation(Traits::wPhaseIdx)); + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + if (satPhaseIdx == Traits::wPhaseIdx) { + values[Traits::wPhaseIdx] = twoPhaseSatDKrw_dSw(params, state.saturation(Traits::wPhaseIdx)); + values[Traits::nPhaseIdx] = 0; + } + else { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = - twoPhaseSatDKrn_dSw(params, 1 - state.saturation(Traits::nPhaseIdx)); + } + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + /*! * \brief Returns the capillary pressure dependend on * the phase saturations. @@ -402,6 +528,9 @@ public: } } + static Scalar twoPhaseSatDpcwn_dSw(const Params ¶ms, Scalar Sw) + { OPM_THROW(std::logic_error, "Not implemented: twoPhaseSatDpcwn_dSw()"); } + /*! * \brief Calculate the wetting phase saturations depending on * the phase pressures. @@ -440,6 +569,9 @@ public: return VanGenuchten::twoPhaseSatKrw(params.mdcParams(), Sw_app); } + static Scalar twoPhaseSatDKrw_dSw(const Params ¶ms, Scalar Sw) + { OPM_THROW(std::logic_error, "Not implemented: twoPhaseSatDKrw_dSw()"); } + /*! * \brief The relative permeability for the non-wetting phase * of the params. @@ -456,6 +588,9 @@ public: return VanGenuchten::twoPhaseSatKrn(params.mdcParams(), Sw_app); } + static Scalar twoPhaseSatDKrn_dSw(const Params ¶ms, Scalar Sw) + { OPM_THROW(std::logic_error, "Not implemented: twoPhaseSatDKrn_dSw()"); } + /*! * \brief Convert an absolute wetting saturation to an apparent one. */ diff --git a/opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp b/opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp index 124290add..2c8f9c61b 100644 --- a/opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp +++ b/opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp @@ -108,7 +108,7 @@ public: * \param values A random access container which stores the * relative pressure of each fluid phase. * \param params The parameter object expressing the coefficients - * required by the van Genuchten law. + * required by the material law. * \param fs The fluid state for which the capillary pressure * ought to be calculated */ @@ -136,7 +136,7 @@ public: * \param values A random access container which stores the * relative permeability of each fluid phase. * \param params The parameter object expressing the coefficients - * required by the van Genuchten law. + * required by the material law. * \param fs The fluid state for which the relative permeabilities * ought to be calculated */ @@ -147,6 +147,132 @@ public: values[Traits::nPhaseIdx] = krn(params, fs); } + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = 0; + if (satPhaseIdx == Traits::wPhaseIdx) + values[Traits::nPhaseIdx] = twoPhaseSatDpcwn_dSw(params, state.saturation(Traits::wPhaseIdx)); + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + if (satPhaseIdx == Traits::wPhaseIdx) { + values[Traits::wPhaseIdx] = twoPhaseSatDKrw_dSw(params, state.saturation(Traits::wPhaseIdx)); + values[Traits::nPhaseIdx] = 0; + } + else { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = - twoPhaseSatDKrn_dSw(params, 1 - state.saturation(Traits::nPhaseIdx)); + } + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + /*! * \brief A regularized Brooks-Corey capillary pressure-saturation * curve. @@ -263,30 +389,30 @@ public: * \brief The derivative of the regularized Brooks-Corey capillary * pressure-saturation curve. */ - static Scalar dpcwn_dSw(const Params ¶ms, Scalar Sw) + static Scalar twoPhaseSatDpcwn_dSw(const Params ¶ms, Scalar Sw) { const Scalar Sthres = params.thresholdSw(); // derivative of the regualarization if (Sw <= Sthres) { // calculate the slope of the straight line used in pcwn() - Scalar m = BrooksCorey::dpcwn_dSw(params, Sthres); + Scalar m = BrooksCorey::twoPhaseSatDpcwn_dSw(params, Sthres); return m; } else if (Sw > 1.0) { // calculate the slope of the straight line used in pcwn() - Scalar m = BrooksCorey::dpcwn_dSw(params, 1.0); + Scalar m = BrooksCorey::twoPhaseSatDpcwn_dSw(params, 1.0); return m; } - return BrooksCorey::dpcwn_dSw(params, Sw); + return BrooksCorey::twoPhaseSatDpcwn_dSw(params, Sw); } /*! * \brief The derivative of the regularized Brooks-Corey * saturation-capillary pressure curve. */ - static Scalar dSw_dpcwn(const Params ¶ms, Scalar pcwn) + static Scalar twoPhaseSatDSw_dpcwn(const Params ¶ms, Scalar pcwn) { const Scalar Sthres = params.thresholdSw(); @@ -341,6 +467,14 @@ public: return BrooksCorey::twoPhaseSatKrw(params, Sw); } + static Scalar twoPhaseSatDKrw_dSw(const Params ¶ms, Scalar Sw) + { + if (Sw <= 0.0 || Sw >= 1.0) + return 0.0; + + return BrooksCorey::twoPhaseSatDKrw_dSw(params, Sw); + } + /*! * \brief Regularized version of the relative permeability of the * non-wetting phase of the Brooks-Corey curves. @@ -357,7 +491,7 @@ public: */ template static Scalar krn(const Params ¶ms, const FluidState &fs) - { return twoPhaseSatKrn(params, fs.saturation(Traits::wPhaseIdx)); } + { return twoPhaseSatKrn(params, 1.0 - fs.saturation(Traits::nPhaseIdx)); } static Scalar twoPhaseSatKrn(const Params ¶ms, Scalar Sw) { @@ -368,6 +502,15 @@ public: return BrooksCorey::twoPhaseSatKrn(params, Sw); } + + static Scalar twoPhaseSatDKrn_dSw(const Params ¶ms, Scalar Sw) + { + if (Sw <= 0.0 || Sw >= 1.0) + return 0.0; + + return BrooksCorey::twoPhaseSatDKrn_dSw(params, Sw); + } + }; } // namespace Opm diff --git a/opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp b/opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp index 936da2ace..e7b780285 100644 --- a/opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp +++ b/opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp @@ -139,6 +139,133 @@ public: values[Traits::nPhaseIdx] = krn(params, fs); } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = 0; + if (satPhaseIdx == Traits::wPhaseIdx) + values[Traits::nPhaseIdx] = twoPhaseSatDpcwn_dSw(params, state.saturation(Traits::wPhaseIdx)); + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + if (satPhaseIdx == Traits::wPhaseIdx) { + values[Traits::wPhaseIdx] = twoPhaseSatDKrw_dSw(params, state.saturation(Traits::wPhaseIdx)); + values[Traits::nPhaseIdx] = 0; + } + else { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = - twoPhaseSatDKrn_dSw(params, 1 - state.saturation(Traits::nPhaseIdx)); + } + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + /*! * \brief A regularized van Genuchten capillary pressure-saturation * curve. @@ -292,10 +419,10 @@ public: else // the slope of the straight line used for the right // side of the capillary pressure function - return params.mHigh(); + return params.pcwnSlopeHigh(); } - return VanGenuchten::dpcwn_dSw(params, Sw); + return VanGenuchten::twoPhaseSatDpcwn_dSw(params, Sw); } /*! @@ -369,6 +496,14 @@ public: return VanGenuchten::twoPhaseSatKrw(params, Sw); } + static Scalar twoPhaseSatDKrw_dSw(const Params ¶ms, Scalar Sw) + { + if (Sw <= 0.0 || Sw >= 1.0) + return 0.0; + + return VanGenuchten::twoPhaseSatDKrw_dSw(params, Sw); + } + /*! * \brief Regularized version of the relative permeability * for the non-wetting phase of @@ -385,7 +520,7 @@ public: */ template static Scalar krn(const Params ¶ms, const FluidState &fs) - { return twoPhaseSatKrn(params, fs.saturation(Traits::wPhaseIdx)); } + { return twoPhaseSatKrn(params, 1.0 - fs.saturation(Traits::nPhaseIdx)); } static Scalar twoPhaseSatKrn(const Params ¶ms, Scalar Sw) { @@ -397,6 +532,14 @@ public: return VanGenuchten::twoPhaseSatKrn(params, Sw); } + + static Scalar twoPhaseSatDKrn_dSw(const Params ¶ms, Scalar Sw) + { + if (Sw <= 0.0 || Sw >= 1.0) + return 0.0; + + return VanGenuchten::twoPhaseSatDKrn_dSw(params, Sw); + } }; } // namespace Opm diff --git a/opm/material/fluidmatrixinteractions/VanGenuchten.hpp b/opm/material/fluidmatrixinteractions/VanGenuchten.hpp index a26f1cac2..0247d79b7 100644 --- a/opm/material/fluidmatrixinteractions/VanGenuchten.hpp +++ b/opm/material/fluidmatrixinteractions/VanGenuchten.hpp @@ -141,6 +141,133 @@ public: values[Traits::nPhaseIdx] = krn(params, fs); } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase saturation. + */ + template + static void dCapillaryPressures_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = 0; + if (satPhaseIdx == Traits::wPhaseIdx) + values[Traits::nPhaseIdx] = dpcwn_dSw(params, state); + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given phase pressure. + */ + template + static void dCapillaryPressures_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * temperature. + */ + template + static void dCapillaryPressures_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all capillary pressures in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dCapillaryPressures_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int pcPhaseIdx = 0; pcPhaseIdx < numPhases; ++pcPhaseIdx) + values[pcPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase saturation. + */ + template + static void dRelativePermeabilities_dSaturation(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int satPhaseIdx) + { + if (satPhaseIdx == Traits::wPhaseIdx) { + values[Traits::wPhaseIdx] = twoPhaseSatDKrw_dSw(params, state.saturation(Traits::wPhaseIdx)); + values[Traits::nPhaseIdx] = 0; + } + else { + values[Traits::wPhaseIdx] = 0; + values[Traits::nPhaseIdx] = - twoPhaseSatDKrn_dSw(params, 1 - state.saturation(Traits::nPhaseIdx)); + } + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given phase pressure. + */ + template + static void dRelativePermeabilities_dPressure(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int pPhaseIdx) + { + // -> not pressure dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * temperature. + */ + template + static void dRelativePermeabilities_dTemperature(ContainerT &values, + const Params ¶ms, + const FluidState &state) + { + // -> not temperature dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + + /*! + * \brief The derivative of all relative permeabilities in regard to + * a given mole fraction of a component in a phase. + */ + template + static void dRelativePermeabilities_dMoleFraction(ContainerT &values, + const Params ¶ms, + const FluidState &state, + int phaseIdx, + int compIdx) + { + // -> not composition dependent + for (int krPhaseIdx = 0; krPhaseIdx < numPhases; ++krPhaseIdx) + values[krPhaseIdx] = 0.0; + } + /*! * \brief The capillary pressure-saturation curve according to van Genuchten. * @@ -280,10 +407,10 @@ public: * ought to be calculated */ template - static Scalar dkrw_dSw(const Params ¶ms, const FluidState &fs) + static Scalar dKrw_dSw(const Params ¶ms, const FluidState &fs) { return twoPhaseSatDkrw_dSw(params, fs.saturation(Traits::wPhaseIdx)); } - static Scalar twoPhaseSatDkrw_dSw(const Params ¶ms, Scalar Sw) + static Scalar twoPhaseSatDKrw_dSw(const Params ¶ms, Scalar Sw) { assert(0 <= Sw && Sw <= 1); @@ -304,7 +431,7 @@ public: */ template static Scalar krn(const Params ¶ms, const FluidState &fs) - { return twoPhaseSatKrn(params, fs.saturation(Traits::wPhaseIdx)); } + { return twoPhaseSatKrn(params, 1.0 - fs.saturation(Traits::nPhaseIdx)); } static Scalar twoPhaseSatKrn(const Params ¶ms, Scalar Sw) { @@ -327,10 +454,10 @@ public: * ought to be calculated */ template - static Scalar dkrn_dSw(const Params ¶ms, const FluidState &fs) + static Scalar dKrn_dSw(const Params ¶ms, const FluidState &fs) { return twoPhaseSatDkrn_dSw(params, fs.saturation(Traits::wPhaseIdx)); } - static Scalar twoPhaseSatDkrn_dSw(const Params ¶ms, Scalar Sw) + static Scalar twoPhaseSatDKrn_dSw(const Params ¶ms, Scalar Sw) { assert(0 <= Sw && Sw <= 1); diff --git a/tests/test_fluidmatrixinteractions.cpp b/tests/test_fluidmatrixinteractions.cpp index 9705f2c63..93c0be77f 100644 --- a/tests/test_fluidmatrixinteractions.cpp +++ b/tests/test_fluidmatrixinteractions.cpp @@ -30,6 +30,7 @@ #include #include #include +#include #include #include #include @@ -121,6 +122,18 @@ void testGenericApi() MaterialLaw::capillaryPressures(destValues, paramsConst, fs); MaterialLaw::saturations(destValues, paramsConst, fs); MaterialLaw::relativePermeabilities(destValues, paramsConst, fs); + + std::array dpc; + MaterialLaw::dCapillaryPressures_dSaturation(dpc, paramsConst, fs, /*phaseIdx=*/0); + MaterialLaw::dCapillaryPressures_dPressure(dpc, paramsConst, fs, /*phaseIdx=*/0); + MaterialLaw::dCapillaryPressures_dTemperature(dpc, paramsConst, fs); + MaterialLaw::dCapillaryPressures_dMoleFraction(dpc, paramsConst, fs, /*phaseIdx=*/0, /*compIdx=*/0); + + std::array dkr; + MaterialLaw::dRelativePermeabilities_dSaturation(dkr, paramsConst, fs, /*phaseIdx=*/0); + MaterialLaw::dRelativePermeabilities_dPressure(dkr, paramsConst, fs, /*phaseIdx=*/0); + MaterialLaw::dRelativePermeabilities_dTemperature(dkr, paramsConst, fs); + MaterialLaw::dRelativePermeabilities_dMoleFraction(dkr, paramsConst, fs, /*phaseIdx=*/0, /*compIdx=*/0); } } @@ -271,6 +284,17 @@ int main(int argc, char **argv) testThreePhaseApi(); //testThreePhaseSatApi(); } + { + typedef Opm::NullMaterial MaterialLaw; + testGenericApi(); + testTwoPhaseApi(); + testTwoPhaseSatApi(); + + typedef Opm::NullMaterial ThreePMaterialLaw; + testGenericApi(); + testThreePhaseApi(); + //testThreePhaseSatApi(); + } { typedef Opm::ParkerLenhard MaterialLaw; testGenericApi(); @@ -296,6 +320,5 @@ int main(int argc, char **argv) testTwoPhaseSatApi(); } - return 0; }