OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Updated integration of VFP tables. Produces almost identical results to bhp control

Browse Source
This commit is contained in:
babrodtk 2015-08-06 15:11:47 +02:00
parent 16a4580219
commit 7eb94caeba
5 changed files with 91 additions and 79 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

5
opm/autodiff/AutoDiffHelpers.hpp
View File

@ -354,7 +354,7 @@ spdiag(const AutoDiffBlock<double>::V& d)
public:
typedef AutoDiffBlock<Scalar> ADB;
enum CriterionForLeftElement { GreaterEqualZero, GreaterZero, Zero, NotEqualZero, LessZero, LessEqualZero };
enum CriterionForLeftElement { GreaterEqualZero, GreaterZero, Zero, NotEqualZero, LessZero, LessEqualZero, NotNaN };
Selector(const typename ADB::V& selection_basis,
CriterionForLeftElement crit = GreaterEqualZero)
@ -385,6 +385,9 @@ spdiag(const AutoDiffBlock<double>::V& d)
case LessEqualZero:
chooseleft = selection_basis[i] <= 0.0;
break;
case NotNaN:
chooseleft = !isnan(selection_basis[i]);
break;
default:
OPM_THROW(std::logic_error, "No such criterion: " << crit);
}

12
opm/autodiff/BlackoilModelBase_impl.hpp
View File

@ -1347,13 +1347,9 @@ namespace detail {
const int np = wells().number_of_phases;
const int nw = wells().number_of_wells;
const ADB& aqua_adb = subset(state.qs, Span(nw, 1, BlackoilPhases::Aqua*nw));
const ADB& liquid_adb = subset(state.qs, Span(nw, 1, BlackoilPhases::Liquid*nw));
const ADB& vapour_adb = subset(state.qs, Span(nw, 1, BlackoilPhases::Vapour*nw));
const V& aqua = aqua_adb.value();
const V& liquid = liquid_adb.value();
const V& vapour = vapour_adb.value();
const ADB& aqua = subset(state.qs, Span(nw, 1, BlackoilPhases::Aqua*nw));
const ADB& liquid = subset(state.qs, Span(nw, 1, BlackoilPhases::Liquid*nw));
const ADB& vapour = subset(state.qs, Span(nw, 1, BlackoilPhases::Vapour*nw));
V bhp_targets = V::Zero(nw);
V rate_targets = V::Zero(nw);
@ -1379,7 +1375,7 @@ namespace detail {
const double& thp = well_controls_iget_target(wc, current);
const double& alq = well_controls_iget_alq(wc, current);
bhp_targets (w) = vfp_properties_->getProd()->bhp(vfp, aqua[w], liquid[w], vapour[w], thp, alq).value;
bhp_targets (w) = vfp_properties_->getProd()->bhp(vfp, aqua.value()[w], liquid.value()[w], vapour.value()[w], thp, alq).value;
rate_targets(w) = -1e100;
}
break;

32
opm/autodiff/VFPProperties.cpp
View File

@ -123,6 +123,7 @@ VFPProdProperties::ADB VFPProdProperties::bhp(int table_id,
const ADB& vapour,
const ADB& thp,
const ADB& alq) const {
const VFPProdTable* table = getProdTable(table_id);
const int nw = thp.size();
assert(aqua.size() == nw);
@ -140,14 +141,22 @@ VFPProdProperties::ADB VFPProdProperties::bhp(int table_id,
dalq.resize(nw);
dflo.resize(nw);
//Find interpolation variables
ADB flo = getFlo(aqua, liquid, vapour, table->getFloType());
ADB wfr = getWFR(aqua, liquid, vapour, table->getWFRType());
ADB gfr = getGFR(aqua, liquid, vapour, table->getGFRType());
//Compute the BHP for each well independently
for (int i=0; i<nw; ++i) {
adb_like bhp_val = bhp(table_id,
aqua.value()[i],
liquid.value()[i],
vapour.value()[i],
thp.value()[i],
alq.value()[i]);
//First, find the values to interpolate between
auto flo_i = find_interp_data(flo.value()[i], table->getFloAxis());
auto thp_i = find_interp_data(thp.value()[i], table->getTHPAxis());
auto wfr_i = find_interp_data(wfr.value()[i], table->getWFRAxis());
auto gfr_i = find_interp_data(gfr.value()[i], table->getGFRAxis());
auto alq_i = find_interp_data(alq.value()[i], table->getALQAxis());
adb_like bhp_val = interpolate(table->getTable(), flo_i, thp_i, wfr_i, gfr_i, alq_i);
value[i] = bhp_val.value;
dthp[i] = bhp_val.dthp;
dwfr[i] = bhp_val.dwfr;
@ -167,10 +176,13 @@ VFPProdProperties::ADB VFPProdProperties::bhp(int table_id,
const int num_blocks = aqua.numBlocks();
std::vector<ADB::M> jacs(num_blocks);
for (int block = 0; block < num_blocks; ++block) {
fastSparseProduct(dthp_diag, aqua.derivative()[block], jacs[block]);
ADB::M temp;
fastSparseProduct(dwfr_diag, liquid.derivative()[block], temp);
jacs[block] += temp;
//Could have used fastSparseProduct and temporary variables
//but may not save too much on that.
jacs[block] = dthp_diag * thp.derivative()[block] +
dwfr_diag * wfr.derivative()[block] +
dgfr_diag * gfr.derivative()[block] +
dalq_diag * alq.derivative()[block] +
dflo_diag * flo.derivative()[block];
}
ADB retval = ADB::function(std::move(value), std::move(jacs));

41
opm/autodiff/VFPProperties.hpp
View File

@ -23,6 +23,7 @@
#include <opm/core/wells.h>
#include <opm/core/utility/ErrorMacros.hpp>
#include <opm/autodiff/AutoDiffBlock.hpp>
#include <opm/autodiff/AutoDiffHelpers.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/VFPProdTable.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/VFPInjTable.hpp>
#include <boost/multi_array.hpp>
@ -235,25 +236,22 @@ public:
template <typename T>
static T getWFR(const T& aqua, const T& liquid, const T& vapour,
const VFPProdTable::WFR_TYPE& type) {
T retval;
switch(type) {
case VFPProdTable::WFR_WOR:
case VFPProdTable::WFR_WOR: {
//Water-oil ratio = water / oil
retval = aqua / liquid;
break;
T wor = aqua / liquid;
return zeroIfNan(wor);
}
case VFPProdTable::WFR_WCT:
//Water cut = water / (water + oil)
retval = aqua / (aqua + liquid);
break;
return zeroIfNan(aqua / (aqua + liquid));
case VFPProdTable::WFR_WGR:
//Water-gas ratio = water / gas
retval = aqua / vapour;
break;
return zeroIfNan(aqua / vapour);
case VFPProdTable::WFR_INVALID: //Intentional fall-through
default:
OPM_THROW(std::logic_error, "Invalid WFR_TYPE: '" << type << "'");
}
return zeroIfNan(retval);
}
/**
@ -263,25 +261,20 @@ public:
template <typename T>
static T getGFR(const T& aqua, const T& liquid, const T& vapour,
const VFPProdTable::GFR_TYPE& type) {
T retval;
switch(type) {
case VFPProdTable::GFR_GOR:
// Gas-oil ratio = gas / oil
retval = vapour / liquid;
break;
return zeroIfNan(vapour / liquid);
case VFPProdTable::GFR_GLR:
// Gas-liquid ratio = gas / (oil + water)
retval = vapour / (liquid + aqua);
break;
return zeroIfNan(vapour / (liquid + aqua));
case VFPProdTable::GFR_OGR:
// Oil-gas ratio = oil / gas
retval = liquid / vapour;
break;
return zeroIfNan(liquid / vapour);
case VFPProdTable::GFR_INVALID: //Intentional fall-through
default:
OPM_THROW(std::logic_error, "Invalid GFR_TYPE: '" << type << "'");
}
return zeroIfNan(retval);
}
@ -339,11 +332,13 @@ private:
return (std::isnan(value)) ? 0.0 : value;
}
static inline ADB::V zeroIfNan(const ADB::V& value) {
ADB::V retval(value.size());
for (int i=0; i<value.size(); ++i) {
retval[i] = zeroIfNan(value[i]);
}
static inline ADB zeroIfNan(const ADB& values) {
Selector<ADB::V::Scalar> not_nan_selector(values.value(), Selector<ADB::V::Scalar>::NotNaN);
const ADB::V z = ADB::V::Zero(values.size());
const ADB zero = ADB::constant(z, values.blockPattern());
ADB retval = not_nan_selector.select(values, zero);
return retval;
}
};
@ -387,6 +382,6 @@ inline VFPProdProperties::adb_like operator*(
return retval;
}
}
} //Namespace
#endif /* OPM_AUTODIFF_VFPPROPERTIES_HPP_ */

80
tests/test_vfpproperties.cpp
View File

@ -336,7 +336,6 @@ BOOST_AUTO_TEST_CASE(InterpolatePlane)
adb_like max_d;
//Check interpolation
int ii=0;
for (int i=0; i<=n; ++i) {
const double thp = i / static_cast<double>(n);
for (int j=1; j<=n; ++j) {
@ -661,16 +660,23 @@ struct ConversionFixture {
ConversionFixture() :
num_wells(5),
aqua(num_wells),
liquid(num_wells),
vapour(num_wells)
aqua(ADB::null()),
liquid(ADB::null()),
vapour(ADB::null())
{
ADB::V aqua_v(num_wells);
ADB::V liquid_v(num_wells);
ADB::V vapour_v(num_wells);
for (int i=0; i<num_wells; ++i) {
aqua[i] = 300+num_wells*15;
liquid[i] = 500+num_wells*15;
vapour[i] = 700+num_wells*15;
aqua_v[i] = 300+num_wells*15;
liquid_v[i] = 500+num_wells*15;
vapour_v[i] = 700+num_wells*15;
}
aqua = ADB::constant(aqua_v);
liquid = ADB::constant(liquid_v);
vapour = ADB::constant(vapour_v);
}
~ConversionFixture() {
@ -679,9 +685,9 @@ struct ConversionFixture {
int num_wells;
ADB::V aqua;
ADB::V liquid;
ADB::V vapour;
ADB aqua;
ADB liquid;
ADB vapour;
};
@ -698,26 +704,26 @@ BOOST_AUTO_TEST_CASE(getFlo)
std::vector<double> ref_flo_liq(num_wells);
std::vector<double> ref_flo_gas(num_wells);
for (int i=0; i<num_wells; ++i) {
ref_flo_oil[i] = liquid[i];
ref_flo_liq[i] = aqua[i] + liquid[i];
ref_flo_gas[i] = vapour[i];
ref_flo_oil[i] = liquid.value()[i];
ref_flo_liq[i] = aqua.value()[i] + liquid.value()[i];
ref_flo_gas[i] = vapour.value()[i];
}
{
ADB::V flo = Opm::VFPProdProperties::getFlo(aqua, liquid, vapour, Opm::VFPProdTable::FLO_OIL);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getFlo(aqua, liquid, vapour, Opm::VFPProdTable::FLO_OIL);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_flo_oil.begin(), ref_flo_oil.end(), computed, computed+num_wells);
}
{
ADB::V flo = Opm::VFPProdProperties::getFlo(aqua, liquid, vapour, Opm::VFPProdTable::FLO_LIQ);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getFlo(aqua, liquid, vapour, Opm::VFPProdTable::FLO_LIQ);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_flo_liq.begin(), ref_flo_liq.end(), computed, computed+num_wells);
}
{
ADB::V flo = Opm::VFPProdProperties::getFlo(aqua, liquid, vapour, Opm::VFPProdTable::FLO_GAS);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getFlo(aqua, liquid, vapour, Opm::VFPProdTable::FLO_GAS);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_flo_gas.begin(), ref_flo_gas.end(), computed, computed+num_wells);
}
}
@ -730,26 +736,26 @@ BOOST_AUTO_TEST_CASE(getWFR)
std::vector<double> ref_wfr_wct(num_wells);
std::vector<double> ref_wfr_wgr(num_wells);
for (int i=0; i<num_wells; ++i) {
ref_wfr_wor[i] = aqua[i] / liquid[i];
ref_wfr_wct[i] = aqua[i] / (aqua[i] + liquid[i]);
ref_wfr_wgr[i] = aqua[i] / vapour[i];
ref_wfr_wor[i] = aqua.value()[i] / liquid.value()[i];
ref_wfr_wct[i] = aqua.value()[i] / (aqua.value()[i] + liquid.value()[i]);
ref_wfr_wgr[i] = aqua.value()[i] / vapour.value()[i];
}
{
ADB::V flo = Opm::VFPProdProperties::getWFR(aqua, liquid, vapour, Opm::VFPProdTable::WFR_WOR);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getWFR(aqua, liquid, vapour, Opm::VFPProdTable::WFR_WOR);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_wfr_wor.begin(), ref_wfr_wor.end(), computed, computed+num_wells);
}
{
ADB::V flo = Opm::VFPProdProperties::getWFR(aqua, liquid, vapour, Opm::VFPProdTable::WFR_WCT);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getWFR(aqua, liquid, vapour, Opm::VFPProdTable::WFR_WCT);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_wfr_wct.begin(), ref_wfr_wct.end(), computed, computed+num_wells);
}
{
ADB::V flo = Opm::VFPProdProperties::getWFR(aqua, liquid, vapour, Opm::VFPProdTable::WFR_WGR);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getWFR(aqua, liquid, vapour, Opm::VFPProdTable::WFR_WGR);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_wfr_wgr.begin(), ref_wfr_wgr.end(), computed, computed+num_wells);
}
}
@ -762,26 +768,26 @@ BOOST_AUTO_TEST_CASE(getGFR)
std::vector<double> ref_gfr_glr(num_wells);
std::vector<double> ref_gfr_ogr(num_wells);
for (int i=0; i<num_wells; ++i) {
ref_gfr_gor[i] = vapour[i] / liquid[i];
ref_gfr_glr[i] = vapour[i] / (liquid[i] + aqua[i]);
ref_gfr_ogr[i] = liquid[i] / vapour[i];
ref_gfr_gor[i] = vapour.value()[i] / liquid.value()[i];
ref_gfr_glr[i] = vapour.value()[i] / (liquid.value()[i] + aqua.value()[i]);
ref_gfr_ogr[i] = liquid.value()[i] / vapour.value()[i];
}
{
ADB::V flo = Opm::VFPProdProperties::getGFR(aqua, liquid, vapour, Opm::VFPProdTable::GFR_GOR);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getGFR(aqua, liquid, vapour, Opm::VFPProdTable::GFR_GOR);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_gfr_gor.begin(), ref_gfr_gor.end(), computed, computed+num_wells);
}
{
ADB::V flo = Opm::VFPProdProperties::getGFR(aqua, liquid, vapour, Opm::VFPProdTable::GFR_GLR);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getGFR(aqua, liquid, vapour, Opm::VFPProdTable::GFR_GLR);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_gfr_glr.begin(), ref_gfr_glr.end(), computed, computed+num_wells);
}
{
ADB::V flo = Opm::VFPProdProperties::getGFR(aqua, liquid, vapour, Opm::VFPProdTable::GFR_OGR);
const double* computed = &flo[0];
ADB flo = Opm::VFPProdProperties::getGFR(aqua, liquid, vapour, Opm::VFPProdTable::GFR_OGR);
const double* computed = &flo.value()[0];
BOOST_CHECK_EQUAL_COLLECTIONS(ref_gfr_ogr.begin(), ref_gfr_ogr.end(), computed, computed+num_wells);
}
}