OilfieldToolsNet/opm-core
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Added class SaturationPropsFromDeck.

Browse Source
This commit is contained in:
Atgeirr Flø Rasmussen 2012-01-05 11:17:52 +01:00
parent 91e7b75a2d
commit 536d58fe4e
3 changed files with 350 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
Makefile.am
View File

@ -28,6 +28,7 @@ opm/core/fluid/blackoil/SinglePvtLiveGas.cpp \
opm/core/fluid/blackoil/SinglePvtLiveOil.cpp \
opm/core/fluid/blackoil/SinglePvtInterface.cpp \
opm/core/fluid/BlackoilPropertiesFromDeck.cpp \
opm/core/fluid/SaturationPropsFromDeck.cpp \
opm/core/utility/MonotCubicInterpolator.cpp \
opm/core/utility/parameters/Parameter.cpp \
opm/core/utility/parameters/ParameterGroup.cpp \
@ -85,6 +86,7 @@ opm/core/fluid/BlackoilFluid.hpp \
opm/core/fluid/BlackoilPropertiesInterface.hpp \
opm/core/fluid/BlackoilPropertiesFromDeck.hpp \
opm/core/fluid/RockFromDeck.hpp \
opm/core/fluid/SaturationPropsFromDeck.hpp \
opm/core/fluid/SimpleFluid2p.hpp \
opm/core/fluid/blackoil/phaseUsageFromDeck.hpp \
opm/core/fluid/blackoil/BlackoilDefs.hpp \

257
opm/core/fluid/SaturationPropsFromDeck.cpp Normal file
View File

@ -0,0 +1,257 @@
/*
Copyright 2010, 2011, 2012 SINTEF ICT, Applied Mathematics.
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 3 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/>.
*/
#include <opm/core/fluid/SaturationPropsFromDeck.hpp>
#include <opm/core/fluid/blackoil/phaseUsageFromDeck.hpp>
#include <opm/core/utility/buildUniformMonotoneTable.hpp>
#include <opm/core/utility/ErrorMacros.hpp>
#include <iostream>
namespace Opm
{
/// Default constructor.
SaturationPropsFromDeck::SaturationPropsFromDeck()
{
}
/// Initialize from deck.
void SaturationPropsFromDeck::init(const Dune::EclipseGridParser& deck)
{
phase_usage_ = phaseUsageFromDeck(deck);
// Extract input data.
// Oil phase should be active.
if (!phase_usage_.phase_used[Liquid]) {
THROW("SaturationPropsFromDeck::init() -- oil phase must be active.");
}
const int samples = 200;
if (phase_usage_.phase_used[Aqua]) {
const Dune::SWOF::table_t& swof_table = deck.getSWOF().swof_;
if (swof_table.size() != 1) {
THROW("We must have exactly one SWOF table.");
}
const std::vector<double>& sw = swof_table[0][0];
const std::vector<double>& krw = swof_table[0][1];
const std::vector<double>& krow = swof_table[0][2];
const std::vector<double>& pcow = swof_table[0][3];
buildUniformMonotoneTable(sw, krw, samples, krw_);
buildUniformMonotoneTable(sw, krow, samples, krow_);
buildUniformMonotoneTable(sw, pcow, samples, pcow_);
krocw_ = krow[0]; // At connate water -> ecl. SWOF
}
if (phase_usage_.phase_used[Vapour]) {
const Dune::SGOF::table_t& sgof_table = deck.getSGOF().sgof_;
if (sgof_table.size() != 1) {
THROW("We must have exactly one SGOF table.");
}
const std::vector<double>& sg = sgof_table[0][0];
const std::vector<double>& krg = sgof_table[0][1];
const std::vector<double>& krog = sgof_table[0][2];
const std::vector<double>& pcog = sgof_table[0][3];
buildUniformMonotoneTable(sg, krg, samples, krg_);
buildUniformMonotoneTable(sg, krog, samples, krog_);
buildUniformMonotoneTable(sg, pcog, samples, pcog_);
}
}
/// Relative permeability.
/// \param[in] n Number of data points.
/// \param[in] s Array of nP saturation values.
/// \param[out] kr Array of nP relperm values, array must be valid before calling.
/// \param[out] dkrds If non-null: array of nP^2 relperm derivative values,
/// array must be valid before calling.
/// The P^2 derivative matrix is
/// m_{ij} = \frac{dkr_i}{ds^j},
/// and is output in Fortran order (m_00 m_10 m_20 m01 ...)
void SaturationPropsFromDeck::relperm(const int n,
const double* s,
double* kr,
double* dkrds) const
{
const int np = phase_usage_.num_phases;
if (dkrds) {
#pragma omp parallel for
for (int i = 0; i < n; ++i) {
evalKrDeriv(s + np*i, kr + np*i, dkrds + np*np*i);
}
} else {
#pragma omp parallel for
for (int i = 0; i < n; ++i) {
evalKr(s + np*i, kr + np*i);
}
}
}
/// Capillary pressure.
/// \param[in] n Number of data points.
/// \param[in] s Array of nP saturation values.
/// \param[out] pc Array of nP capillary pressure values, array must be valid before calling.
/// \param[out] dpcds If non-null: array of nP^2 derivative values,
/// array must be valid before calling.
/// The P^2 derivative matrix is
/// m_{ij} = \frac{dpc_i}{ds^j},
/// and is output in Fortran order (m_00 m_10 m_20 m01 ...)
void SaturationPropsFromDeck::capPress(const int n,
const double* s,
double* pc,
double* dpcds) const
{
const int np = phase_usage_.num_phases;
if (dpcds) {
#pragma omp parallel for
for (int i = 0; i < n; ++i) {
evalPcDeriv(s + np*i, pc + np*i, dpcds + np*np*i);
}
} else {
#pragma omp parallel for
for (int i = 0; i < n; ++i) {
evalPc(s + np*i, pc + np*i);
}
}
}
void SaturationPropsFromDeck::evalKr(const double* s, double* kr) const
{
if (phase_usage_.num_phases == 3) {
// Stone-II relative permeability model.
double sw = s[Aqua];
double sg = s[Vapour];
double krw = krw_(sw);
double krg = krg_(sg);
double krow = krow_(sw + sg); // = 1 - so
double krog = krog_(sg); // = 1 - so - sw
double krocw = krocw_;
kr[Aqua] = krw;
kr[Vapour] = krg;
kr[Liquid] = krocw*((krow/krocw + krw)*(krog/krocw + krg) - krw - krg);
if (kr[Liquid] < 0.0) {
kr[Liquid] = 0.0;
}
return;
}
// We have a two-phase situation. We know that oil is active.
if (phase_usage_.phase_used[Aqua]) {
int wpos = phase_usage_.phase_pos[Aqua];
int opos = phase_usage_.phase_pos[Liquid];
double sw = s[wpos];
double krw = krw_(sw);
double krow = krow_(sw);
kr[wpos] = krw;
kr[opos] = krow;
} else {
ASSERT(phase_usage_.phase_used[Vapour]);
int gpos = phase_usage_.phase_pos[Vapour];
int opos = phase_usage_.phase_pos[Liquid];
double sg = s[gpos];
double krg = krg_(sg);
double krog = krog_(sg);
kr[gpos] = krg;
kr[opos] = krog;
}
}
void SaturationPropsFromDeck::evalKrDeriv(const double* s, double* kr, double* dkrds) const
{
const int np = phase_usage_.num_phases;
std::fill(dkrds, dkrds + np*np, 0.0);
if (np == 3) {
// Stone-II relative permeability model.
double sw = s[Aqua];
double sg = s[Vapour];
double krw = krw_(sw);
double dkrww = krw_.derivative(sw);
double krg = krg_(sg);
double dkrgg = krg_.derivative(sg);
double krow = krow_(sw + sg);
double dkrow = krow_.derivative(sw + sg);
double krog = krog_(sg);
double dkrog = krog_.derivative(sg);
double krocw = krocw_;
kr[Aqua] = krw;
kr[Vapour] = krg;
kr[Liquid] = krocw*((krow/krocw + krw)*(krog/krocw + krg) - krw - krg);
if (kr[Liquid] < 0.0) {
kr[Liquid] = 0.0;
}
dkrds[Aqua + Aqua*np] = dkrww;
dkrds[Vapour + Vapour*np] = dkrgg;
dkrds[Liquid + Aqua*np] = krocw*((dkrow/krocw + dkrww)*(krog/krocw + krg) - dkrww);
dkrds[Liquid + Vapour*np] = krocw*((krow/krocw + krw)*(dkrog/krocw + dkrgg) - dkrgg)
+ krocw*((dkrow/krocw + krw)*(krog/krocw + krg) - dkrgg);
return;
}
// We have a two-phase situation. We know that oil is active.
if (phase_usage_.phase_used[Aqua]) {
int wpos = phase_usage_.phase_pos[Aqua];
int opos = phase_usage_.phase_pos[Liquid];
double sw = s[wpos];
double krw = krw_(sw);
double dkrww = krw_.derivative(sw);
double krow = krow_(sw);
double dkrow = krow_.derivative(sw);
kr[wpos] = krw;
kr[opos] = krow;
dkrds[wpos + wpos*np] = dkrww;
dkrds[opos + wpos*np] = dkrow; // Row opos, column wpos, fortran order.
} else {
ASSERT(phase_usage_.phase_used[Vapour]);
int gpos = phase_usage_.phase_pos[Vapour];
int opos = phase_usage_.phase_pos[Liquid];
double sg = s[gpos];
double krg = krg_(sg);
double dkrgg = krg_.derivative(sg);
double krog = krog_(sg);
double dkrog = krog_.derivative(sg);
kr[gpos] = krg;
kr[opos] = krog;
dkrds[gpos + gpos*np] = dkrgg;
dkrds[opos + gpos*np] = dkrog;
}
}
void SaturationPropsFromDeck::evalPc(const double* s, double* pc) const
{
pc[phase_usage_.phase_pos[Liquid]] = 0.0;
if (phase_usage_.phase_used[Aqua]) {
int pos = phase_usage_.phase_pos[Aqua];
pc[pos] = pcow_(s[pos]);
}
if (phase_usage_.phase_used[Vapour]) {
int pos = phase_usage_.phase_pos[Vapour];
pc[pos] = pcog_(s[pos]);
}
}
void SaturationPropsFromDeck::evalPcDeriv(const double* /*s*/, double* /*pc*/, double* /*dpcds*/) const
{
THROW("Evaluation of capillary pressure derivatives not yet implemented.");
}
} // namespace Opm

91
opm/core/fluid/SaturationPropsFromDeck.hpp Normal file
View File

@ -0,0 +1,91 @@
/*
Copyright 2010, 2011, 2012 SINTEF ICT, Applied Mathematics.
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 3 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/>.
*/
#ifndef OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED
#define OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED
#include <opm/core/eclipse/EclipseGridParser.hpp>
#include <opm/core/utility/UniformTableLinear.hpp>
#include <opm/core/fluid/blackoil/BlackoilPhases.hpp>
namespace Opm
{
class SaturationPropsFromDeck : public BlackoilPhases
{
public:
/// Default constructor.
SaturationPropsFromDeck();
/// Initialize from deck.
void init(const Dune::EclipseGridParser& deck);
/// Relative permeability.
/// \param[in] n Number of data points.
/// \param[in] s Array of nP saturation values.
/// \param[out] kr Array of nP relperm values, array must be valid before calling.
/// \param[out] dkrds If non-null: array of nP^2 relperm derivative values,
/// array must be valid before calling.
/// The P^2 derivative matrix is
/// m_{ij} = \frac{dkr_i}{ds^j},
/// and is output in Fortran order (m_00 m_10 m_20 m01 ...)
void relperm(const int n,
const double* s,
double* kr,
double* dkrds) const;
/// Capillary pressure.
/// \param[in] n Number of data points.
/// \param[in] s Array of nP saturation values.
/// \param[out] pc Array of nP capillary pressure values, array must be valid before calling.
/// \param[out] dpcds If non-null: array of nP^2 derivative values,
/// array must be valid before calling.
/// The P^2 derivative matrix is
/// m_{ij} = \frac{dpc_i}{ds^j},
/// and is output in Fortran order (m_00 m_10 m_20 m01 ...)
void capPress(const int n,
const double* s,
double* pv,
double* dpcds) const;
private:
void evalKr(const double* s, double* kr) const;
void evalKrDeriv(const double* s, double* kr, double* dkrds) const;
void evalPc(const double* s, double* pc) const;
void evalPcDeriv(const double* s, double* pc, double* dpcds) const;
private:
PhaseUsage phase_usage_;
Dune::utils::UniformTableLinear<double> krw_;
Dune::utils::UniformTableLinear<double> krow_;
Dune::utils::UniformTableLinear<double> pcow_;
Dune::utils::UniformTableLinear<double> krg_;
Dune::utils::UniformTableLinear<double> krog_;
Dune::utils::UniformTableLinear<double> pcog_;
double krocw_; // = krow_(s_wc)
};
} // namespace Opm
#endif // OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED