Support for kw VAPPARS.

When this kw is active, BlackoilPropsAdFromDeck now modifies rvSat
and rsSat curves cell-wise by a power of (sat_oil_cell /
sat_oil_cell_historical_max).   Currently, the associated jacobians do
not reflect terms of type d/d_sat_oil, but code for doing this is given
as comments to BlackoilPropsAdFromDeck::applyVap(ADB& r, ...).
This commit is contained in:
osae 2014-07-05 15:06:12 +02:00
parent bb12bdd1fd
commit 4ce61b7c7c
8 changed files with 375 additions and 11 deletions

View File

@ -609,6 +609,23 @@ namespace Opm
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V BlackoilPropsAd::rsSat(const V& po,
const V& so,
const Cells& cells) const
{
// Suppress warning about "unused parameters".
static_cast<void>(po);
static_cast<void>(so);
static_cast<void>(cells);
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -623,6 +640,23 @@ namespace Opm
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB BlackoilPropsAd::rsSat(const ADB& po,
const ADB& so,
const Cells& cells) const
{
// Suppress warning about "unused parameters".
static_cast<void>(po);
static_cast<void>(so);
static_cast<void>(cells);
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
// ------ Rs bubble point curve ------
/// Bubble point curve for Rs as function of oil pressure.
@ -639,6 +673,23 @@ namespace Opm
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V BlackoilPropsAd::rvSat(const V& po,
const V& so,
const Cells& cells) const
{
// Suppress warning about "unused parameters".
static_cast<void>(po);
static_cast<void>(so);
static_cast<void>(cells);
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -653,6 +704,23 @@ namespace Opm
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB BlackoilPropsAd::rvSat(const ADB& po,
const ADB& so,
const Cells& cells) const
{
// Suppress warning about "unused parameters".
static_cast<void>(po);
static_cast<void>(so);
static_cast<void>(cells);
OPM_THROW(std::runtime_error, "Method rsMax() not implemented.");
}
// ------ Relative permeability ------
/// Relative permeabilities for all phases.
@ -828,7 +896,12 @@ namespace Opm
{
OPM_THROW(std::logic_error, "BlackoilPropsAd class does not support hysteresis.");
}
/// Update for max oil saturation.
void BlackoilPropsAd::updateSatOilMax(const std::vector<double>& saturation)
{
OPM_THROW(std::logic_error, "BlackoilPropsAd class does not support this functionality.");
}
} // namespace Opm

View File

@ -254,6 +254,15 @@ namespace Opm
V rsSat(const V& po,
const Cells& cells) const;
/// Solution gas/oil ratio and its derivatives at saturated condition as a function of p.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V rsSat(const V& po,
const V& so,
const Cells& cells) const;
/// Solution gas/oil ratio and its derivatives at saturated condition as a function of p.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -261,6 +270,15 @@ namespace Opm
ADB rsSat(const ADB& po,
const Cells& cells) const;
/// Solution gas/oil ratio and its derivatives at saturated condition as a function of p.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB rsSat(const ADB& po,
const ADB& so,
const Cells& cells) const;
// ------ Rv condensation curve ------
/// Vapor oil/gas ratio and its derivatives at saturated conditions as a function of p.
@ -270,6 +288,15 @@ namespace Opm
V rvSat(const V& po,
const Cells& cells) const;
/// Vapor oil/gas ratio and its derivatives at saturated conditions as a function of p.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V rvSat(const V& po,
const V& so,
const Cells& cells) const;
/// Vapor oil/gas ratio and its derivatives at saturated conditions as a function of p.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -277,6 +304,15 @@ namespace Opm
ADB rvSat(const ADB& po,
const Cells& cells) const;
/// Vapor oil/gas ratio and its derivatives at saturated conditions as a function of p.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB rvSat(const ADB& po,
const ADB& so,
const Cells& cells) const;
// ------ Relative permeability ------
/// Relative permeabilities for all phases.
@ -320,6 +356,11 @@ namespace Opm
/// \param[in] cells Array of n cell indices to be associated with the saturation values.
void updateSatHyst(const std::vector<double>& saturation,
const std::vector<int>& cells);
/// Update for max oil saturation.
void updateSatOilMax(const std::vector<double>& saturation);
private:
const BlackoilPropertiesInterface& props_;
PhaseUsage pu_;

View File

@ -188,6 +188,15 @@ namespace Opm
OPM_THROW(std::runtime_error, "Input is missing PVDG or PVTG\n");
}
}
// Oil vaporization controls (kw VAPPARS)
vap1_ = vap2_ = 0.0;
if (deck->hasKeyword("VAPPARS") && deck->hasKeyword("VAPOIL") && deck->hasKeyword("DISGAS")) {
vap1_ = deck->getKeyword("VAPPARS")->getRecord(0)->getItem(0)->getRawDouble(0);
vap2_ = deck->getKeyword("VAPPARS")->getRecord(0)->getItem(1)->getRawDouble(0);
satOilMax_.resize(number_of_cells, 0.0);
} else if (deck->hasKeyword("VAPPARS")) {
OPM_THROW(std::runtime_error, "Input has VAPPARS, but missing VAPOIL and/or DISGAS\n");
}
SaturationPropsFromDeck<SatFuncGwsegNonuniform>* ptr
= new SaturationPropsFromDeck<SatFuncGwsegNonuniform>();
@ -746,6 +755,20 @@ namespace Opm
return rbub;
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V BlackoilPropsAdFromDeck::rsSat(const V& po,
const V& so,
const Cells& cells) const
{
V rs = rsSat(po, cells);
applyVap(rs, so, cells, vap2_);
return rs;
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -770,6 +793,20 @@ namespace Opm
return ADB::function(rbub, jacs);
}
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB BlackoilPropsAdFromDeck::rsSat(const ADB& po,
const ADB& so,
const Cells& cells) const
{
ADB rs = rsSat(po, cells);
applyVap(rs, so, cells, vap2_);
return rs;
}
// ------ Condensation curve ------
/// Condensation curve for Rv as function of oil pressure.
@ -790,6 +827,20 @@ namespace Opm
return rv;
}
/// Condensation curve for Rv as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V BlackoilPropsAdFromDeck::rvSat(const V& po,
const V& so,
const Cells& cells) const
{
V rv = rvSat(po, cells);
applyVap(rv, so, cells, vap1_);
return rv;
}
/// Condensation curve for Rv as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -814,6 +865,20 @@ namespace Opm
return ADB::function(rv, jacs);
}
/// Condensation curve for Rv as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB BlackoilPropsAdFromDeck::rvSat(const ADB& po,
const ADB& so,
const Cells& cells) const
{
ADB rv = rvSat(po, cells);
applyVap(rv, so, cells, vap1_);
return rv;
}
// ------ Relative permeability ------
/// Relative permeabilities for all phases.
@ -987,6 +1052,74 @@ namespace Opm
const int n = cells.size();
satprops_->updateSatHyst(n, cells.data(), saturation.data());
}
/// Update for max oil saturation.
void BlackoilPropsAdFromDeck::updateSatOilMax(const std::vector<double>& saturation)
{
if (!satOilMax_.empty()) {
const int n = satOilMax_.size();
const int np = phase_usage_.num_phases;
const int posOil = phase_usage_.phase_pos[Oil];
const double* s = saturation.data();
for (int i=0; i<n; ++i) {
if (satOilMax_[i] < *(s+np*i+posOil)) satOilMax_[i] = *(s+np*i+posOil);
}
}
}
/// Apply correction to rs/rv according to kw VAPPARS
/// \param[in/out] r Array of n rs/rv values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the r and so values.
/// \param[in] vap Correction parameter.
void BlackoilPropsAdFromDeck::applyVap(V& r,
const V& so,
const std::vector<int>& cells,
const double vap) const
{
if (!satOilMax_.empty() && vap > 0.0) {
const int n = cells.size();
V factor = V::Ones(n, 1);
for (int i=0; i<n; ++i) {
if (satOilMax_[cells[i]] > 0.01 && so[i] < satOilMax_[cells[i]]) {
factor[i] = std::pow(so[i]/satOilMax_[cells[i]], vap);
}
}
r = factor*r;
}
}
/// Apply correction to rs/rv according to kw VAPPARS
/// \param[in/out] r Array of n rs/rv values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the r and so values.
/// \param[in] vap Correction parameter.
void BlackoilPropsAdFromDeck::applyVap(ADB& r,
const ADB& so,
const std::vector<int>& cells,
const double vap) const
{
if (!satOilMax_.empty() && vap > 0.0) {
const int n = cells.size();
V factor = V::Ones(n, 1);
//V dfactor_dso = V::Zero(n, 1); TODO: Consider effect of complete jacobian (including so-derivatives)
for (int i=0; i<n; ++i) {
if (satOilMax_[cells[i]] > 0.01 && so.value()[i] < satOilMax_[cells[i]]) {
factor[i] = std::pow(so.value()[i]/satOilMax_[cells[i]], vap);
//dfactor_dso[i] = vap*std::pow(so.value()[i]/satOilMax_[cells[i]], vap-1.0)/satOilMax_[cells[i]];
}
}
//ADB::M dfactor_dso_diag = spdiag(dfactor_dso);
//const int num_blocks = so.numBlocks();
//std::vector<ADB::M> jacs(num_blocks);
//for (int block = 0; block < num_blocks; ++block) {
// jacs[block] = dfactor_dso_diag * so.derivative()[block];
//}
//r = ADB::function(factor, jacs)*r;
r = factor*r;
}
}
} // namespace Opm

View File

@ -277,6 +277,15 @@ namespace Opm
V rsSat(const V& po,
const Cells& cells) const;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V rsSat(const V& po,
const V& so,
const Cells& cells) const;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -284,6 +293,15 @@ namespace Opm
ADB rsSat(const ADB& po,
const Cells& cells) const;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB rsSat(const ADB& po,
const ADB& so,
const Cells& cells) const;
// ------ Rv condensation curve ------
/// Condensation curve for Rv as function of oil pressure.
@ -293,6 +311,15 @@ namespace Opm
V rvSat(const V& po,
const Cells& cells) const;
/// Condensation curve for Rv as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
V rvSat(const V& po,
const V& so,
const Cells& cells) const;
/// Condensation curve for Rv as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -300,6 +327,15 @@ namespace Opm
ADB rvSat(const ADB& po,
const Cells& cells) const;
/// Condensation curve for Rv as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
ADB rvSat(const ADB& po,
const ADB& so,
const Cells& cells) const;
// ------ Relative permeability ------
/// Relative permeabilities for all phases.
@ -344,6 +380,9 @@ namespace Opm
void updateSatHyst(const std::vector<double>& saturation,
const std::vector<int>& cells);
/// Update for max oil saturation.
void updateSatOilMax(const std::vector<double>& saturation);
private:
/// Initializes the properties.
template <class CentroidIterator>
@ -356,6 +395,17 @@ namespace Opm
int dimension,
const bool init_rock);
/// Correction to rs/rv according to kw VAPPARS
void applyVap(V& r,
const V& so,
const std::vector<int>& cells,
const double vap) const;
void applyVap(ADB& r,
const ADB& so,
const std::vector<int>& cells,
const double vap) const;
RockFromDeck rock_;
std::unique_ptr<SaturationPropsInterface> satprops_;
@ -380,6 +430,12 @@ namespace Opm
std::vector<int> pvtTableIdx_;
std::vector<std::array<double, BlackoilPhases::MaxNumPhases> > densities_;
// VAPPARS
double vap1_;
double vap2_;
std::vector<double> satOilMax_;
};
} // namespace Opm

View File

@ -247,6 +247,16 @@ namespace Opm
V rsSat(const V& po,
const Cells& cells) const = 0;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
virtual
V rsSat(const V& po,
const V& so,
const Cells& cells) const = 0;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -255,6 +265,16 @@ namespace Opm
ADB rsSat(const ADB& po,
const Cells& cells) const = 0;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
virtual
ADB rsSat(const ADB& po,
const ADB& so,
const Cells& cells) const = 0;
// ------ Rs bubble point curve ------
/// Bubble point curve for Rs as function of oil pressure.
@ -265,6 +285,16 @@ namespace Opm
V rvSat(const V& po,
const Cells& cells) const = 0;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
virtual
V rvSat(const V& po,
const V& so,
const Cells& cells) const = 0;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
@ -273,6 +303,16 @@ namespace Opm
ADB rvSat(const ADB& po,
const Cells& cells) const = 0;
/// Bubble point curve for Rs as function of oil pressure.
/// \param[in] po Array of n oil pressure values.
/// \param[in] so Array of n oil saturation values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n bubble point values for Rs.
virtual
ADB rvSat(const ADB& po,
const ADB& so,
const Cells& cells) const = 0;
// ------ Relative permeability ------
/// Relative permeabilities for all phases.
@ -321,6 +361,10 @@ namespace Opm
virtual
void updateSatHyst(const std::vector<double>& saturation,
const std::vector<int>& cells) = 0;
/// Update for max oil saturation.
virtual
void updateSatOilMax(const std::vector<double>& saturation) = 0;
};
} // namespace Opm

View File

@ -254,18 +254,22 @@ namespace Opm {
V
fluidRsSat(const V& p,
const V& so,
const std::vector<int>& cells) const;
ADB
fluidRsSat(const ADB& p,
const ADB& so,
const std::vector<int>& cells) const;
V
fluidRvSat(const V& p,
const V& so,
const std::vector<int>& cells) const;
ADB
fluidRvSat(const ADB& p,
const ADB& so,
const std::vector<int>& cells) const;
ADB

View File

@ -495,6 +495,10 @@ namespace {
state.saturation[ pu.phase_pos[ Gas ] ] = sg;
so = so - sg;
std::vector<int> all_cells = buildAllCells(nc);
ADB rsSat = fluidRsSat(state.pressure, so, all_cells);
ADB rvSat = fluidRvSat(state.pressure, so, all_cells);
if (has_disgas_) {
state.rs = (1-isRs) * rsSat + isRs*xvar;
} else {
@ -595,18 +599,20 @@ namespace {
const ADB bw = fluid_.bWat(perf_press, well_cells);
b.col(pu.phase_pos[BlackoilPhases::Aqua]) = bw.value();
}
assert(active_[Oil]);
const ADB perf_so = subset(state.saturation[pu.phase_pos[Oil]], well_cells);
if (pu.phase_used[BlackoilPhases::Liquid]) {
const ADB perf_rs = subset(state.rs, well_cells);
const ADB bo = fluid_.bOil(perf_press, perf_rs, perf_cond, well_cells);
b.col(pu.phase_pos[BlackoilPhases::Liquid]) = bo.value();
const V rssat = fluidRsSat(perf_press.value(), well_cells);
const V rssat = fluidRsSat(perf_press.value(), perf_so.value(), well_cells);
rssat_perf.assign(rssat.data(), rssat.data() + nperf);
}
if (pu.phase_used[BlackoilPhases::Vapour]) {
const ADB perf_rv = subset(state.rv, well_cells);
const ADB bg = fluid_.bGas(perf_press, perf_rv, perf_cond, well_cells);
b.col(pu.phase_pos[BlackoilPhases::Vapour]) = bg.value();
const V rvsat = fluidRvSat(perf_press.value(), well_cells);
const V rvsat = fluidRvSat(perf_press.value(), perf_so.value(), well_cells);
rvsat_perf.assign(rvsat.data(), rvsat.data() + nperf);
}
// b is row major, so can just copy data.
@ -1285,8 +1291,8 @@ namespace {
// phase translation sg <-> rs
const V rsSat0 = fluidRsSat(p_old, cells_);
const V rsSat = fluidRsSat(p, cells_);
const V rsSat0 = fluidRsSat(p_old, s_old.col(pu.phase_pos[Oil]), cells_);
const V rsSat = fluidRsSat(p, so, cells_);
std::fill(primalVariable_.begin(), primalVariable_.end(), PrimalVariables::Sg);
@ -1310,8 +1316,8 @@ namespace {
}
// phase transitions so <-> rv
const V rvSat0 = fluidRvSat(p_old, cells_);
const V rvSat = fluidRvSat(p, cells_);
const V rvSat0 = fluidRvSat(p_old, s_old.col(pu.phase_pos[Oil]), cells_);
const V rvSat = fluidRvSat(p, so, cells_);
if (has_vapoil_) {
// The obvious case
@ -1863,9 +1869,10 @@ namespace {
template<class T>
V
FullyImplicitBlackoilSolver<T>::fluidRsSat(const V& p,
const V& satOil,
const std::vector<int>& cells) const
{
return fluid_.rsSat(p, cells);
return fluid_.rsSat(p, satOil, cells);
}
@ -1875,17 +1882,19 @@ namespace {
template<class T>
ADB
FullyImplicitBlackoilSolver<T>::fluidRsSat(const ADB& p,
const ADB& satOil,
const std::vector<int>& cells) const
{
return fluid_.rsSat(p, cells);
return fluid_.rsSat(p, satOil, cells);
}
template<class T>
V
FullyImplicitBlackoilSolver<T>::fluidRvSat(const V& p,
const V& satOil,
const std::vector<int>& cells) const
{
return fluid_.rvSat(p, cells);
return fluid_.rvSat(p, satOil, cells);
}
@ -1895,9 +1904,10 @@ namespace {
template<class T>
ADB
FullyImplicitBlackoilSolver<T>::fluidRvSat(const ADB& p,
const ADB& satOil,
const std::vector<int>& cells) const
{
return fluid_.rvSat(p, cells);
return fluid_.rvSat(p, satOil, cells);
}

View File

@ -289,6 +289,9 @@ namespace Opm
SimulatorReport sreport;
// Max oil saturation
props_.updateSatOilMax(state.saturation());
// Solve pressure equation.
// if (check_well_controls_) {
// computeFractionalFlow(props_, allcells_,