411 lines
15 KiB
C++
411 lines
15 KiB
C++
/*
|
|
Copyright (C) 2008-2013 by Andreas Lauser
|
|
|
|
This file is part of the Open Porous Media project (OPM).
|
|
|
|
OPM 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.
|
|
|
|
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
|
|
*/
|
|
/*!
|
|
* \file
|
|
* \copydoc Opm::EclDefaultMaterial
|
|
*/
|
|
#ifndef OPM_ECL_DEFAULT_MATERIAL_HPP
|
|
#define OPM_ECL_DEFAULT_MATERIAL_HPP
|
|
|
|
#include "EclDefaultMaterialParams.hpp"
|
|
|
|
#include <opm/material/Valgrind.hpp>
|
|
|
|
#include <opm/core/utility/Exceptions.hpp>
|
|
#include <opm/core/utility/ErrorMacros.hpp>
|
|
|
|
#include <algorithm>
|
|
|
|
namespace Opm {
|
|
|
|
/*!
|
|
* \ingroup material
|
|
*
|
|
* \brief Implements the default three phase capillary pressure law
|
|
* used by the ECLipse simulator.
|
|
*
|
|
* This material law is valid for three fluid phases and only depends
|
|
* on the saturations.
|
|
*
|
|
* The required two-phase relations are supplied by means of template
|
|
* arguments and can be an arbitrary other material laws. (Provided
|
|
* that these only depend on saturation.)
|
|
*/
|
|
template <class TraitsT,
|
|
class GasOilMaterialLawT,
|
|
class OilWaterMaterialLawT,
|
|
class ParamsT = EclDefaultMaterialParams<TraitsT,
|
|
typename GasOilMaterialLawT::Params,
|
|
typename OilWaterMaterialLawT::Params> >
|
|
class EclDefaultMaterial : public TraitsT
|
|
{
|
|
public:
|
|
typedef GasOilMaterialLawT GasOilMaterialLaw;
|
|
typedef OilWaterMaterialLawT OilWaterMaterialLaw;
|
|
|
|
// some safety checks
|
|
static_assert(TraitsT::numPhases == 3,
|
|
"The number of phases considered by this capillary pressure "
|
|
"law is always three!");
|
|
static_assert(GasOilMaterialLaw::numPhases == 2,
|
|
"The number of phases considered by the gas-oil capillary "
|
|
"pressure law must be two!");
|
|
static_assert(OilWaterMaterialLaw::numPhases == 2,
|
|
"The number of phases considered by the oil-water capillary "
|
|
"pressure law must be two!");
|
|
static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
|
|
typename OilWaterMaterialLaw::Scalar>::value,
|
|
"The two two-phase capillary pressure laws must use the same "
|
|
"type of floating point values.");
|
|
|
|
static_assert(GasOilMaterialLaw::implementsTwoPhaseSatApi,
|
|
"The gas-oil material law must implement the two-phase saturation "
|
|
"only API to for the default Ecl capillary pressure law!");
|
|
static_assert(OilWaterMaterialLaw::implementsTwoPhaseSatApi,
|
|
"The oil-water material law must implement the two-phase saturation "
|
|
"only API to for the default Ecl capillary pressure law!");
|
|
|
|
typedef TraitsT Traits;
|
|
typedef ParamsT Params;
|
|
typedef typename Traits::Scalar Scalar;
|
|
|
|
static const int numPhases = 3;
|
|
static const int waterPhaseIdx = Traits::wettingPhaseIdx;
|
|
static const int nonWettingPhaseIdx = Traits::nonWettingPhaseIdx;
|
|
static const int oilPhaseIdx = Traits::nonWettingPhaseIdx;
|
|
static const int gasPhaseIdx = Traits::gasPhaseIdx;
|
|
|
|
//! Specify whether this material law implements the two-phase
|
|
//! convenience API
|
|
static const bool implementsTwoPhaseApi = false;
|
|
|
|
//! Specify whether this material law implements the two-phase
|
|
//! convenience API which only depends on the phase saturations
|
|
static const bool implementsTwoPhaseSatApi = false;
|
|
|
|
//! Specify whether the quantities defined by this material law
|
|
//! are saturation dependent
|
|
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 Implements the default three phase capillary pressure law
|
|
* used by the ECLipse simulator.
|
|
*
|
|
* This material law is valid for three fluid phases and only
|
|
* depends on the saturations.
|
|
*
|
|
* The required two-phase relations are supplied by means of template
|
|
* arguments and can be an arbitrary other material laws.
|
|
*
|
|
* \param values Container for the return values
|
|
* \param params Parameters
|
|
* \param state The fluid state
|
|
*/
|
|
template <class ContainerT, class FluidState>
|
|
static void capillaryPressures(ContainerT &values,
|
|
const Params ¶ms,
|
|
const FluidState &state)
|
|
{
|
|
values[gasPhaseIdx] = pcgn(params, state);
|
|
values[oilPhaseIdx] = 0;
|
|
values[waterPhaseIdx] = - pcnw(params, state);
|
|
}
|
|
|
|
/*!
|
|
* \brief Capillary pressure between the gas and the non-wetting
|
|
* liquid (i.e., oil) phase.
|
|
*
|
|
* This is defined as
|
|
* \f[
|
|
* p_{c,gn} = p_g - p_n
|
|
* \f]
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar pcgn(const Params ¶ms,
|
|
const FluidState &fs)
|
|
{
|
|
Scalar Sw = 1 - fs.saturation(gasPhaseIdx);
|
|
return GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), Sw);
|
|
}
|
|
|
|
/*!
|
|
* \brief Capillary pressure between the non-wetting liquid (i.e.,
|
|
* oil) and the wetting liquid (i.e., water) phase.
|
|
*
|
|
* This is defined as
|
|
* \f[
|
|
* p_{c,nw} = p_n - p_w
|
|
* \f]
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar pcnw(const Params ¶ms,
|
|
const FluidState &fs)
|
|
{
|
|
Scalar Sw = fs.saturation(waterPhaseIdx);
|
|
return OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(), Sw);
|
|
}
|
|
|
|
/*!
|
|
* \brief The inverse of the capillary pressure
|
|
*/
|
|
template <class ContainerT, class FluidState>
|
|
static void saturations(ContainerT &values,
|
|
const Params ¶ms,
|
|
const FluidState &fs)
|
|
{
|
|
OPM_THROW(std::logic_error, "Not implemented: saturations()");
|
|
}
|
|
|
|
/*!
|
|
* \brief The saturation of the gas phase.
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar Sg(const Params ¶ms,
|
|
const FluidState &fluidState)
|
|
{
|
|
OPM_THROW(std::logic_error, "Not implemented: Sg()");
|
|
}
|
|
|
|
/*!
|
|
* \brief The saturation of the non-wetting (i.e., oil) phase.
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar Sn(const Params ¶ms,
|
|
const FluidState &fluidState)
|
|
{
|
|
OPM_THROW(std::logic_error, "Not implemented: Sn()");
|
|
}
|
|
|
|
/*!
|
|
* \brief The saturation of the wetting (i.e., water) phase.
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar Sw(const Params ¶ms,
|
|
const FluidState &fluidState)
|
|
{
|
|
OPM_THROW(std::logic_error, "Not implemented: Sw()");
|
|
}
|
|
|
|
/*!
|
|
* \brief The relative permeability of all phases.
|
|
*
|
|
* The relative permeability of the water phase it uses the same
|
|
* value as the relative permeability for water in the water-oil
|
|
* law with \f$S_o = 1 - S_w\f$. The gas relative permebility is
|
|
* taken from the gas-oil material law, but with \f$S_o = 1 -
|
|
* S_g\f$. The relative permeability of the oil phase is
|
|
* calculated using the relative permeabilities of the oil phase
|
|
* in the two two-phase systems.
|
|
*
|
|
* A more detailed description can be found in the "Three phase
|
|
* oil relative permeability models" section of the ECLipse
|
|
* technical description.
|
|
*/
|
|
template <class ContainerT, class FluidState>
|
|
static void relativePermeabilities(ContainerT &values,
|
|
const Params ¶ms,
|
|
const FluidState &fluidState)
|
|
{
|
|
values[waterPhaseIdx] = krw(params, fluidState);
|
|
values[nonWettingPhaseIdx] = krn(params, fluidState);
|
|
values[gasPhaseIdx] = krg(params, fluidState);
|
|
}
|
|
|
|
/*!
|
|
* \brief The relative permeability of the gas phase.
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar krg(const Params ¶ms,
|
|
const FluidState &fluidState)
|
|
{
|
|
Scalar Sw = 1 - fluidState.saturation(gasPhaseIdx);
|
|
return GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), Sw);
|
|
}
|
|
|
|
/*!
|
|
* \brief The relative permeability of the wetting phase.
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar krw(const Params ¶ms,
|
|
const FluidState &fluidState)
|
|
{
|
|
Scalar Sw = fluidState.saturation(waterPhaseIdx);
|
|
return OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
|
|
}
|
|
|
|
/*!
|
|
* \brief The relative permeability of the non-wetting (i.e., oil) phase.
|
|
*/
|
|
template <class FluidState>
|
|
static Scalar krn(const Params ¶ms,
|
|
const FluidState &fluidState)
|
|
{
|
|
Scalar Sw = std::min(1.0, std::max(0.0, fluidState.saturation(waterPhaseIdx)));
|
|
Scalar So = std::min(1.0, std::max(0.0, fluidState.saturation(oilPhaseIdx)));
|
|
Scalar Sg = std::min(1.0, std::max(0.0, fluidState.saturation(gasPhaseIdx)));
|
|
|
|
// connate water. According to the Eclipse TD, this is
|
|
// probably only relevant if hysteresis is enabled...
|
|
Scalar Swco = 0; // todo!
|
|
|
|
Scalar krog = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So + Swco);
|
|
Scalar krow = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), 1 - So);
|
|
|
|
if (Sg + Sw - Swco < 1e-30)
|
|
return 1.0; // avoid division by zero
|
|
else {
|
|
Scalar tmp = (Sg*krog + (Sw - Swco)*krow) / (Sg + Sw - Swco);
|
|
return std::min(1.0, std::max(0.0, tmp));
|
|
}
|
|
}
|
|
|
|
/*!
|
|
* \brief The derivative of all capillary pressures in regard to
|
|
* a given phase saturation.
|
|
*/
|
|
template <class ContainerT, class FluidState>
|
|
static void dCapillaryPressures_dSaturation(ContainerT &values,
|
|
const Params ¶ms,
|
|
const FluidState &fluidState,
|
|
int satPhaseIdx)
|
|
{
|
|
OPM_THROW(std::logic_error,
|
|
"Not implemented: dCapillaryPressure_dSaturation()");
|
|
}
|
|
|
|
/*!
|
|
* \brief The derivative of all capillary pressures in regard to
|
|
* a given phase pressure.
|
|
*/
|
|
template <class ContainerT, class FluidState>
|
|
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 <class ContainerT, class FluidState>
|
|
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 <class ContainerT, class FluidState>
|
|
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 <class ContainerT, class FluidState>
|
|
static void dRelativePermeabilities_dSaturation(ContainerT &values,
|
|
const Params ¶ms,
|
|
const FluidState &state,
|
|
int satPhaseIdx)
|
|
{
|
|
OPM_THROW(std::logic_error,
|
|
"Not implemented: dRelativePermeabilities_dSaturation()");
|
|
}
|
|
|
|
/*!
|
|
* \brief The derivative of all relative permeabilities in regard to
|
|
* a given phase pressure.
|
|
*/
|
|
template <class ContainerT, class FluidState>
|
|
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 <class ContainerT, class FluidState>
|
|
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 <class ContainerT, class FluidState>
|
|
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;
|
|
}
|
|
};
|
|
} // namespace Opm
|
|
|
|
#endif
|